Fibre Beam Element

Download as pdf or txt
Download as pdf or txt
You are on page 1of 55

A FiBRE/TIMOSHENKO BEAM ELEMENT

IN
CASTEM 2000
J. Guedes, P. Pegon and A.V. Pinto
Applied Mechanics Unit, Safety Technology Institute,
Joint Research Centre, European Commission,
I-21020 Ispra (VA) Italy

Special Publication Nr. I.94.31


July 1994

The work presented in this report has been realized in the framework of the collaboration contract
num. 5090-92-11 TS ISP F between CEA/DMT/SMTS and JRC/IST/AMU for the development of
the computer code CASTEM 2000.
1 of 55

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

2 of 55

CONTENTS
1 INTRODUCTION 5
2 THE BASIC ASSUMPTIONS OF TIMOSHENKO BEAM THEORY AND
INTRODUCTION OF THE FIBRE MODELLING 7
2.1 Assumed displacement field 7
2.2 Strain-displacement relationships 8
2.3 Virtual work expression 8
2.4 The elastic case 9
2.5 Generalized fibre modelling 10
2.6 Usual fibre modelling 10
2.7 Cross-sectional warping 11
3 SOME TYPICAL CONSTITUTIVE RELATIONS FOR USUAL REINFORCED
CONCRETE FIBRE MODELLING 13
3.1 Concrete Constitutive Laws 13
3.1.1 Monotonic loading 13
3.1.2 Cyclic loading 16
3.2 Steel constitutive law 18
3.2.1 Monotonic loading 18
3.2.2 Cyclic loading 19
4 FINITE ELEMENT IDEALIZATION FOR TIMOSHENKO BEAMS 23
4.1 Displacement and strain representation 23
4.2 Stiffness matrix evaluation 24
4.3 Element stress resultant and nodal resultant 24
4.4 The mass matrix 25
5 THE TIMO ELEMENT WITH FIBRES IN CASTEM 2000 27
5.1 Overview 27
5.2 Short remarks about CASTEM 2000 28
5.3 Global and local referentials for structural beam elements 29
5.4 The two-level approach 29
5.4.1 Description of the section. 29
5.4.2 Description of the beam 30
6 IMPLEMENTATION NOTES 33
6.1 Definition of the elements TRIS and QUAS 33
6.2 Models 33
6.3 Identification subroutines 34
6.4 Complex operators working for the beam (POUTRE) formulation 35
6.4.2 Element elastic stress calculation (SIGM operators) 35
6.4.3 Stiffness calculation (RIGI operator) 36
6.4.4 Mass matrix calculation (MASS operator) 36
6.4.5 Initialisation (ZERO operator) 37
6.4.6 Non linear flow evaluation (ECOU operator) 37
6.4.7 Tangent Stiffness (KTAN operator) 38
A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

3 of 55

6.5 Operators working for the SECTION formulation 38


6.5.1 Biaxial bending of a section (MOCU operator) 38
6.6 How to implement a new fibre model 39
7 NUMERICAL APPLICATIONS 41
7.1 Introduction 41
7.2 Monotonic tests 41
7.2.1 Doubly reinforced rectangular cross sections 41
7.2.2 Influence of the concrete tensile strength 42
7.3 Cyclic tests 43
7.3.1 Uniaxial Loading on a Cantilever T Beam 43
7.3.2 Biaxial Loading on a R/C Square Section 45
8 SUMMARY AND CONCLUSION 49
ACKNOWLEDGEMENTS 49
REFERENCES 50
ANNEX 1- A SAMPLE GIBIANE INPUT FOR FIBRE APPLICATION 51

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

4 of 55

INTRODUCTION

1. INTRODUCTION
Seismic analysis of Reinforced Concrete (R/C) structures requires realistic models able to represent
non-linear stiffness and strength characteristics of members under cyclic loading.
Many analytical models have been proposed in the past for the different structural elements. They
can be classified into two main groups: a) detailed models based on the mechanics of solids and the
consequent description of the local material behaviour (microscopic approach) and b) models based
on a simplified idealization of the overall behaviour (macroscopic approach).
Examples of the second group are the heuristic Takeda-type models able to represent uniaxial bending behaviour of R/C beams. In the first group are included the classical finite element models (continuum element support) and the semi-finite element models (discretization at section level) - fibre
type models (beam element support).
While the finite element models are a powerful tool for simulating the non-linear behaviour of complex structural parts (e.g. joints) their application to the common structural assemblages are impractical due to the large computational and memory requirements. On the other hand, the fibre type
modelling (see Fig. 1), taking advantage of the simplified kinematic hypotheses of the Euler-Bernoulli or Timoshenko theories offers a reliable and practical solution for the non-linear analysis of
composite structural elements such as building columns and walls and bridge piers. Moreover, with
this intermediate modelling approach it is in general possible to overcome the numerical difficulties
associated with the non-linear behaviour of concrete.
This report presents a fibre Timoshenko beam element implemented in the general purpose computer
code CASTEM 2000 [1]. This code allows the user to manipulate meaningful data-structures
(objects) by mean of a high level language (GIBIANE). The consequence of this unusual approach

Reinforced concrete
column

Beam assumption
and fibre modelling

Finite element
discretization

Figure 1- The scope of the fibre modelling

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

5 of 55

INTRODUCTION

(at least in the development and use of large structural analysis codes) for the implementation of
models has been investigated in [2].
The Timoshenko beam theory, assuming that plane sections remain plane after deformation but not
necessarily normal to the beam axial axis, has been considered in order to allow the development of
the so called generalized fibre model in the sense that very complex interactions between normal
force, bending moment and shear can be taken into account. A suitable constitutive model is still to
be derived, however the general formulation and the corresponding numerical data structure have
been created. In the mean time, uniaxial material constitutive relationships (usual fibre modelling)
for concrete and steel have been implemented.
Concrete behaviour is represented by a parabolic curve up to the peak stress point followed by a
straight line in the softening zone. Confinement is taken into account by the modification of the
plane concrete curve and including an additional plateau zone. Cyclic behaviour accounts for stiffness degradation and crack closing phenomena. Tensile resistance has also been considered.
Steel behaviour is represented by a modified Menegoto-Pinto model with a three-stage monotonic
curve (linear, plateau and hardening). Bauschinger effects are taken into account and buckling effects
can be simulated.
It is to be noticed that the present report intends to be a general document able to both describe the
work developed (finite element, constitutive models, etc.) and to assist in the implementation of new
models (implementation manual) in CASTEM 2000. In order to facilitate the access to the document, modelling and implementation aspects have been clearly separated as follows:
section 2, 3 and 4 are dedicated to the Timoshenko beam element and to the fibre modelling
including the material constitutive relations for reinforced concrete elements,
sections 5 and 6 deal with the specific aspects of the implementation in CASTEM 2000. A short
description of the code features is also included,
section 7 presents some numerical applications and Annex 1 one of the GIBIANE input,
finally the main conclusions are drawn and needs for further improvement are identified.

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

6 of 55

THE BASIC ASSUMPTIONS OF TIMOSHENKO BEAM THEORY AND INTRODUCTION OF THE FIBRE MODELLING

2. THE BASIC ASSUMPTIONS OF TIMOSHENKO BEAM


THEORY AND INTRODUCTION OF THE FIBRE MODELLING
2.1 Assumed displacement field
Consider the beam of Fig. 2. It is assumed that planes which are normal to the axis Ox before deformation remain plane, but not necessarily normal to the Ox axis after deformation. This implies that
the displacement field u = ( u x, u y, u z ) can be expressed in terms of the generalized displacement
of any section, that is the displacement U = ( Ux, U y, U z ) of the Ox axis and the rotations
= ( x, y, z ) of the plane.
U ( P ) = U ( R ) + RP
or, using components:
u x ( x, y, z ) = U x ( x ) y z ( x ) + z y ( x )
u y ( x, y, z ) = U y ( x ) z x ( x )

(1)

u z ( x, y, z ) = U z ( x ) + y x ( x )
Note that the planar normal rotations z and y are equal (except for the sign) to the slope of the
Ox axis corrected by the rotations y and z which are due to the transverse shear deformation:
z =

dU y
y
dx

y = z

dU z
dx

(2)

Let U denote the vector of generalized displacement


U = ( U, )

(3)

y
z
x

Figure 2- Timoshenko beam

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

7 of 55

THE BASIC ASSUMPTIONS OF TIMOSHENKO BEAM THEORY AND INTRODUCTION OF THE FIBRE MODELLING

2.2 Strain-displacement relationships


Assuming small deflection theory, the axial and shear components of the strain tensor are obtained
from (1)-(2) as follows:
x =
xy =
xz =

d y
dz
u x
dU x
=
+z
y
dx
dx
x
dx

d x
d x
u x u y
dU y
+
=
z z
= y z
dx
dx
y x
dx

(4)

d x
d x
dU z
u x uz
+
=
+ y + y
= z + y
dx
dx
dx
z x

Evaluating the appropriate partial derivatives, it is found that the three other components y , z and
yz of the strain tensor are apparently null. In fact, no provision was made at the displacement field
level (1) to account for these lateral deformations since the lateral stresses y , z
and yz are
assumed to be null; the values of these deformations can be eliminated from the expression of the
constitutive law.

2.3 Virtual work expression


Consider a Timoshenko beam of length L , of reference axis Ox and section S ( x ) . The beam is subjected to a lateral distributed loading q of components qy and q z . If the beam undergoes a set of virtual displacements u x , u y and u z and rotations x , y and z , the virtual work expression
is:
L

( x x + xy xy + xz xz )dS dx ( Uy qy + Uz qz )dx
0S

= 0

(5)

where the first term takes the following form using (2) and (4):
L

d x Ux + z d x y y d x z x + y xy + z xz + d x x ( zxy + yxz ) dS dx
d

0S

Introducing now the following quantities:


normal force:
Fx =

x dS

(6.1)

shear forces:
Fy =

xy dS
S

Fz =

xz dS

(6.2)

twisting moment:
Mx =

( zxy + yxz )dS

(6.3)

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

8 of 55

THE BASIC ASSUMPTIONS OF TIMOSHENKO BEAM THEORY AND INTRODUCTION OF THE FIBRE MODELLING

bending moments:
My =

zx dS

M z = y x dS

(6.4)

the virtual work takes the expression:


L

Fx d x Ux + Fy y + Fz z + Mx d x x + My d x y + Mz d x z dx
d

(7)

( U y q y + U z q z )dx = 0
0

It is thus possible to associate to the components of the generalized stress vector


F = ( F x, F y, F z, M x, M y, M z ) the components of the generalized strain field E defined by
dU x
d x d y d z
E = ( E x, E y, E z, C x, C y, C z ) =
, y, z,
,
,
dx
dx dx dx

(8)

2.4 The elastic case


For an isotropic elastic material for which E is the Young modulus and the Poisson ratio, the constitutive law takes the simple form:
x
y
z

E
= -------------------------------------( 1 + ) ( 1 2 )

x
1

1 y

1
z

(9.1)

and:
xy
xz
yz

xy
1 0 0
E
= -------------------- 0 1 0 xz
2(1 + )
0 0 1
yz

(9.2)

From the condition y = z = yz = 0 , the following expression are obtained for y , z and yz :
y = z = x

yz = 0

The stress-strain relationships (9.1)-(9.2) then simplify under the form:


x = E x

xy = G xy

xz = G xz

(10)

where G = E 2 ( 1 + ) is the shear modulus


Assume now that these material properties are not homogeneous on the section S and denote by E s
and G s the elastic properties of any fibre of the section. Using (10) and (2)-(4), the components of

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

9 of 55

THE BASIC ASSUMPTIONS OF TIMOSHENKO BEAM THEORY AND INTRODUCTION OF THE FIBRE MODELLING

the generalized stress vector F = ( F x, F y, F z, M x, M y, M z ) , namely the normal force (6.1), the
shear forces (6.2), the twisting moment (6.3) and the bending moments (6.4) take the expression:

Es dS Es zdS
S

Fx

Es zdS Es z

My =

Mz

Es y

Cz
dS
(11)

G s zdS
S

Fy

Ex

dS Es yzdS C y

E s ydS E s yzdS

Gs dS
Mx

Fz

E s ydS

Gs dS

Gs ydS

Gs zdS G s ydS G s ( y + z )dS


2

Ey
Ez
Cx

where no assumption about the position of the intersection R of the reference axis of the beam and
the section (see Fig. 2) nor on the orientation of the local axis Ry and Rz have been made.

2.5 Generalized fibre modelling


The elastic constitutive relation (9.1)-(9.2) may be considered as the reversible part of a more complex constitutive relation of the generic form
= f ( , ... )
In particular, any constitutive law developed for performing refined analysis (e.g. plasticity or damage modelling) can be used at any point of the section. Though it may be questionable to use such a
kind of 3-D constitutive modelling in a structural beam element (and not in a 3-D continuum element), a generalized fibre modelling is sometimes used when very complex interactions between
normal force, bending and transverse shear are expected (see for instance [3]). The Mazars model is
currently used in CASTEM 2000 in order to evaluate the potential of this approach [4].

2.6 Usual fibre modelling


The usual fibre modelling, which accurately accounts only for interaction between normal force and
bending, uses (9.1)-(9.2) as the elastic part of a model in which a complex relationship is introduced
only between the axial component x of the stress field (the fibre stress) and the axial deformation
x (the fibre strain)
x = f ( x, ... )

(12)

In this report we considered only usual fibre modelling. Examples of such modelling for a reinforced
concrete beam are presented in the next section.

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

10 of 55

THE BASIC ASSUMPTIONS OF TIMOSHENKO BEAM THEORY AND INTRODUCTION OF THE FIBRE MODELLING

2.7 Cross-sectional warping


It should be noted that, in order to account for possible cross-sectional warping, the expressions
(6.2)-(6.3) of the shear forces and the twisting moment are changed according to
Fy =

s xy dS
y

s xz dS
z

Fz =

Mx =
y

y
s z xy

(13)

z
s y xz )dS

where s and s are parameters lower or equal to 1.


In elasticity, the expression (11) of the shear forces and the twisting moment is then changed to:

s Gs dS
y

Mx

Fy
Fz

s G s zdS

z
s G s dS

y
s G s zdS

z
s G s ydS

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

z
s G s ydS

z 2
Gs ( s y

y 2
s z )dS

Ey
Ez

(14)

Cx

11 of 55

THE BASIC ASSUMPTIONS OF TIMOSHENKO BEAM THEORY AND INTRODUCTION OF THE FIBRE MODELLING

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

12 of 55

SOME TYPICAL CONSTITUTIVE RELATIONS FOR USUAL REINFORCED CONCRETE FIBRE MODELLING

3. SOME TYPICAL CONSTITUTIVE RELATIONS FOR USUAL


REINFORCED CONCRETE FIBRE MODELLING
3.1 Concrete Constitutive Laws
3.1.1 Monotonic loading
3.1.1.1 Compression stresses
As the concrete is subjected to loading of increasing intensity, it undergoes different phases of damage, from micro-cracking up to ultimate failure. Fig. 3 shows a typical response curve of the concrete
under increasing deformation. After an almost linear zone up to about one-half of the compression
strength, the concrete presents a degradation of the stiffness due to micro-cracking and consequently
nonlinear behaviour is pronounced.
The peak stress zone in the stress-strain curve, where transversal strains increase rapidly owing to
internal cracking parallel to the loading direction, is relatively sharp for high-strength concrete and
has a flat top for low-strength concrete. The strain at the maximum stress is approximately
( = 0.2 %). After the peak, the concrete can still carry on compression stresses although it presents
high decreasing strength ratios for increasing strain.
The evolution of the post-peak curve depends mainly on the degree of confinement (e.g. presence of
stirrups on R/C structural elements). This phenomenon is explained by the triaxial behaviour of the
concrete: due to the Poisson ratio, during an uniaxial compression loading the transversal strains
increase. Although the volume of the specimen starts to decrease, at stresses near the compression
strength, the transversal strains become so high that the volume of the specimen will actually
increase. If transversal bars exist, the expansion of the concrete core is prevented, subjecting the
material to lateral compression forces that improve the ductile performances of the concrete.
In this case, the post-peak curve presents not only a lower decreasing strength ratio but also a residual strength for important strains, reflecting the coupling between the confinement effect and the
interlocking of the aggregate.
Analytical modeling of concrete behaviour can be performed using several approaches ranging from
the classical Plasticity Models to the recent Continuum Damage Mechanics Models. However, due
to the complex behaviour under cyclic loading, a more flexible approach has been adopted herein.
Thus, a two branches law of Hognestad type [5] was considered for the non-linear behaviour of the
concrete under monotonic compression loading. The first branch, a parabolic function, defines the
ascending part of the curve and goes from zero to the maximum compression stress point, ( co, co).
The second branch, a descending straight line, represents the concrete softening behaviour after
maximum strength until failure:

-------- = ------- 2.0 -------


co
co
co

(15)

-------- = 1.0 + Z ( co )
co

(16)

for 0 < < co ,

for co till failure.


The slope of the second branch, Z , depends on the degree of confinement of the concrete: the stirrups, preventing the normal expansion of the concrete core, determine a triaxial compression state in

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

13 of 55

SOME TYPICAL CONSTITUTIVE RELATIONS FOR USUAL REINFORCED CONCRETE FIBRE MODELLING

co
co

Unconfined
(Z)

co

Confined
(Z*)

co*

Figure 3- Confinement effect in the concrete.

the concrete. This phenomenon increases the peak value of the compression strength modifying its
position during the deformation process, as illustrated in Fig. 3. Notice that all the strain and stress
values presented in this report follow the rule: negative sign for compression and positive sign for
tension.
According to the Eurocode specifications [6] and Tassios [7], this effect can be taken into account in
the model through a confinement parameter depending on the section characteristics:
co = co

2
co = co

(17)

0.85
Z = ------------------------------------------------------------------ ( 0.1 + 0.0035 + co )

(18)

in which
= min ( 1 + 2.5 , 1.125 + 1.25 )
where represents the mechanical volumetric ratio of the stirrups equal to:
A s y ( l s )
= -----------------------------------------b c h c co
with A s and l denoting the cross-section of the leg and the total length of the stirrups, s the distance between stirrups along the member axis, b c and h c the dimensions of the confined concrete
core measured from the centre-line of the stirrups and y the yield strength of the stirrups. The
coefficient expresses the effect of both the longitudinal bars and the density of the stirrups on the
degree of confinement of the concrete core:
8
s
s
= 1 ------ 1 -------- 1 --------
3n
2b c
2h c
being n the number of longitudinal bars on the perimeter of the cross-section, engaged in the angle
of a stirrup.
With this parameter it is possible to represent, in a proper way, both the confined and unconfined
concrete, allowing a good description of phenomena like the spalling of concrete.

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

14 of 55

SOME TYPICAL CONSTITUTIVE RELATIONS FOR USUAL REINFORCED CONCRETE FIBRE MODELLING

Tension law

co*

tm

pt
co*

Compression law

Figure 4- Concrete response curves for monotonic loading.

In order to improve the model for confined concrete, a third branch is considered after the softening
compression branch and before failure: a zero slope straight line defining a compression plateau (see
Fig. 4). This additional condition accounts for the residual strength of the concrete core for important
axial post-peak deformations.
The strength for this new branch is a parameter of the model, pt . Park and Priestley[8] suggest
pt = 0.2 co .
3.1.1.2 Tensile stresses
Tensile experimental tests on concrete specimens are not easy to perform and indirect methods may
be considered in measuring the strength of the material, usually less than 20% of the compressive
strength.
For maximum tensile stress, the concrete open cracks and the strength reduces suddenly to almost
zero strength. Because of this non ductile behaviour together with low-strength values, the tensile
strength of the concrete is usually ignored. However, in order to obtain a more realistic response for
structures under earthquake type loading, the tensile behaviour shall be taken into account. The
model herein adopted assumes a tensile stress-strain curve idealized by a straight line, with a slope
equal to the initial compression Young modulus, up to the tensile strength.
When R/C structural elements are considered, a new phenomenon is present due to the bond between
the longitudinal bars and the concrete: the tension-stiffening effect, i.e., the ability of intact concrete between adjacent cracks to carry on tensile stresses, giving to the concrete a fictitious supplementary post-cracking strength (smeared cracking approach).
Following the above considerations, a bilinear stress-strain curve has been adopted. From zero to
maximum tensile strength, point ( t, t) , the model presents a linear elastic behaviour with a slope
equal to the initial compression Young modulus, ( E o = 2.0 ( co co ) ) . The second branch, allowing a softening behaviour after cracking, i.e., after the maximum tensile strength is reached, presents
a straight line from point ( t, t) to a zero stress point at tm , as indicated in Fig. 4.
The constitutive equations are then given by:
= Eo

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

(19)

15 of 55

SOME TYPICAL CONSTITUTIVE RELATIONS FOR USUAL REINFORCED CONCRETE FIBRE MODELLING

for 0 < < t , and


r ( t )
= t -----------------------
r1

(20)

r = tm t

(21)

for t < tm , with

Numerical tests performed on R/C elements show that during the initial loading steps tensile stresses
may assume an important role on the structures global behaviour and that a non-smooth response
curve is obtained if sudden cracking is assumed, ( tm = t ) . The second branch of the curve, representing a post-cracking softening behaviour, not only accounts for the tension-stiffening effect but
also contributes to give better convergence characteristics to the solution algorithm.
3.1.2 Cyclic loading
Under dynamic earthquake type loading conditions, the concrete is submitted to loading and unloading processes. Thus, the non-linear model presented in Section 3.1.1, applicable only to monotonic
stress-strain type loading, must be extended to cyclic loading.
Experimental results show that the envelope for the concrete stress-strain curve under cyclic loading
may be considered unique and identical to the monotonic one. Unloading from the envelope reflects
the stiffness degradation due to monotonic damage and the reloading curve from zero stress tends to
the envelope evidencing some strength degradation. The intersection point between the two curves
defines the so-called common point. Stresses above it produce additional plastic strains, while cycles
below the line defined by the locus of the common points result in stable loops.
A small hysteresis effect is present during the cyclic process, as well as a decrease in strength, i.e.,
for the maximum strain ever reached during the loading history the reloading curve assumes a stress
lower than the one given by the previous unloading branch. Along with this phenomenon, the concrete presents, for both the unloading and the reloading branches, decreasing stiffness values as the
maximum strain increases on the stress-strain curve (see Fig. 5).
Under tensile cyclic loading, the envelope for the concrete may also be considered unique and identical to the monotonic curve. Before reaching the tensile strength the concrete presents an almost linear elastic behaviour. After that point, experimental results on R/C elements show that the unloading
and reloading branches follow a straight line and present increasing stiffness degradation for increasing maximum strains on the envelope.

Common points

Envelope

Figure 5- Concrete response under compression cyclic loading.

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

16 of 55

SOME TYPICAL CONSTITUTIVE RELATIONS FOR USUAL REINFORCED CONCRETE FIBRE MODELLING

No experimental tests have yet been performed subjecting the concrete simultaneously to compressive and tensile stresses under cyclic loading. This lack of information prevents a good knowledge of
the crack closing phenomenon, i.e., the way the concrete re-starts to work in compression when, after
crack opening, there is a reversedl strain loading. Nevertheless, some considerations on this subject
will be presented later (see Section 7.3.1).
The cyclic behaviour law adopted in the present model is based on experimental results. The compression monotonic curve will be the envelope of the concrete behaviour law under compression
cyclic loading.
In order to simplify the model, unloading from the envelope follows a law similar to the one proposed by Mercer [9], i.e., a straight line with a slope depending on the maximum strain ever reached
during the loading history, max :
2

( max co )

-2
E d = E o 1.0 ---------------------------------------------------------------------------
1.0 + ( max co ) + ( max co )

(22)

The decrease of the slope with the increase of the maximum strain (equation (22)) tries to represent,
in a simple way, the degradation of the material stiffness. The reloading compression curve, also a
straight line (see Fig. 6), goes from the zero stress point at the plastic strain, pl , to the last loading
point reached on the envelope. No strength degradation is considered.
Concerning tensile stresses, the model considers an envelope going from zero stress point at co (see
Fig. 7) to the tensile strength point on the monotonic curve, ( t, t) . Loading in tension follows a
straight line so that the maximum stress is reached for ( tr1 = pl + t ) . The second branch, a
straight line defining the post-cracking behaviour, presents a slope such that the intersection with the
zero stress line is achieved at ( tr2 = pl + tm ), as shown in Fig. 7.
Strain tm is also a parameter of the model. Barzegar-Jamshidi and Schnobrich [10] suggest that
when the reinforcing steel intersects the cracks at right angles, an appropriate value for tm is the
yield strain of the longitudinal reinforcement steel bars.
Any unloading from the tensile curve follows a slope equal to the secant modulus of the concrete at
the reversal stress point.
With this model, strength and stiffness degradation of the concrete under tensile stresses is implicitly
taken into account. Unloading after tr2 follows a zero stress path. When the maximum tensile

3
7
9

= f()

Figure 6- Numerical model for the concrete behaviour under cyclic loading.

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

17 of 55

SOME TYPICAL CONSTITUTIVE RELATIONS FOR USUAL REINFORCED CONCRETE FIBRE MODELLING

Envelope

= f()
2

pl

3
4

co

tr1

1
7

tm

tr2
t

tm

Figure 7- Numerical model for the concrete under tensile stresses.

strength is reached, no more tensile stresses can be supported by the concrete in subsequent unloadings from compression.
If no tensile stresses are considered in the model, or its maximum value has already been reached in
previous cycles, the plastic strain is determined by the zero stress point from the compression
unloading curve. Otherwise, if during the unloading no cracks were present in the concrete and at the
end the strain goes beyond tr2 , the plastic strain defined by the compression unloading curve is corrected by t , ( pl = pl + t ). Fig. 7 illustrates a possible strain history starting from an unloading
compression curve.

3.2 Steel constitutive law


3.2.1 Monotonic loading
Typical monotonic stress-strain curves for steel in general come from bars loaded in tension. After a
linear elastic zone, the steel exhibits a yielding plateau, i.e., a point beyond which the stress is almost
constant. In the plateau zone the material suffers important and irreversible changes in its internal
matrix, and the Poisson ratio, generally in the range 0.2 to 0.3 , becomes equal to 0.5 , meaning that
during this stage the steel bars do not change their volume.
The stress defining the plateau, referred to as the yield stress, may present two different values: the
upper and the lower yield stress, the last one being the value usually adopted to represent it. The
existence of the upper value is determined not only by the geometric characteristics of the bar but
also by the speed of the test. The length of the plateau is mainly determined by the strength and internal characteristics of the steel: steel bars with higher strength and carbon levels, usually present a
shorter plateau.
After this stage, the steel presents a strain-hardening range up to the maximum strength, followed by
a strain-softening curve pointing to failure. Either by using a cold worked or a high-strength and carbon steel, the tensile monotonic curve may not present the yielding plateau, going from the linear
elastic zone immediately to the strain-hardening curve.
Although the behaviour of steel bars under compression forces can be considered identical to the tensile one, if the tie spacing is not small enough steel reinforcing bars in R/C elements can evidence
inelastic buckling. Experimental tests [11] show that when the ratio between the length and the diameter of the bar exceeds 5 , this phenomenon occurs. In this case, the yielding plateau and the strainhardening zones may even disappear from the response, being substituted by a softening behaviour
curve immediately after the linear elastic zone.

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

18 of 55

SOME TYPICAL CONSTITUTIVE RELATIONS FOR USUAL REINFORCED CONCRETE FIBRE MODELLING

In terms of modelling, the monotonic tensile stress-strain curve proposed for the reinforcing steel
bars is a five parameters law, including the following three zones: the linear elastic, defined by the
Young modulus, E , and the yield strain, sy ,
= E

(23)

the yielding plateau from the yield strain, sy , up to the hardening strain, sh , and defined by the
yield stress, sy ,
= sy

(24)

and finally the hardening curve up to the ultimate stress-strain point, ( su, su) , represented by a
fourth degree polynomial:
su 4
= su ( su sy ) --------------------
su sh

(25)

The same constitutive law is adopted for compression stresses if no buckling effect is considered for
the reinforcing bars. Otherwise, some modifications must be introduced after the linear elastic zone.
This issue will be referred to later on when presentating the reinforcing steel cyclic behaviour.
3.2.2 Cyclic loading
3.2.2.1 basic model
According to experimental results, steel bars under cyclic loading present an unloading straight line
pattern with a slope equal to the initial Young modulus. If after loading beyond the yield strain
reversed stresses are applied to the bars, the stress-strain response curve becomes nonlinear at a
stress much lower than the initial yield stress. This phenomenon, the so-called Bauschinger effect,
depends on the strain history.
Following the experimental results, the monotonic curve together with the explicit formulation proposed by Giuffr and Pinto, later on implemented by Menegotto [12], were chosen to represent the
behaviour of the reinforcing steel bars under cyclic loading.
If small unloading cycles are imposed to the bars, no hysteretic effect is evidenced and a straight
line, with a slope equal to the elastic Young modulus, can be adopted for both unloading and reloading curves. Under virgin strains, reloading follows the monotonic curve.
Otherwise, if after unloading from a post-yielding position the strain satisfies the relation:
sy
max < --------3.0
the Giuffr curve is activated and an important hysteretic effect may occur. Any subsequent loading
or unloading curve follows this formulation.
The selection of this model to represent the steel cyclic behaviour was based on its simplicity,
numerical efficiency and good agreement with experimental results from cyclic tests on reinforcing
steel bars. The simplicity is well expressed by the fact that a single equation is enough to represent
both loading and unloading stress-strain curves:

(1 b)
= b + ----------------------------------
1

R
R
( 1 + ( ) )

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

(26)

19 of 55

SOME TYPICAL CONSTITUTIVE RELATIONS FOR USUAL REINFORCED CONCRETE FIBRE MODELLING

where
s r
= -------------o r

s r
= ----------------o r

a1
R = R o -------------a2 +

As illustrated in Fig. 8, this equation defines a family of transition curves between two asymptotes
with slopes E o and E h and having ( o, o) as a common point. The pair of values ( r, r) are the
coordinates of the last reversal loading point. On the other hand, the b factor represents the slope
hardening ratio, i.e., the ratio between the hardening slope, E h , and the initial slope, E o , and R is a
parameter defining the shape of the transition branch of the curve. This parameter allows a good representation of the Bauschinger effect and its value depends on , i.e., the difference between the
maximum value of the strain ever reached on the loading direction and o , normalized by ( o r ) .
Numerical tests performed on several structures showed that when the tangent stiffness is considered
in the resolution algorithm the convergence process can be very sensitive to the adopted slope hardening ratio. A small b factor, together with a softening behaviour in some of the concrete fibres, may
introduce singularities in the global stiffness matrix.
Nevertheless, the hardening slope can be estimate by:
su sy
E h = -------------------- su sy

(27)

Parameters a 1 and a 2 , together with R o , should be obtained from experimental results. For all the
numerical tests performed, the values suggested by Menegotto [12] were adopted:
R o = 20.0

a 1 = 18.5

a 2 = 0.15

To start the cyclic loading model, the stress and strain values for both reversal points are initialized
by the yield parameters.

2x(o-r)2
(o,o)2

(r,r)1

(o,o)

Eh

= f()

R(2)

Eo

Eo

R(1)

(r,r)2

Eh
(o,o)1

(r,r)
1x(o-r)1

Figure 8- Numerical model for the steel under cyclic loading.

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

20 of 55

SOME TYPICAL CONSTITUTIVE RELATIONS FOR USUAL REINFORCED CONCRETE FIBRE MODELLING

3.2.2.2 Inelastic buckling


The proposed model also accounts for the inelastic buckling effect. Tests made by Monti and Nuti
[11] on several steel bars showed that, for ratios between the length L and the diameter D less or
equal to 5 , the compression curve is similar to the tensile one, i.e., no buckling effect is observed.
According to the same author, if L D ratios are greater than 5 , modifications should be introduced
in the model, namely the unloading curve from tensile stresses and the reloading curve after reversal
from compression (see Fig. 9). In this case the post-yield compression zone has a softening behaviour and the new b factor is computed by:

b c = a ( 5.0 L D )e

Eo
b------------------ sy

(28)

sy
in which = 4.0 ----------- and is the major difference between the reversal strain values found in
LD
both directions during the loading history.
On the other hand, during reloading after reversal from compression, a reduced Young modulus is
adopted for the first asymptote:
E r = E ( a 5 + ( 1.0 a 5 )e

( a 6 )
2

(29)

in which a 5 = 1.0 + ( 5.0 L D ) 7.5 and = ( o r ) , and consequently a new b factor for
tensile stresses must be adopted,
b t = ( bE ) E r

(30)

If 5 is the strain at which the compression curve diverges more than 5% from the tensile one
towards lower values, an additional stress given by
b bc
s = s E o -----------------1.0 b c

(31)

= f()
E

s
bcxE

bxE
Er

Figure 9- Inelastic buckling on the steel cyclic behaviour.

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

21 of 55

SOME TYPICAL CONSTITUTIVE RELATIONS FOR USUAL REINFORCED CONCRETE FIBRE MODELLING

where
11.0 L D
s = o 5 = ----------------------------------------c(L D)
10 ( e
1.0 )

(32)

must be added to the yield value so that the new hardening asymptote could be correctly positioned. This shift allows the response curve to follow first the unbuckled curve and then the buckling
compression branch (see Fig. 9).
Expressions (27) to (31) are valid for ( 5.0 < L D 11.0 ) . Parameters a , c and a 6 are experimental
values to be defined for each steel grade. In all the numerical tests presented in this work the values
suggested by Monti and Nuti [11] were adopted:
a = 0.006

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

c = 0.500

a 6 = 620.0

22 of 55

FINITE ELEMENT IDEALIZATION FOR TIMOSHENKO BEAMS

4. FINITE ELEMENT IDEALIZATION FOR TIMOSHENKO BEAMS


4.1 Displacement and strain representation
In the linear beam element, the displacement field U = ( Ux, U y, U z ) is approximated by the relationship:
1

U ( x ) = N ( x )U + N ( x )U
1

(33)
2

where U and U are the nodal displacements and N and N are the conventional linear shape
functions:
x2 x
1
N ( x ) = ------------l

x x1
2
N ( x ) = ------------l

(34)

and l = x 2 x 1 is the length of the element.


Similarly, the rotations = ( x, y, z ) are represented by:
1

( x ) = N ( x ) + N ( x )
1

(35)

where and are the nodal rotations.


1

Denoting by the set of nodal variables ( U , , U , ) and using (2), (8), (33), (34) and (35)
the expressions for the generalized strains in the element can be readily evaluated:

Ex =
Ey =

dUx
1 1 1 2
s
= --- U x + --- Ux = B x
l
l
dx

dU y
1 1 1 1 1 2 1 2
s
z = --- U y --- z + --- U y --- z = B y
l
2
l
2
dx

(36)

dU z
1 1 1 1 1 2 1 2
s
Ez =
+ y = --- U z + --- + --- U z + --- y = B z
2
l
2 y l
dx
C =

d
1 1 1 2
b
= --- + --- = B
l
l
dx

= x, y, z

In fact, the expressions of the generalized transverse shear strains, say E y and E z have been modified (elimination of the linear terms according to Donea and Lamain[13]), in order to a priori remove
the problem of shear locking which occurs for equal approximation elements. A discussion of this
problem has been given in a previous report dealing with the usual linear Timoshenko element [14].
The equations (36) take the following symbolic form
E = B

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

(37)

23 of 55

FINITE ELEMENT IDEALIZATION FOR TIMOSHENKO BEAMS

4.2 Stiffness matrix evaluation


In the elastic case, the introduction in the virtual work principle (7) of the expressions (11) and (14)
relating the generalized stresses and the generalized strains (8) leads to:
L

E
K
Edx

(
U
q
+
U
q
)dx
=
0
y y
z z

(38)

where the section constitutive matrix K is given by:


K 11 =

Es dS

K 15 =

K 22 =

Es zdS

K 16 = E s ydS

s Gs dS

K 24 = s G s zdS

K 33 =

s Gs dS
z

K 34 =

K 42 =

y
s G s zdS

s Gs ydS
z

K 43 =

z
s G s ydS

K 51 =

Es zdS

(39)

K 44 =

z 2
ES ( s y

K 55 =

Es z

dS

K 56 = E s yzdS

K61 = E s ydS

y 2
s z )dS

K 65 = E s yzdS

K 66 =

Es y

dS

the other elements being null.


Using now a discretization of the domain [ 0, L ] in elements and introducing the linear approximation of the type (33), (34) and (35), (38) is equivalent to the solution of a linear system.
K = f

(40)

The contribution K of the element to the global stiffness K takes as usual the following form:
x2
e

K =

[B]

K B dx

(41)

x1

which can be evaluated exactly using a one-point Gauss-Legendre rule.


Integrating now the second term in (38), the consistent nodal force vector for the element is derived:
qy l qz l
qy l qz l
e
f = 0, -------, -------, 0, 0, 0, 0, -------, -------, 0, 0, 0
2 2

2 2

(42)

which, unlike the Euler-Bernouilli cubic Hermitian element, only has lateral nodal point forces and
no nodal moments.

4.3 Element stress resultant and nodal resultant


Using the expressions (36) of the generalized strain and, in the elastic case, the relations (11) and
(14), it is easy to observe that the axial force F x and the moments M x , M y and M z are constant. The

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

24 of 55

FINITE ELEMENT IDEALIZATION FOR TIMOSHENKO BEAMS

transverse shear forces F y and F z are also constant in the element according to the expression (36)
of the transverse generalized strains E y and E z .
Knowing this generalized stress state ( F x, F y, F z, M x, M y, M z ) it is possible to determine the contri1
1
2
2
bution of the element to the nodal resultant ( F , M , F , M ) of the element, using again the virtual
work principle (7), the discretization (36) of the generalized strain and a one-point Gauss-Legendre
rule:
1

F = F

= x, y, z

F = F

Mx = Mx

Mx = Mx

l
1
M y = --- F z M y
2
l
1
M z = --- F y M z
2

(43)

l
2
M y = --- F z + M y
2
l
2
M z = --- F y + M z
2

4.4 The mass matrix


Since the Timoshenko element is to be used also for performing dynamic analysis, it may be interesting to derive its element mass matrix.
The virtual work W I of the inertial forces which have to be added to (5) takes the form:
2

d ux
d uy
d u z

W I = u x 2 + u y 2 + u z 2 dSdx

dt
dt
dt
S

(44)

using the expression (1) of the assumed displacement field and introducing the generalized displacement vector U defined by (3), W I in (44) may be rewritten under the form:
L

WI =

dU
dt

(45)

dx

where the section mass matrix M is given by:


M 11 =

s dS

M15 =

M 22 =

s zdS
S

s s dS

M 24 = s s zdS

M 33 =

s s dS
z

M 34 =

M 42 =

y
s s zdS

M 51 =

M 43 =

z
s s ydS

M 55 =

M 61 = s ydS

s s ydS
z

s zdS

M 16 = s ydS

M 44 =

z 2
S ( s y

s z

dS

M 56 = s yzdS

M 65 = s yzdS

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

(46)
y 2
s z )dS

M 66 =

s y

dS

25 of 55

FINITE ELEMENT IDEALIZATION FOR TIMOSHENKO BEAMS

Introducing in (45) the discretization (33)-(34)-(35) leads to the addition of the linear term
2
2
M ( d dt ) to the first member of the discrete equilibrium equation (40). Performing exact integration for all terms:
x2

x2

x2

l
( N N )dx = ( N N )dx = --3-

(N

x1

x1

x1

2
l
N )dx = --6

the element mass matrix M may be readily derived.

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

26 of 55

THE TIMO ELEMENT WITH FIBRES IN CASTEM 2000

5. THE TIMO ELEMENT WITH FIBRES IN CASTEM 2000


5.1 Overview
The implementation of the TIMOshenko element with fibres could be realized following a 2-D
multi-layer approach, for which selected integration points are used for performing the integration
over the section (or part of the section defined by a material). This approach would generalize the
one which has been introduced in CASTEM 2000 for treating thick shells. Though being natural at a
first glance, this approach suffers from severe drawbacks:
the difficulty to account for the a priory unknown complexity of the section that have to be
described,
the resulting cryptic character of any solution aiming to have some generality,
the difficulty to post-treat the relevant informations at the level of the section.
In order to overcome these difficulties, a rather different approach is adopted here. In fact, the sections are described by a model with various sub-zones corresponding to various materials (say steel
fibre, concrete fibre,...) and associated with the meshes of these various components. This description allows in turn to define the beam modelling. The description of the sections (shape, material
characteristics) is then general, easy to generate and to post-treat.
The evaluation of the stress resultant for each beam element proceeds as follow (see illustration in
Fig. 10):

( U, ) E

F ( F, M )

BEAM level

SECTION level

Figure 10- The fibre element

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

27 of 55

THE TIMO ELEMENT WITH FIBRES IN CASTEM 2000

Evaluation of the generalized strain E at the integration point of each beam element, from the
1
1
2
2
nodal generalized displacement ( U , , U , ) .
Use of the beam model in order to evaluate the strain tensor , and in particular its normal component x at the level of each fibre, located at the Gauss integration points of the elements describing the section.
Use of the constitutive relationship in order to evaluate the stress tensor at the level of each
fibre, and in particular its normal component x .
Integration over the section of the relevant stress components in order to compute the generalized
stress F for the section (the use of Gauss numerical integration rules allows to obtain precise evaluations for these quantities).
1

Computation of the stress resultant ( F , M , F , M ) for the beam element.

5.2 Short remarks about CASTEM 2000


In order to understand better the next sections, it may be interesting to present an outline of the modelling concepts of CASTEM 2000. In-depth descriptions may be found in [1] or [2].
CASTEM 2000 uses an internal language (GIBIANE) which allows the user to create, manipulate
and destroy objects belonging to a relevant class by means of operators. The implementation of the
object classes is realized by means of generic data-structure types, the implementation of the operators by means of subroutines.
The geometry of a problem is specified creating objects of the type mesh. The corresponding datastructure is subdivided in zones corresponding to an identical parent geometric element. A geometric
element is different from a finite element since the same geometry may be associated for example
withcontinuum or structural elements. It is said that the geometrical element does not have the same
formulation.
The modelling of a mechanical problem is specified creating objects of the type model. The corresponding data-structure is subdivided in zones corresponding to identical modelling and identical
parent finite elements. They are called the sub-zones of the model.
CASTEM 2000 being a displacement based finite-element code, the quantities important for the
modelling (say material characteristics, strain, stress, internal variables, etc...) have to be known at
the integration points of the mesh. These informations are hold by objects of the type field-by-element. These objects have specific types (characteristics, strain, stress, internal variable, etc...)
depending on their content. Since almost all the operators creating an object of the type field-by-element involve a model, a field-by-element is subdivided in sub-zones according to a model.
The operations involving a model and the related field-by-elements are generally implemented by
means of 3 nested loops: the first loop on the sub-zone of the specified model, the second on the elements of the zone and the third on the integration points of the element.
In terms of the data-structure, a field-by-element on each sub-zone is a collection of two-entry arrays
containing the values of the field, associated to each components of the field. The two entries are the
number of elements in the sub-zone and the number of integration points. The values are either real
numbers or more complex typed data-structures. The number and the name of the components
depends on the specific type of the field-by-element, on the name of the model and on the formulation of the finite element. For example, consider a field-by-element of the specific type stress in a
sub-zone of parent finite element QUA8 (continuum 2-D plane strain element corresponding to a
quadratic 8-node quadrilateral geometric element with 3x3 Gauss-integration points). It is a collection of two-entry arrays corresponding to the active components xx (component SMXX),
yy (component SMYY), zz (component SMZZ) and xy (component SMXY) of the stress

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

28 of 55

THE TIMO ELEMENT WITH FIBRES IN CASTEM 2000

tensor. For the same geometrical sub-zone, the situation is quite different when changing the formulation of the element. Consider for instance the finite element COQ8, associated with the same 8node geometric element but corresponding now to a 3-D structural thick shell formulation: not only
the number (and the name) of the active components are different (5 components for the stress tensor
in the local axes), but the integration points are arranged in such a way that an integration over the
thickness of the shell is also allowed.
To create the field-by-elements and to further check the internal data consistency, CASTEM 2000
uses internal identification routines, one for each specific field-by-element type: these routines
return, for a given model and a given finite-element formulation, the number and the name of the
components as well as the type of the value.
CASTEM 2000 was developed using the ESOPE language. ESOPE [15] is in fact an extension of
FORTRAN 77, allowing the manipulation (creation, re-sizing, activation, inactivation, destruction)
of data-structures.

5.3 Global and local referentials for structural beam elements


In order to be able to compute the equilibrium of a frame structure, CASTEM 2000, as many codes
do, distinguishes between a global referential ( OXYZ ) in which the equilibrium (or the stiffness or
the mass) is assembled and a local referential ( Oxyz ) in which the computations are mainly performed. A local referential is associated to each element in the following way:
O position of the first node of the element;
Ox axis passing through the two nodes of the element and oriented from the first to the second
node.
Oy axis orthogonal to Ox , in the plane Oxv where v is a vector supplied by the user (default
OY ) and verifying Oy v > 0 .
Oz axis defined by Oz = Ox Oy , for which ( Oxyz ) is a direct orthogonal referential.
Then for all computations at element level, CASTEM 2000 computes the local referential ( Oxyz ) ,
eventually projects global nodal information (say the generalized displacement field UX, UY,
UZ, RX, RY, RZ for U x , U y , U z , x , y , z ) expressed in ( OXYZ ) ) in the local referential, performs element computations in the local referential and eventually projects back the results in
the global one before assembling.
The same approach is of course used for the TIMO element with fibres so that many functions used
in the implementation of the usual Timoshenko element (TIMO) or the Hermitian element
(POUT) can be reused.

5.4 The two-level approach


As outlined in the previous overview section, a two-level approach which uses the existing data
structures and operators of CASTEM 2000 in order to specify the behaviour of each individual section and in order to post-treat them has been adopted: the first level is the description of the section,
the second is the description of the beam.
5.4.1 Description of the section.
Each section is described by a object of the type model, associated to mesh(es) of simple 2-D linear
elements (the 3-node triangle TRIS and the 4-node quadrilateral QUAS) associated with a new
formulation called SECTION. The local referential yOz of the section is assumed to coincide with
the global referential XOY of the mesh. This assumption allows to generate the mesh of the various
sections using the 2-D options of CASTEM 2000, and to switch to the 3-D level only when the beam

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

29 of 55

THE TIMO ELEMENT WITH FIBRES IN CASTEM 2000

mesh is to be generated. Note also that the degree of eccentricity of a given section is simply controlled by the position of the section with respect to the origin of the XOY reference.
The fibres state (e.g. the local stress ) is described by standard objects of the type field-by-element associated with the object of the type model.
A new formulation and accordingly new finite elements are required for at least two reasons:
The number of active components of the stress tensor ( x , xy and xz ) on the section reduces to
3.
the fibres being located at the Gauss integration points of the elements, the useful integrations
over the section, say for deriving the generalized stress state or the section stiffness matrix, are
performed numerically looping over the elements of each sub-zone of the section model and using
the standard Gauss integration rules. Due to the nature of the expression to be integrated (see e.g.
(39)), the full integration (that is, an integration which gives exact results for an elastic case),
should be performed with a number of integration points which may differ from the number used
for instance in the continuum element associated with the same geometric element. This is in particular the case for the TRIS element which requires more than the one integration point in consistent use for the TRI3 (continuum) element.
Only two new elements are useful for the following two reasons: first it is possible to mesh any 2-D
section by mean of triangular or quadrilateral elements, and second neither nodal quantities nor
shape functions being really used, unless for computing the coordinates of the integration points,
there is no need to introduce any element of order higher than one (linear element).
The effective specification of the modelling of the section is completed after the specification of the
material and geometrical characteristics (object of the type field-by-element of the specific type
characteristics). The material characteristics are the usual elastic parameters ( E , and ) and
eventually the non-linear parameters of the constitutive laws of the fibres. The geometrical characteristics are the two cross-sectional warping parameters of the fibre.
In order to check the physical consistency of a given section, an operator called MOCU has been
developed. This operator uses as input the model and the characteristics of the section which defines
a beam, and three objects of the type list-of-real-numbers of same length associated with two
sequences of curvatures C y and Cz and one of normal force F x . These sequences define a bi-axial
loading path for the considered beam. The main output of MOCU are the related two sequences of
bending moments M y and, M z , and the sequence of normal deformation E x .
5.4.2 Description of the beam
A frame of beams with complex sections is described by an object of the type model, associated to
a mesh of 3-D linear structural beam element TIMO[14]. This element, associated with the beam
(POUTRE) formulation, involves one integration point corresponding to the mid-section of the element.
The field-by-elements relevant for the modelling, namely the material or geometric characteristics,
strain, stress and internal variables field are, as usual, a collection of arrays associated to each component of the each field. These arrays contain information relative to the section of integration.
Instead of being simple objects belonging to the real-number, vector, curve or list-of-real-numbers classes, the content of the array may be the complex objects of the type model or field-byelement which effectively describe the section of integration.
The resort to these more complex informations are useful for the TIMO element with fibres:
the material characteristics of each element of a sub-zone are the model (component MODS)
and the field-by-element of characteristics (component MATS) of the section of integration,

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

30 of 55

THE TIMO ELEMENT WITH FIBRES IN CASTEM 2000

the internal variables of each TIMO element are the objects of the type field-by-element and of
the specific type stress (component CONS) and internal variables (component VAIS) on
the section.
These informations are extensively used when a sub-zone of the beam model is non-linear. In order
to limit the computational effort in the case of elastic sections, an optional component have been
added to the field-by-element of material characteristic: the component name is MAHO and the
associated value is a list-of-real-number containing the section constitutive matrix K (39) (Hook
matrix). This matrix can be computed using the HOOK operator on an equivalent beam element
and by extracting from the resulting field-by-element of the specific type Hook-matrix the value
of the component MAHO at the integration point of this element.
It would have been possible to address the problem of the elastic section in a maybe more general
way by introducing a new beam element and consequently a new formulation. This solution has not
been retained in order to limit the duplication of coding effort between this new element (and its
associated formulation) and the already existing TIMO element with beam (POUTRE) formulation.

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

31 of 55

THE TIMO ELEMENT WITH FIBRES IN CASTEM 2000

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

32 of 55

IMPLEMENTATION NOTES

6. IMPLEMENTATION NOTES
6.1 Definition of the elements TRIS and QUAS
The name of the elements TRIS and QUAS with the name of the related formulation SECTION
are defined in BDATA. TRIS and QUAS are declared in the NOMTP array, whereas SECTION is
declared on the array NOMFR.
The new relations element/formulation and element/geometry are introduced in the two subroutines
NUMMFR and NUMGEO. The value of the index MELE associated with the element TRIS and
QUAS are respectively 104 and 105 and the value of the index MFR associated with the SECTION formulation is 47.
The other informations relative to the elements (e.g. number of integration points, associated geometric element) are introduced, as usual, in the subroutine ELQUOI. The array INFELE which is
constructed by this subroutine is consulted before every computation at element level.
Recall that for the TIMO element: MELE=84 and MFR=7.

6.2 Models
For the SECTION formulation, the models are assumed to be at least elastic isotropic (ELASTIQUE ISOTROPE). The non linear models for the section are introduced by the word PLASTIQUE followed by the name of the non-linear fibre model. The two models already presented (see
Section 3) are available under the name ACIER_UNI for the steel and BETON_UNI for the concrete.
For the beam (POUTRE) formulation, the models are assumed to be at least elastic section
(ELASTIQUE SECTION). Only one non linear model is obviously available: its name is PLASTIQUE SECTION.
At the implementation level, the highest entry point for the modelling, directly associated with the
operator MODE, is the subroutine MODELI. This routine would be informed of all the default in
use. MODELI calls the subroutine MODEL2 which is responsible for the mechanical (MECANIQUE) modelling in use here. The names of the new elements that can be associated with such
modelling should be added here. This routine checks for the existence of elastic modelling using the
subroutine MODELA and then for a non-linear modelling using the subroutine MODNLI, which recognizes the PLASTIQUE behaviour, and eventually branches to the routine responsible for such a
behaviour, MODPLA, which analyzes the possible names of the PLASTIQUE models.
At the end of MODEL2 the name of the elements TRIS and QUAS has been added in the array
LESTEF,
In MODELA a new slot has been added for holding the elastic SECTION behaviour (relative
index IELAS=7).
In MODPLA 3 new slots have been added for the 3 new plastic constitutive behaviours:
BETON_UNI (relative index INPLAS=20), ACIER_UNI (INPLAS=21) and SECTION
(INPLAS=22).
According to these modifications/additions, the manual of the operator MODE has been also
changed.
The routine NOMATE associates to each sub-zone of a given model a character string CMATE, characteristics of the elastic part (say isotropic, anisotropic or... SECTION), and an absolute index
INATU for the non linear material. This routine has been expanded in order to account for the previous modelling: at the SECTION level CMATE=ISOTROPE (isotropic) and INATU=0 (elastic

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

33 of 55

IMPLEMENTATION NOTES

case), INATU=39 (BETON_UNI) or INATU=40 (ACIER_UNI); at the beam level


CMATE=SECTION and INATU=0 (elastic case) or INATU=41 (SECTION).
Finally, to perform the computations at the SECTION level, it is interesting to define a relative
index NUFIBR ordering the fibre material. A new subroutine TEMANF associates to the absolute
index INATU the relative index NUFIBR: NUFIBR=1 for BETON_UNI, NUFIBR=2 for
ACIER_UNI and NUFIBR=0 in the other cases and in particular the elastic one.

6.3 Identification subroutines


As usual, when a new formulation is introduced, the number, the name and the type of values of the
components of the field by elements that can be associated with this formulation should be declared
in the related identification routines.
For the SECTION formulation, three specific field-by-elements are required for operational reasons: CHARACTERISTICS for the material and geometrical characteristics, and CONTRAINTE
(stress) and VARIABLES INTERNES (internal variables) for the non-linear sections.
The corresponding modifications/additions have been performed:
Material characteristics: the usual elastic material characteristics (YOUN= E , NU= ,
RHO= ) are specified by IDMATR, whereas the non-linear material characteristics of the two
available plastic models (ACIER_UNI and BETON_UNI) are specified as usual in the satellite of IDMATR called IDPLAS.
The name of the non linear components for the BETON_UNI fibres are (see Section 3.1):
STFC (maximum compression stress co ), EZER (strain at the maximum compression stress
co ), STFT (maximum tensile stress t ), ALF1 ( factor), OME1 ( factor), ZETA
(modulus Z of the descending part of the compression curve), ST85 (stress value pt defining
the compression plateau), TRAF (descending traction curve factor r = tm t ).
The name of the non linear components for the ACIER_UNI fibres are (see Section 3.2):
STSY (yield stress sy ), STSU (ultimate stress su ), EPSH (hardening strain sh ), EPSU
(ultimate strain su ), ROFA ( R o factor), BFAC ( b factor), A1FA ( a 1 factor), A2FA ( a 2
factor), FALD ( L D factor), A6FA ( a 6 factor), CFAC ( c factor), AFAC ( a factor).
y

Geometric characteristics: the names of the two cross-sectional warpings s (ALPY) and s
(ALPZ) are introduced in IDCARB.
Stress: the names of the active stress components x (SMXX), xy (SMXY) and xz
(SMXZ) are introduced in IDCONT.
Internal variables: the names of the internal variables of the two fibre models are introduced in
IDVARI: 6 internal variables for BETON_UNI and 12 for ACIER_UNI. In both models two
internal variables are TANG (current uniaxial stiffness) and EPSO (current uniaxial strain).
The values of all these fields belong to the real-numbers class.
For the beam (POUTRE) formulation some new features are introduced for the following fieldby-elements:
Material characteristics: The material characteristics are defined in the subroutine IDMATR.
They are the section model (MODS, type model), the section characteristics (MATF, type
field-by-element) and, as an optional component, the section constitutive matrix K (39)
(MAHO, type list-of-real-number)
Internal variables: the internal variables for a non linear section (IDVARI) are the stress state of
the section (CONS, type field-by-element) and the associated internal variables (VAIS, type
field-by-element).

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

34 of 55

IMPLEMENTATION NOTES

The manual of the operator MATE has been changed according to the new possible characteristics.

6.4 Complex operators working for the beam (POUTRE) formulation


The simple operators working for the beam (POUTRE) formulation and the TIMO element
with fibres are those ones which do not use the model of the section. These operators does not have
to be modified when fibre modelling is introduced. This is for example the case for the EPSI operator (which computes the field-by-element of strain) and the operator BSIG (which computes the
field-by-node of internal nodal force).
Of course, all the operators using the information at the SECTION level should be modified: this is
the case of the operators HOOK, SIGI, RIGI and MASS used for performing elastic and/or
dynamic analysis and of the operators ZERO, ECOU and KTAN used for performing non-linear
analysis. All these complex operators involve the particular field-by-element defined at the section
level, introduced for performing the non-linear SECTION computations.
6.4.1 Hook matrix calculation (HOOK operator)
The HOOK operator computes the field-by-element of (elastic) Hook matrix at the Gauss points of
the mesh associated with the sub-zones of a given model with associated material characteristics.The
Hook matrix which is constant for each element having the same material and geometric characteristic, is, for the TIMO element with fibres, the 6x6 section constitutive matrix K whose expression is
given by (39).
Due to the particular role played by this operator for building-up efficient elastic calculations (see
Section 5.4.2) it is assumed that the component MAHO of the field-by-element of material characteristics is not (yet) existing. As a consequence the HOOK operator works at both beam and
SECTION level.
As usual in this type of computation, a driver routine (here HOOK2P) performs a loop on each subzone of the input model, collecting and checking all the useful informations for each sub-zone and
then, inside the loop, calls the routine (here HOOK2D) which achieves the work at the element level.
HOOK2P has been modified in order to allow the TIMO element to have complex material characteristics (IF-THEN-ELSE structure on the value of CMATE). The subroutine HOOK2D has been
modified using the same structure, and, when required, calls the new routine FRIGIE which returns
two arrays of 12 integrals allowing to compute the Hook matrix and the section mass matrix.
In fact, FRIGIE is a driver routine at the SECTION level, which performs a loop on each sub-zone
of the section model and which calls a new routine FRIGI2.FRIGI2 loops on the elements/Gauss
points of the sub-zone, and computes the numerical contribution to the integrals which appear in the
expressions of K (39) and M (46).
When available at the beam level, the array relative to the Hook matrix is elaborated (still in
HOOK2D) by a new routine: DOHTIF.
6.4.2 Element elastic stress calculation (SIGM operators)
Three possibilities exist in CASTEM 2000 for evaluating the field-by-element of (local generalized)
stress at the Gauss point of a given mesh. These possibilities have to be compatible.
Compute the field-by-element of the (elastic) Hook matrix (HOOK operator) at the Gauss points
of the mesh, compute the deformation using the EPSI operator and make the product of the two
results.
Use the operator SIGM with the field-by-element of Hook matrix and the field-by-node of (generalized) displacement.

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

35 of 55

IMPLEMENTATION NOTES

Use the SIGM operator with the field-by-element of geometrical and material characteristic and
the field-by-node of (generalized) displacement.
Note that the last syntax is the most widely used for instance when performing the elastic predictor
step in the non-linear model computations. However, even in this case, the Hook matrix (which is
not known as input) is internally computed, which means, for the TIMO element with fibres and as
shown in the previous subsection, the need for a work at the SECTION level. It was then important
to design a possibility to escape this call and this was the reason why the MAHO component was
optionally allowed in the field-by-element of material characteristics.
At the beam level the driver routine HOOK2P has been modified for the TIMO element with fibres in
order to account for complex material characteristics. The driver calls HOOK2D. When the first
syntax of SIGM is used or when the value of the component MAHO is available, the Hook matrix
is directly loaded. In the other case, a work at the section level is performed (call to FRIGIE and
FRIGI2) and the matrix is loaded using DOHTIF.
Since in contrast with the usual (isotropic) TIMO element, the TIMO element with fibres leads
generally to an Hook matrix which is not diagonal, two routines TIFSTR (driving the transformation
global/local for the displacement) and TIFOLO (computing the local element stress resultant) have
been added. These routines are called in SIGMA3 when CMATE=SECTION.
Recall that the components of the field-by-element of generalized stress have the following name:
EFFX= F x , EFFY= Fy , EFFZ= F z , MOMX= M x , MOMY= M y and MOMZ= M z .
6.4.3 Stiffness calculation (RIGI operator)
The operator RIGI calculates the elastic contribution of all the elements of a mesh to the global
stiffness matrix.
As for the SIGM operator, RIGI works with two syntaxes: the first one using a field-by-element
of Hook matrix and the second one a field-by-element of material characteristics. Although this is
obviously less strategic than for the SIGM operator (usually the elastic stiffness is computed only
one time when a model is set), it is possible to prevent RIGI from activating calls at the SECTION level in the case of the second syntax when the component MAHO is available.
As usual the driver routine (RIGI1) is modified in order to account for complex material characteristics. RIGI1 calls RIGI4 which either directly loads the Hook matrix or computes it (FRIGIE,
FRIGI2 and DOHTIF). The element stiffness matrix is then computed by mean of two new routines: TIFILO (computing the local element stiffness matrix) and TIFRIG (driving the transformation local axis/global axis).
Note that the product of the stiffness matrix computed with RIGI and any global displacement field
should be equal to the global resultant, computed by BSIG, of the element elastic stress, computed
before using SIGM and associated with the same displacement field.
6.4.4 Mass matrix calculation (MASS operator)
The operator MASS computes the contribution of each element to the global mass matrix. This
operator is generally used at most one time after having setting up a model.
For the TIMO element with fibres, the driver routine MASS1 is modified (complex material characteristics). The routine MASS3 then activates the SECTION level (call to FRIGIE and FRIGI2)
in order to determined the section mass matrix which is loaded by DOHTIF. Two routines TIFALO
(computing the local element mass matrix) and TIFMAS (driving the transformation local axis/global axis) have been added.

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

36 of 55

IMPLEMENTATION NOTES

6.4.5 Initialisation (ZERO operator)


The ZERO operator is used, between other things, for making the initialization of the field-by-element of internal variables at the beginning of any non linear calculation. In particular, this is the case
of the two components (CONS and VAIS) associated to the beam formulation. As already discussed, these components are field-by-elements over the section. Instead of performing the initialization at the two levels, beam and SECTION, which would require the knowledge of the
SECTION model defined as beam material characteristics, only the beam level is initialized
according to the following trick: the two component fields are created, are correctly typed, but are
associated to a void object. The further modelling routines, used for the non-linear TIMO element with fibres, are allowed to recognize a void object as being a null object and to perform, if
required, the initialization at the SECTION level (calling usual slave routines of the operator
ZERO).
The trick used for the initialization at the beam level is introduced in the routine ZEROP.
6.4.6 Non linear flow evaluation (ECOU operator)
The ECOU operator is used for performing all the incremental evaluation of the non-linear constitutive laws in CASTEM 2000 written under the form:

new

old

+ (

old

old

, )

new

(47)
old

where the stress


, at the end of the increment, is a function of the stress and of the internal
old
e
variables at the beginning of the increment, and of the increment of elastic stress .
The two levels of modelling are treated in a separate way.
The beam (POUTRE) level is implemented as usual [2]: introduction of a new model number
(NOMATE), modification of ECOUL1 for accounting to the special value type of the non linear
components of the material characteristics and of the internal variables (the default being real
number value), addition of a slot in ECOUL2 and call to the new routine which is devoted to the
fibre modelling: BIFLEX.
The subroutine BIFLEX, called for each integration point of each TIMO element with fibres, is
the interface between the two levels. In particular, the data coming from the beam level are rearold
ranged in such a way the model, the characteristics, the initial stress
and the initial internal
old
variable
are accessible as usual. The behaviour of each fibre being specified kinematically
(Eq. (4)), the increment of elastic generalized stress is inverted using the section stiffness. In consequence, the SECTION models are generally formulated as

new

old

+ (

old

old

, )

(48)

The SECTION level is then implemented using the subroutines FCOUL1, FCOUL2 and specialized routines for each kind of fibre. In fact FCOUL1 and FCOUL2 work for the section as
ECOUL1 and ECOUL2 work in general. FCOUL1 loops on the various sub-zones of the section
model and prepares the constant data for the zone, FCOUL2 loops on the elements and on the integration point(s) of each zone. The main difference between FCOUL2 and ECOUL2 comes from
the introduction of the coordinates of the integration points in order to retrieve, using (4) and the
generalized strain increment, the strain increment at the current point. In the same way, after having computed the new local section stress state using a specific constitutive model, the numerical
contribution to the generalized stress integrals (6.1)-(6.4) is derived using the standard Gauss formula. A computed GOTO in FCOUL2 allows to branch to each specific model depending on the
relative section model index INFIBR (ranging from 0). The transformation between the absolute
model index (set by NOMATE) and the relative section model is performed by the new subroutine

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

37 of 55

IMPLEMENTATION NOTES

TEMANF. The relative numbering of the non linear section models should be compliant with the
one introduced for the model name definition (MODPLA), the material characteristics (IDPLAS)
and the internal variables (IDVARI). As in ECOUL2, the slot 0 is reserved to the elastic model.
Three specific section models are available. The first one is purely elastic and the corresponding
subroutine is FIBELA. The second model concerns the concrete behaviour and the corresponding
subroutine is FIBETO which calls in turn NEWBET. NEWBET is a pure FORTRAN subroutine
which implements the 1-D concrete fibre model described in Section 3.1. The third model deals
with the steel behaviour. The associated subroutine is FIBSTE which calls in turn STEEL1.
STEEL1 is a pure FORTRAN subroutine which implements the 1-D steel fibre model described in
Section 3.2. FIBETO and FIBSTE are not only simple interface routines with the CASTEM 2000
data-structures and the dedicated routines which perform the modelling, they also complete elastically the behaviour with respect to the shear components xy and xz of the stress (see the discussion about generalized and simple fibre modelling in Section 2.5 and Section 2.6)
6.4.7 Tangent Stiffness (KTAN operator)
Due to the complex nature of the global section behaviour, the convergence of the algorithms introduced for performing the step-by-step non linear analysis requires a frequent evaluation of the tangent stiffness matrix. This evaluation is to be performed by the KTAN operator.
Due to the particular nature of the fibre models herein introduced (usual fibre modelling of chapter
2.6), the tangent constitutive matrix in the current configuration depends only on the tangent value
tang
tang
E
of the constitutive curve of all the fibres in this configuration. Substituting E by E
in (39),
tang
it is possible to derive the tangent section matrix K
which is then used for deriving the global tangent stiffness.
As in the implementation of the ECOU operator, the two levels of modelling are treated in a separate way.
The beam level requires the usual modifications of the subroutine KTANGA: convenient treatment of the special kind of component of the field-by-element of the specific type characteristics and internal variables for the TIMO element with fibres. In the same routine are
performed the nested loops on the elements and Gauss points of each TIMO with fibres subzone. The work inside the loops is arranged as follow: First is computed (FRIGTA) the set of the
12 section integrals which allows to derive (DOHTIF) the tangent section constitutive matrix
tang
K
(tangent Hook matrix). The computation then proceed as for the evaluation of the elastic
stiffness matrix (use of TIFRIG and TIFILO as in Section 6.4.3).
The SECTION level is implemented using the subroutine FRIGTA and FRIGT2. FRIGTA performs the loop on the various sub-zones of the section models and prepares the constant data for
the zone. FRIGT2 loops on the elements and on the integration points of each zone. Depending
tang
on the model, the value of the internal variable containing E
is retrieved and the numerical
contributions to the section integrals (39) are computed.

6.5 Operators working for the SECTION formulation


Letting aside the operators MODE, MATE and CARA which work directly at the section level
after having modified the model-related routines (see Section 6.2) and the characteristics-related
identification routines (see Section 6.3), a new operator named MOCU have been specifically
developed for dealing with the section modelling.
6.5.1 Biaxial bending of a section (MOCU operator)
The MOCU operator (see chapter 5.4.1) allows to derive the bi-axial moment/curvature curves of
any section described by a SECTION model and the associated characteristics. Normal forces can
A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

38 of 55

IMPLEMENTATION NOTES

also be considered. At the implementation level, the subroutine MOCUR first reads the input, computes the elastic constitutive matrix (FRIGIE and FRIGI2) and then equilibrates the section for
each point of the loading path using a simple secant method and the evaluation of the generalized
stress on the section (FCOUL1 and FCOUL2). A manual of the operator MOCU is also available.

6.6 How to implement a new fibre model


The method to follow in order to implement a new SECTION non linear model is more or less the
one indicated in [2], with some slight differences due to the two-level implementation. A short
guideline is here indicated.
Model: in the subroutine MODPLA add 1 to the number NMOD of possible models and add a new
position to the array MOMODL by setting the name of the model. The index position of the model is
the relative plastic model index INPLAS.
Model index: an absolute model index (INATU) is associated to the new model in the subroutine
NOMATE. The correspondence between this absolute number and the relative SECTION model
index INFIBR should be established in the subroutine TEMANF.
Material characteristics: The names of the required non linear material characteristics are specified in the subroutine IDMATR by adding a new slot (ELSEIF structure for the local index
IPLAC) corresponding to the new relative plastic model index. The names are stored in the array
MOTS.
Internal variables: The names of the required internal variables are specified in the subroutine
IDVARI by adding a new slot (ELSEIF structure for the local index MATPLA) corresponding to
the new relative plastic model index. The names are stored in the array LESOBL. A position for
tang
the component TANG (containing E ) should be presented. A position for the component
EPSO (containing the total elongation x of each fibre) is strongly advised for performing the
post processing.
Interface with FCOUL2 (ECOU): In FCOUL2 a new slot (ELSEIF structure for the relative
fibre index INFIBR) should be added. A call to a subroutine, say FIXXXX, is then required. This
routine is used as an interface between FCOUL2 and, say MOXXXX, which implements effectively
tang
the new model. It is also the role of the interface to collect the value of E
and possibly the one
of x .
Tangent modulus in FRIGT2 (KTAN): in FRIGT2 a new slot (ELSEIF structure for the local
index INFIBR) should be added. The local variable YOUNGT is set according to the value of the
internal variable associated to the component TANG. This value is found in the array XVAR at
the position decided by the developer.

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

39 of 55

IMPLEMENTATION NOTES

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

40 of 55

NUMERICAL APPLICATIONS

7. NUMERICAL APPLICATIONS
7.1 Introduction
In order to validate the fibre model, as well as the presented fibre constitutive laws, and to illustrate
its potential, a set of monotonic and cyclic numerical tests on R/C structural elements has been performed. According to the kind of loading, two main groups are considered: monotonic and cyclic.
Two tests are included in the first group:
The first one simulates the response of a rectangular section for seven different reinforcement
ratios. The concrete and steel stress-strain behaviour laws are defined according to the theoretical
curves assumed by Park and Paulay [16]. It must be noticed that the aim of this test is not to validate the constitutive laws but to check the fibre algorithm procedures in CASTEM 2000.
The second test, also performed on a rectangular R/C structural element, intends to show the
potential of the fibre model in the sense that, if a proper constitutive law is adopted, phenomena
like cracking can be well represented. The importance of the concrete tensile strength in the global response of R/C elements under flexural moments is also analysed.
In the second group, where almost all the cases could be included, two numerical tests are presented:
The first one, performed on a R/C cantilever T beam with B500 steel, tries to simulate the LNECs
experiment on specimen S1V3 as described by Pipa and Carvalho [17]. From the analysis of the
results, some considerations concerning the crack closing phenomenon are presented.
Classical global behaviour models usually fail when biaxial type loading is considered. In order to
validate the fibre model under this type of loading, the numerical simulation of one of the tests
referred by Bousias [18] is presented as second test. Some improvements are then suggested in the
Giuffr-Pinto model so that strength degradation at high strain level could be properly taken into
account.

7.2 Monotonic tests


7.2.1 Doubly reinforced rectangular cross sections
The first numerical test is performed on a R/C rectangular section. Different volumetric steel ratios
are considered at both sides of the beam. The loading, an increasing curvature path, is imposed to the
structure and moments are computed by the fibre algorithm.
The response curves (see Fig. 11) are to be compared with the theoretical results obtained by another
algorithm but supposing the same behaviour (stress-strain laws) for the materials.
The monotonic curve for the concrete considers the compression strength ( co = 27.6MPa ) , for
( co = 0.2%). The softening behaviour is described by ( Z = 100.0 ). No compression plateau nor
tensile stresses are considered.
For the steel, an elastic perfectly plastic law is adopted. This two parameter law is defined by a yield
stress ( sy = 276.0MPa ) and a Young modulus ( E = 207.5GPa ).
Fig. 11 shows the fibre discretization considered in the section as well as the moment-curvature
results for the different reinforcement ratios of Table 1.
A very good agreement with the theoretical results presented by Park and Paulay [16] was achieved.

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

41 of 55

NUMERICAL APPLICATIONS

Table 1: Steel reinforcement ratios (see Fig. 11)


Beam

0.0375

0.0375

0.0375

0.0250

0.0250

0.0125

0.0125

0.0250

0.0125

0.0125

0.0125

Moment/(b*d2) (MN/m2)
11.0
As

d/10
1
2

8.8
3

6.6

4
5
As

d/10

4.4

b
6
7

2.2

As
= -----bd

Curvature (1/m)
0.0
.00

.56

1.12

1.68

2.24

As
= ------bd

2.80

X1.E-2

Figure 11- Fibre results according to Park and Paulay tests.

7.2.2 Influence of the concrete tensile strength


The analysis of the influence of the concrete tensile strength in the global response of R/C structural
elements under flexural moments is presented. The structure under analysis, a square doubly symmetrical section, is illustrated in Fig. 12 together with the fibre discretization.
In this test, a compressive strength ( co = 30.0MPa ) , for ( co = 0.2 %) and a tensile strength
( t = 3.0MPa ) , were adopted for the concrete. For the steel, a Young modulus ( E o = 203.0GPa )
and a yield stress ( sy = 460.0MPa ) were considered.
The concrete cover, ( c = 0.19m ), was modeled by two external fibres, one on the top and the other
on the bottom of the section. For all the other fibres representing the concrete core, a confinement
factor ( = 1.331 ) was adopted. Since the aim of this test is the analysis of the response curve just
before and after the cracking process, all the other parameters defining both the steel and the concrete constitutive laws are irrelevant to the response and so are not presented.
Four numerical tests were performed and the associated results are presented in Fig. 12. In each one
a different cracking model was adopted. For the first test a sudden rupture was chosen, ( tm = t ). In
the second, third and fourth tests the ratio ( r = tm t ) takes values equal to 1.5 , 2.0 , and 22.7 ,
respectively. The dashed line illustrates the response of the section when zero tensile strength is considered. The fourth result follows the suggestion of Barzegar [10].

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

42 of 55

NUMERICAL APPLICATIONS

Moment (KNm)
20.0
4 -> r = 22.7

1 -> r = 1.0
2 -> r = 1.5

16.0

3 -> r = 2.0
3

12.0

316mm

2
.031m

1
8.0

.25m
.019m

4.0

.25m

Curvature (1/m)
.0
.00

1.00

2.00

3.00

4.00

5.00

6.00
X1.E-3

Figure 12- Influence of the tensile parameters in the section global response.

Thus, when a section is subjected to an increasing curvature loading path, the moment-curvature
response before cracking is represented by an almost straight line. During this first step the R/C
behaves as linear elastic and the three results are identical. When the first set of the fibres Gauss
points reaches the tensile strength, cracks start opening and the curve inflects tending to the response
curve of the zero tensile strength case.
These results illustrate the importance of adopting a tensile softening behaviour after cracking, since
a non ductile behaviour law, ( tm = t ), gives the discontinuous response represented by curve 1.
This curve reflects the fact that when cracks open the tensile strength goes immediately to zero,
changing instantaneously the internal equilibrium of the section. Part of the tensile force is then
transferred to other concrete and steel fibres and another part is compensated by a movement of the
neutral axis, allowing a decrease of the global compression forces.
This problem, affecting the convergence of the nonlinear algorithm, can, of course, be minimized by
considering a greater number of fibre elements, but this will also increase the computation costs.
Thus, post-cracking ductile behaviour has been adopted in the model, allowing not only a more
smooth response, as presented by curves 2, 3 and even 4 in Fig. 12, but also the possibility to account
for phenomena such as tension-stiffening in the final response.
It should be noticed that Barzegar [10] suggestion produces a ratio ( r = 22.7 ) which seems to be too
high.

7.3 Cyclic tests


7.3.1 Uniaxial Loading on a Cantilever T Beam
A series of tests on R/C cantilever T beams were performed at LNEC within Task D of the Cooperative Research Program on Seismic Response of R/C Structures (2nd Phase). The main aim of these
experimental tests was to provide information on the nonlinear behaviour of the beams of a four storey full scale building tested in ELSA using the PSeudo Dynamic (PSD) test method.
A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

43 of 55

NUMERICAL APPLICATIONS

412mm

46mm

0.10m
312mm

6mm // 7mm

0.20m

0.20m

0.40m

Figure 13- Section of the T beam.

The test performed on the specimen S1V3 is simulated by the fibre algorithm. All the concrete was
considered confined with a compressive strength ( co = 35.0 MPa ), at ( co = 0.2 %), and a slope
for the descending part of the envelope ( Z = 100.0 ) . A confinement factor ( = 1.33 ) and a Poisson ratio ( = 0.25 ) were also adopted. No compression plateau nor tensile strength was considered
in the concrete model.
For the steel we adopted the following characteristics: yield stress ( sy = 550.0MPa ), Young modulus ( E = 203.0GPa ), Poisson ratio ( = 0.25 ) , hardening and ultimate strain equal to 0.22 % and
9.0 %, respectively, ultimate stress ( su = 610.0MPa ) and hardening ratio ( b = 0.03 ).
Fig. 13 shows the cross section of the cantilever T beam and Fig. 14 the discretization adopted for
the numerical application. The structure was divided into 7 elements with increasing length from the
base to the free top, as follows: 0.15, 0.15, 0.15, 0.15, 0.30, 0.30 and 0.30m . The loading history, a
horizontal displacement, was imposed at the free top of the beam. No axial compression forces were
considered.

Force Y (KN)
60.0
40.0

Fibre
Experimental

20.0
.0
-20.0
-40.0
-60.0
Disp. UY (m)
-80.0
-.15

-.12

-.09

-.06

-.03

.00

.03

.06

.09

.12

.15

Figure 14- Response curves of the cantilever T beam.

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

44 of 55

NUMERICAL APPLICATIONS

31.0mm
0.25m

816mm

Y
19.0mm

15.0mm

8mm // 70mm
0.25m

Figure 15- Column section tested by Bousias.

The response curve from the fibre model, illustrated in Fig. 14 by the free top force-displacement
diagram, presents a rather good agreement with the experimental results. The discontinuities
observed in the dashed line curve result from the buckling of some of the bottom steel bars caused by
the opening of the stirrups.
The major differences between both diagrams appear during the crack closing. The pinching effect is
in the case of the fibre model much more pronounced. Although some modifications could be performed in the future in order to improve the behaviour of the concrete under cyclic loading, the
present model seems to fail on representing in a realistic way the crack closing phenomenon.
7.3.2 Biaxial Loading on a R/C Square Section
The structure tested by Bousias [18], a R/C column element with 1.49m length, framing into an
enlarged end block of reinforced concrete, intended to simulate part of a column between the foundation and the point of inflection. From the twelve specimens tested, the specimen S9 was chosen. The
cross section, a square doubly symmetrical one, is presented in Fig. 15.
For the loading history, the axial force and the horizontal displacements imposed at the free top of
the column during the experimental test were adopted. These values were read directly from a data
base and are presented in Fig. 16. The displacement path is what is usually called a shrinking path,
since, after an important displacement imposed to the structure, there is a progressive movement, in
both directions, towards the initial zero displacement position. The loading follows a four squares

Disp. UZ (mm)

Force X (KN)
110.

140.
2
4

90.

70.

1
.0

50.
7

30.

5
Disp. UX (mm)

Disp. UY (mm)

10.

-140.
-2.0

.0

2.0

4.0

-140.

.0

140.

Figure 16- Loading path: axial force and horizontal displacements.

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

45 of 55

NUMERICAL APPLICATIONS

displacement path centered at the origin with half-side lengths equal to 100, 80, 60 and 40mm in
each direction.
The characteristics adopted for the materials are the mean values found from the laboratory tests performed on concrete and steel specimens.
Around all the section a concrete cover of ( c = 19mm ) thickness with a compressive strength
( co = 30.0MPa ) , at ( co = 0.2 %), and a slope for the descending part of the envelope
( Z = 100.0 ), was considered. For the concrete core a confinement factor ( = 1.33 ) was introduced. For the tensile parameters, a strength ( t = 3.0MPa ) and a ratio ( r = 3.0 ) were used. The
model was completed assuming a Poisson ratio ( = 0.25 ) and no compression plateau.
For the steel we adopted the following average characteristics: yield stress ( sy = 440.0MPa ) ,
Young modulus ( E = 203.0GPa ), Poisson ratio ( = 0.25 ), hardening and ultimate strain equal to
0.8 % and 13.0 %, respectively, and ultimate stress ( su = 760.0MPa ) . For the hardening ratio
( b = 0.13 ) was adopted.
In terms of finite element discretization, the structure was divided into 7 elements with increasing
length from the base to the free top, as follows: 0.05, 0.10, 0.15, 0.20, 0.25, 0.35 and 0.39m. The discretization adopted for the section is represented in Fig. 15.
The response curves from the fibre model, showed in Fig. 17 and Fig. 18, present a good agreement
with the experimental results. The reduction of strength along the cycles, not evidenced in the fibre
results, is due to the degradation of the concrete and steel properties. After the first square path
almost all the concrete is destroyed and the response curve is mainly controlled by the steel.
The high strength degradation verified in the section, which results from the high strain levels
imposed to the column, is not well represented by the fibre model; Since no strength degradation is
considered in the steel model, the fibre response curves from the second and third square loading
cycles tend to the same asymptote. In what concerns the present test , this characteristic of the Giuffr-Pinto model is responsible for the most important differences found between the experimental
and the numerical results.

Force Y (KN)
60.0
Fibre
Experimental

40.0
20.0
.0
-20.0
-40.0

Disp. UY (m)
-60.0
-.10

-.08

-.06

-.04

-.02

.0

.02

.04

.06

.08

.10

Figure 17- Free top force-displacement response curves on Y direction.

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

46 of 55

NUMERICAL APPLICATIONS

Force X (KN)
80.0
Fibre
Experimental

60.0
40.0
20.0
.0
-20.0
-40.0

Disp. UX (m)
-60.0
-.10

-.08

-.06

-.04

-.02

.0

.02

.04

.06

.08

.10

Figure 18- Free top force-displacement response curves on Y direction.

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

47 of 55

NUMERICAL APPLICATIONS

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

48 of 55

SUMMARY AND CONCLUSION

8. SUMMARY AND CONCLUSION


A fibre Timoshenko beam element has been implemented in CASTEM 2000. The Timoshenko beam
theory has been considered in order to allow the development of the so-called generalized fibre
model in the sense that very complex interactions between normal force, bending moment and shear
can be potentially taken into account. Although a suitable constitutive model is still to be derived, the
most general formulation and the corresponding numerical data structure have been created.
In the mean time, uniaxial material constitutive relationships (usual fibre modelling) for concrete and
steel have been implemented.
Concrete behaviour is represented by a parabolic curve up to the peak stress followed by a straight
line in the softening zone. Confinement is taken into account by a modification of the plane concrete curve and including an additional plateau zone. Cyclic behaviour accounts for stiffness degradation and crack closing phenomena. Tensile resistance has been also considered.
Steel behaviour is represented by a modified Menegoto-Pinto model with a three stage monotonic
curve (linear, plateau and hardening). The Bauschinger effect is taken into account and buckling
can be simulated.
The set of numerical tests included in the report shows the potential of this uniaxial modelling for the
study of complex structural elements and loading paths. However some aspects need further
improvements; in particular, the modelling of the concrete crack closing.
Moreover, a suitable model able to represent the non-linear behaviour of structural elements dominated by shear (e.g. short bridge piers) should be developed. A series of tests programmed for the
ELSA laboratory on bridges (bridge piers) will constitute the basis for this development.
The implementation of the beam fibre-modelling in CASTEM 2000 has been performed following
the standard procedures for model implementation. However, in order to be general enough and to
benefit at most from the existing capabilities of the code, a two-level approach has been adopted.
The structural beam level allows to discretize complex boundary value problems, the Gauss integration points of each beam element being described at a lowest level: the section level.
The section level allows not only to specify the geometry and the material characteristics of any
section, but also to easily analyse the distribution of stress, strain and internal variables resulting
from the straining of this section.
As a result of this two-level approach, it is very easy to implement other fibre constitutive laws.

ACKNOWLEDGEMENTS
The authors would like to thank Alain Millard and Laurent Le Ber for their contribution (and hopefully a touch of rational elegance) to the last step of this work: the effective implementation of the
fibre approach in the official version of CASTEM 2000.

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

49 of 55

REFERENCES

REFERENCES
[1] CASTEM 2000, Guide dutilisation, CEA, France, 1990.
[2] Pegon, P., Model Implementation in CASTEM 2000: Some General Considerations and a Simple Practical Realization, JRC-Technical Note No. I.93.06, 1993, JRC-Special Publication No
I.94.05, 1994 and CEA-Report DMT/94-195, 1994.
[3] Oller, S., Barbat, A.H., Oate, E., Hanganu, A., A Damage model for the seismic analysis of
building structures, Proc. 10th World Conf. on Earthquake Engng., Ed. A.A. Balkema, Rotterdam, Vol. 5, pp 2593-2598, 1992.
[4] Combescure, D., Pegon, P., A Fiber model accounting for transverse shear in CASTEM 2000,
JRC-Special Publication, to appear in 1994.
[5] Hognestad, E., A Study of Combined Bending and Axial Load in Reinforced Concrete, Bulletin Series 339, Univ. of Illinois Exp. Sta., November 1951.
[6] Eurocode No 8 - Structures in Seismic Regions - Design, Part 1, General and Building, May
1988.
[7] Tassios, T. and Lefas, N., Ductility of confined reinforced concrete elements, Report NTUA (in
Greek), Athens, Greece, 1983.
[8] Park, R., Priestley, M., Ductility of Square-Confined Concrete Columns, Journal of the Structural Division, ASCE, Vol. 108, No ST4, pp. 929-950, April 1982.
[9] Mercer, C. and Martin, J., A Beam Element for Cyclically Loaded Reinforced Concrete Structures, Thecnical Report No. 98, FRD/UCT Centre for Research in Computational and Applied
Mechanics, University of Cape Town, November 1987.
[10] Barzegar-Jamshidi, F. and Schnobrich, W., Nonlinear Finite Element Analysis of Reinforced
Concrete Under Short Term Monotonic Loading, Civil Engrg. Studies, SRS 530, Univ. of Illinois at Urbana-Champaign, Urbana, Illinois, November 1986.
[11] Monti, G. and Nuti, C., Nonlinear Cyclic Behaviour of Reinforcing Bars Including Buckling,
Journal of Structural Engineering, Vol. 118, No 12, December 1992.
[12] Menegotto, M. and Pinto, P., Method of Analysis for Cyclically Loaded Reinforced Concrete
Plane Frames Including Changes in Geometry and Nonelastic Behaviour of Elements under
Combined Normal Force and Bending, IABSE Symposium on Resistance and Ultimate
Deformability of Structures Acted On by Well-Defined Repeated Loads, Final Report, Lisbon,
1973.
0

[13] Donea, J. and Lamain, L.G., A modified representation of transverse shear in C quadrilateral
plate elements, Int. J. Num. Meth. Engng., Vol. 63, pp 183-207, 1987.
[14] Pegon, P., A Timoshenko simple beam element in CASTEM 2000, JRC-Technical Note No.
I.93.05, 1993, JRC-Special Publication No I.94.04, 1994 and CEA-Report DMT/94-200, 1994.
[15] ESOPE/GEMAT, MANUEL dutilisation, CEA, France, 1989.
[16] Park, R. and Paulay, T., Reinforced Concrete Structures, John Wiley & Sons, 1975.
[17] Pipa, M., Carvalho, E., Cyclic Tests on R/C Cantilever T Beams with B500 Steel - Cooperative
Research Program on Seismic Response of R/C Structures, LNEC, Lisboa, June 1993.
[18] Bousias, S., Load Path Effects in Reinforced Concrete Column - Biaxial Bending with Axial
Force, Technical Note No I.92.63, Commission of the European Communities, Joint Research
Centre, Ispra Site, Italy, June 1992.

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

50 of 55

ANNEX 1- A SAMPLE GIBIANE INPUT FOR FIBRE APPLICATION

ANNEX 1- A SAMPLE GIBIANE INPUT FOR FIBRE


APPLICATION
*-------------------------------------------*
*
EXAMPLE OF A "T" R.C. SECTION WITH STEEL ON
*
THE TOP AND BOTTOM OF THE BEAM
*
(L.N.E.C. TESTS)
*
*-----------------------------------------*
Quadrangular Elements
*-----------------------------------------opti dime 2 elem qua4;
*-----------------------------------------*
Triangular Elements (commented out)
*-----------------------------------------* opti dime 2 elem tri3;
*-----------------------------------------*
DEFINITION OF THE STEEL GEOMETRY (SECTION)
*-----------------------------------------ps1 = -1.41372d-2 -17.6d-2;
ps2 = -1.41372d-2 -16.4d-2;
ps3 = 2.82743d-2
0.0
;
cars1 = ps1 d 1 ps2 tran 1 ps3;
*
ps4 = -1.88496d-2 7.6d-2;
ps5 = -1.88496d-2 6.4d-2;
ps6 = 3.76991d-2 0.0
;
cars2 = ps4 d 1 ps5 tran 1 ps6;
*
ps7 = -30.9425d-2 7.0d-2;
ps8 = -30.9425d-2 6.4d-2;
ps9 = 1.88496d-2 0.0
;
cars3 = ps7 d 1 ps8 tran 1 ps9;
*
cars4 = cars3 plus (60.0d-2
0.0);
*
mesf = coul (cars1 et cars2 et cars3 et cars4) roug;
*-----------------------------------------*
DEFINITION OF THE UNCONFINED CONCRETE GEOMETRY (SECTION)
*-----------------------------------------pc1 = -50.0d-2 10.0d-2;
pc2 = -50.0d-2
0.0d-2;
pc3 = 40.0d-2
0.0d-2;
carc1 = pc1 d 4 pc2 tran 1 pc3;
carc2 = carc1 plus (60.0d-2 0.0);
*
pc7 = -10.0d-2 10.0d-2;
pc8 = 10.0d-2 10.0d-2;
pc9 =
0.0d-2 -2.1d-2;
carc3 = pc7 d 1 pc8 tran 1 pc9;
*
carc4 = carc3 plus (0.0 -27.9d-2);
*
mesu = coul (carc1 et carc2 et carc3 et carc4) turq;
*-----------------------------------------*
DEFINITION OF THE CONFINED CONCRETE GEOMETRY (SECTION)
*-----------------------------------------pc10 = -10.0d-2
7.9d-2;
pc11 = -10.0d-2 -17.9d-2;

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

51 of 55

ANNEX 1- A SAMPLE GIBIANE INPUT FOR FIBRE APPLICATION

pc12 = 20.0d-2
0.0d-2;
carc5 = pc10 d 10 pc11 tran 1 pc12;
*
mesc = coul carc5 bleu;
*
elim (mesf et mesu et mesc) 0.001;
*-----------------------------------------*
TOTAL MESH (SECTION)
*-----------------------------------------*
titre LNEC TESTS - "T" SECTION;
*
trac (mesf et mesu et mesc);
*
*-----------------------------------------*
DEFINITION OF THE TESTED BEAM (BEAM)
*-----------------------------------------opti dime 3;
*
pp0 =
0.0d-2 0. 0.;
pp1 = 15.0d-2 0. 0.;
pp2 = 40.0d-2 0. 0.;
pp3 = 75.0d-2 0. 0.;
pp4 = 110.0d-2 0. 0.;
pp5 = 150.0d-2 0. 0.;
*
llb = pp0 d 1 pp1 d 1 pp2 d 1 pp3 d 1 pp4;
llp = pp4 d 1 pp5;
*-----------------------------------------*
DEFINITION OF THE AUXILIARY BEAM (FOR THE MAHO COMPONENT)
*-----------------------------------------pa0= 0 0 0; pa1=1 0 0; lla=pa0 d 1 pa1;
*-----------------------------------------*
CARACTERIZATION OF THE STEEL AND CONCRETE MODELS (SECTION LEVEL)
*-----------------------------------------*-----------------------------------------*
Quadrangular Elements QUAS
*-----------------------------------------modf=mode mesf mecanique elastique plastique acier_uni quas;
modu=mode mesu mecanique elastique plastique beton_uni quas;
modc=mode mesc mecanique elastique plastique beton_uni quas;
*-----------------------------------------*
Triangular Elements TRIS (commented out)
*-----------------------------------------* modf=mode mesf mecanique elastique plastique acier_uni tris;
* modu=mode mesu mecanique elastique plastique beton_uni tris;
* modc=mode mesc mecanique elastique plastique beton_uni tris;
*-----------------------------------------*
Steel
*-----------------------------------------matf=mate modf YOUN 2.06e5 NU
0.30
STSY 548.0 EPSU .098
STSU 625.0 EPSH 0.023
FALD 4.900 A6FA 620.0 CFAC 0.5 AFAC 0.008
ROFA 20.0
BFAC 0.0116 A1FA 18.5 A2FA 0.15;
carf=carb modf ALPY 1.0
ALPZ 1.0;
*-----------------------------------------*
Unconfined concrete
*-----------------------------------------matu=mate modu YOUN 0.35e5 NU
.25
STFC 35.0
EZER .002
STFT 0.0
ALF1 .22429 OME1 .12953

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

52 of 55

ANNEX 1- A SAMPLE GIBIANE INPUT FOR FIBRE APPLICATION

ZETA 100.0
ST85 .00
TRAF 2.0;
caru=carb modu ALPY .66
ALPZ .66;
*-----------------------------------------*
Confined concrete
*-----------------------------------------*
Initial concrete Young modulus = 2 * STIFC / ( BETA * EZERO )
*-----------------------------------------matc=mate modc YOUN 0.33e5 NU
.25
STFC 35.0
EZER .002
STFT 0.0
ALF1 .22429 OME1 .12953
ZETA 0.0
ST85 10.0
TRAF 2.0;
carc=carb modc ALPY .66
ALPZ .66;
*-----------------------------------------* Complete model and material at the SECTION level
*-----------------------------------------modq=modf et modu et modc;
macq=matf et matu et matc et carf et caru et carc;
*-----------------------------------------* Determination of the hook matrix at the BEAM level
*-----------------------------------------moa=mode lla mecanique elastique section timo;
maa=mate moa MODS modq MATS macq;
hoa=hook moa maa;
lrig=extr hoa MAHO 1 1 1;
*-----------------------------------------*
CHARACTERIZATION OF THE BEAM
*-----------------------------------------modfi=mode llb mecanique elastique section plastique section timo;
matfi=mate modfi MAHO lrig MODS modq MATS macq;
*
modpo=mode llp mecanique elastique section timo;
modpo=mode llp mecanique elastique section plastique section timo;
matpo=mate modpo MAHO lrig MODS modq MATS macq;
*
modto=modfi et modpo;
matto=matfi et matpo;
*-----------------------------------------*
STIFFNESS MATRIX
*-----------------------------------------rigi1 = rigi modto matto;
*
rigb1 = bloq DEPL pp0;
rigb2 = bloq ROTA pp0;
rigdz = bloq UZ pp5;
*
rigt = rigi1 et rigb1 et rigb2 et rigdz;
*-----------------------------------------*
LOADING (CYCLIC DISPLACEMENT ON THE TOP OF THE CANTILIVER BEAM)
*-----------------------------------------*
Loading program
*-----------------------------------------dhoz = prog -.0002 pas -.002 -10.0e-3
pas -.002 -30.0e-3
pas .002
30.0e-3
pas -.002 -40.0e-3;
npdy = dime dhoz;
xxxx = prog 1. pas 1. npas (npdy-1);
evolz = evol manu xxxx dhoz;
*-----------------------------------------*
Horizontal displacement on the top of the beam
*------------------------------------------

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

53 of 55

ANNEX 1- A SAMPLE GIBIANE INPUT FOR FIBRE APPLICATION

deslo = 1.0;
fdepz = depi rigdz deslo;
*
cargz = char fdepz evolz;
chart = cargz;
*-----------------------------------------*
RESOLUTION THE NON-LINEAR PROBLEM (USE OF NONLIN)
*-----------------------------------------tab1 = table;
tab1.PLASTIQUE = vrai;
tab1.MOVA = RIEN;
tab1.MAXITERATION = 20;
tab1.KTA = vrai;
*
nonlin rigt matto chart xxxx modto tab1;
*-----------------------------------------*
ANALYSIS OF THE RESULTS
*-----------------------------------------*
Curve (shear top force)/(horizontal top displacement)
*-----------------------------------------recont=tab1.RESUCONT;
redepl=tab1.RESUDEPL;
indj=index redepl;
ftop=prog 0;dtop=prog 0;
j=0; repe loopt (dime indj); j=j+1;
timj=indj.j;
ftop=inse ftop (extr (bsig modto recont.timj matto) pp5 FZ) (j+1);
dtop=inse dtop (extr redepl.timj pp5 UZ) (j+1);
fin loopt;
*
tt = table;
tt.1 = MARQ CARR;
tt.2 = MARQ CROI;
fdtop=evol manu DISPL. (m) dtop FORCE (MN) ftop;
*-----------------------------------------*
Reference curve (for checking the results)
*-----------------------------------------dtopref=prog
.00000000000000E+00 -2.00000000000000E-04 -2.16000000000000E-03
-4.12000000000000E-03 -6.08000000000000E-03 -8.04000000000001E-03
-1.00000000000000E-02 -1.20000000000000E-02 -1.40000000000000E-02
-1.60000000000000E-02 -1.80000000000000E-02 -2.00000000000000E-02
-2.20000000000000E-02 -2.40000000000000E-02 -2.60000000000000E-02
-2.80000000000000E-02 -3.00000000000000E-02 -2.80000000000000E-02
-2.60000000000000E-02 -2.40000000000000E-02 -2.20000000000000E-02;
dtopref=dtopref et (prog
-2.00000000000000E-02 -1.80000000000000E-02 -1.60000000000000E-02
-1.40000000000000E-02 -1.20000000000000E-02 -1.00000000000000E-02
-8.00000000000000E-03 -6.00000000000000E-03 -4.00000000000000E-03
-2.00000000000000E-03 6.07153216591883E-18 2.00000000000001E-03
4.00000000000001E-03 6.00000000000001E-03 8.00000000000001E-03
1.00000000000000E-02 1.20000000000000E-02 1.40000000000000E-02
1.60000000000000E-02 1.80000000000000E-02 2.00000000000000E-02);
dtopref=dtopref et (prog
2.20000000000000E-02 2.40000000000000E-02 2.60000000000000E-02
2.80000000000000E-02 3.00000000000000E-02 2.80000000000000E-02
2.60000000000000E-02 2.40000000000000E-02 2.20000000000000E-02
2.00000000000000E-02 1.80000000000000E-02 1.60000000000000E-02
1.40000000000000E-02 1.20000000000000E-02 1.00000000000000E-02
8.00000000000000E-03 6.00000000000000E-03 4.00000000000000E-03

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

54 of 55

ANNEX 1- A SAMPLE GIBIANE INPUT FOR FIBRE APPLICATION

2.00000000000000E-03 -6.07153216591883E-18 -2.00000000000001E-03);


dtopref=dtopref et (prog
-4.00000000000001E-03 -6.00000000000001E-03 -8.00000000000001E-03
-1.00000000000000E-02 -1.20000000000000E-02 -1.40000000000000E-02
-1.60000000000000E-02 -1.80000000000000E-02 -2.00000000000000E-02
-2.20000000000000E-02 -2.40000000000000E-02 -2.60000000000000E-02
-2.80000000000000E-02 -3.00000000000000E-02 -3.20000000000000E-02
-3.40000000000000E-02 -3.60000000000000E-02 -3.80000000000000E-02);
dtopref=dtopref et (prog -4.00000000000000E-02);
ftopref=prog
.00000000000000E+00 -1.14481311676523E-03 -1.22901793761370E-02
-2.32933890701765E-02 -3.41409060368847E-02 -4.48163634737750E-02
-5.52993760011572E-02 -6.25736672822648E-02 -6.34682582020739E-02
-6.37997368461639E-02 -6.39516852406453E-02 -6.40648825527967E-02
-6.41051047602156E-02 -6.40881950555459E-02 -6.40683885171247E-02
-6.40457149947723E-02 -6.40198575442582E-02 -5.31916838295068E-02
-4.24307708413968E-02 -3.18940241565422E-02 -2.15933027921714E-02;
ftopref=ftopref et (prog
-1.14642689013442E-02 -1.95490782397935E-03 4.34527534900837E-03
1.05039345010342E-02 1.64943367735540E-02 2.22738649284235E-02
2.77888455743349E-02 3.11954965063240E-02 3.11996112204543E-02
3.12036823442664E-02 3.12076823925978E-02 3.12115835624900E-02
3.12155450724221E-02 3.12193553958296E-02 3.12232037374745E-02
3.12269455176290E-02 3.13243163279959E-02 3.17866134984408E-02
3.22108792714725E-02 3.26000050114313E-02 3.29579270685599E-02);
ftopref=ftopref et (prog
3.32842841187497E-02 3.35810040717798E-02 3.38496367881653E-02
3.40918486154093E-02 3.43092715455886E-02 2.76024463612365E-02
2.10585068426046E-02 1.46524661635717E-02 8.49910700480509E-03
2.69258722007945E-03 -2.62926487360908E-03 -8.01317688346363E-03
-1.25038980801999E-02 -1.60322003928289E-02 -1.87262001823851E-02
-2.07658260777817E-02 -2.23232577557867E-02 -2.35408842921363E-02
-2.45206665833014E-02 -2.53333845014858E-02 -2.60273227929928E-02);
ftopref=ftopref et (prog
-2.66354492379779E-02 -2.71804925526341E-02 -2.76783464385084E-02
-2.81403048598080E-02 -2.85745259272785E-02 -2.89869974898484E-02
-2.93821816783630E-02 -2.97634513185043E-02 -3.01333896091121E-02
-3.42706858331460E-02 -4.38448863358707E-02 -5.31853566475303E-02
-6.18463597099409E-02 -6.72796758934449E-02 -6.94003595649139E-02
-7.03019196590381E-02 -7.08358037573039E-02 -7.12371916595880E-02);
ftopref=ftopref et (prog -7.16045522313159E-02);
fdtopr=evol rouge manu DISPL. (m) dtopref FORCE (MN) ftopref;
titre Computed curve (white) and reference curve (red);
dess (fdtop et fdtopr) tt;

A FiBRE/TIMOSHENKO BEAM ELEMENT IN CASTEM 2000

55 of 55

You might also like