Safety formats for non-linear finite
element analyses of reinforced concrete beams loaded to shear failure
Anton Ekesi¨o¨o Andreas Ekhamre
Master Thesis in Concrete Structures, June 2018 TRITA-ABE-MTB-18333
Department of Civil and Architectural Engineering Division of Concrete Structures
Abstract
There exists several different methods that can be used to implement a level of safety when performing non-linear finite element analysis of a structure. These methods are called safety formats and they estimate safety by different means and formulas which are partly discussed further in this thesis.
The aim of this master thesis is to evaluate a model uncertainty factor for one safety format method called the estimation of coefficient of variation method (ECOV) since it is suggested to be included in the next version of Eurocode. The ECOV method will also be compared with the most common and widely used safety format which is the partial factor method (PF).
The first part of this thesis presents the different safety formats more thoroughly followed by a theoretical part. The theory part aims to provide a deeper knowledge for the finite element method and non-linear finite element analysis together with some beam theory that explains shear mechanism in different beam types.
The study was conducted on six beams in total, three deep beams and three slender beams. The deep beams were previously tested in the 1970s and the slender beams were previously tested in the 1990s, both test series were performed in a laboratory. All beams failed due to shear in the experimental tests. A detailed description of the beams are presented in the thesis. The simulations of the beams were all performed in the FEM- programme ATENA 2D to obtain high resemblance to the experimental test.
In the results from the simulations it could be observed that the ECOV method generally got a higher capacity than the PF method. For the slender beams both methods received rather high design capacities with a mean of about 82% of the experimental capacity.
For the deep beams both method reached low design capacities with a mean of around 46% of the experimental capacity. The results regarding the model uncertainty factor showed that the mean value for slender beams should be around 1.06 and for deep beams it should be around 1.25.
Keywords: nonlinear finite element analysis, estimation of coefficient of variation, par- tial factor, model uncertainty, deep beams, shear
Sammanfattning
Det finns flera olika metoder som kan anv¨andas f¨or att implementera en s¨akerhetsaspekt vid utf¨orande av en icke-linj¨ar finit element analys p˚a en konstruktion. Dessa metoder kallas s¨akerhetsformat och de uppskattar s¨akerheten med olika medel och formler som delvis diskuteras i denna avhandling.
Syftet med detta examensarbete ¨ar att utv¨ardera en modellos¨akerhetsfaktor f¨or ett s¨akerhetsformat som kallas estimerad variationskoefficient-metoden (ECOV), eftersom det ligger som
f¨orslag att den ska ing˚a i n¨asta version av Eurocode. ECOV-metoden kommer ocks˚a att j¨amf¨oras med det vanligaste och mest anv¨anda s¨akerhetsformatet vilket ¨ar partialfak- tormetoden (PF).
Den f¨orsta delen av denna avhandling presenterar de olika s¨akerhetsformaten mer nog- grant f¨oljt av en teoretisk del. Teoridelen syftar till att ge en djupare kunskap om den finita elementmetoden och icke-linj¨ar finit elementsanalys tillsammans med en del balk- teori som f¨orklarar skjuvmekanismer i olika balktyper.
Studien genomf¨ordes p˚a totalt sex balkar, tre h¨oga balkar och tre slanka balkar. De h¨oga balkarna vart tidigare pr¨ovade under 1970-talet och de slanka balkarna vart tidigare pr¨ovade p˚a 1990-talet, b˚ada provningarna utf¨ordes i ett laboratorium. Alla balkar gick s¨onder p˚a grund av skjuvning under de tidigare provningarna. En detaljerad beskrivn- ing av balkarna presenteras i avhandlingen. Simuleringarna av balkarna utf¨ordes alla i FEM-programmet ATENA 2D f¨or att uppn˚a h¨og likhet med de tidigare provningarna.
I resultaten fr˚an simuleringarna kan det konstateras att ECOV-metoden generellt sett har en h¨ogre kapacitet ¨an PF-metoden. F¨or de smala balkarna uppn˚adde b˚ada metoderna generellt ganska h¨oga dimensionerande lastkapaciteter, vilka i medeltal var omkring 82%
av den experimentella kapaciteten. F¨or de h¨oga balkarna uppn˚adde b˚ada metoderna generellt l˚ag dimensionerande lastkapacitet, vilka i medeltal var omkring 46% av den experimentella kapaciteten. Resultaten avseende modellos¨akerhetsfaktorn visade p˚a ett medelv¨arde omkring 1.06 f¨or slanka balkar och omkring 1.25 f¨or h¨oga balkar.
Nyckelord: Icke-linj¨ar finit elementanalys, estimering av variationskoefficient, partial- faktor, modellos¨akerhet, h¨oga balkar, skjuvning
Preface
This master thesis was conducted through a collaboration between the Division of Con- crete Structures from the Department of Civil and Architectural Engineering at the Royal Institute of Technology (KTH) and Tyr´ens AB. The thesis was conducted during the spring semester of 2018 at the headquarters of Tyr´ens AB in Stockholm.
We would like to send our deepest gratitude to our supervisor and examiner Adjunct Professor Mikael Hallgren for his guidance, support and for always taking the time to assist throughout the duration of this thesis. We would also send a special thanks to Dr.
Richard Malm for taking his time to provide us with necessary information, which this thesis would have been hard to conduct without. We would also like to thank Professor Anders Ansell for his effort to find us a suitable topic as well as for providing us with a contact at Tyr´ens and to Professor Anders Ericsson for his helpful and valuable tips.
Finally, we thank Tyr´ens AB and especially department “Konstruktion 5” for providing us with a space at their department, with a warm, including and positive atmosphere and for providing us with all the tools necessary to conduct this master thesis.
Stockholm, June 2018
Anton Ekesi¨o¨o Andreas Ekhamre
Contents
Abstract iii
Sammanfattning v
Preface vii
Symbols xi
Abbreviations xv
1 Introduction 1
1.1 Background . . . . 1
1.2 Aim . . . . 2
1.3 Limitations . . . . 2
1.4 Outline . . . . 2
2 Safety formats 5 2.1 Partial safety factor method (PF) . . . . 5
2.2 Global resistance method with estimation of coefficient of variation of re- sistance (ECOV) . . . . 5
2.3 Other methods . . . . 6
2.3.1 Global resistance factor method . . . . 6
2.3.2 Full probabilistic method . . . . 7
2.4 Model uncertainty factor . . . . 7
3 Finite element method 9 3.1 Finite element analysis . . . . 9
3.2 Non-linear finite element analysis . . . . 11
3.2.1 Iterative procedure . . . . 12
3.2.2 Non-linear material behaviour . . . . 13
3.2.2.1 Non-linear behaviour of concrete . . . . 13
3.2.2.1.1 Smeared crack models . . . . 14
3.2.2.2 Non-linear behaviour of reinforcement . . . . 15
4 Shear beam theory 17 4.1 Deep beams . . . . 17
4.2 Slender beams . . . . 18
5 Description of the beams 21 5.1 Beams by Holmgren and Nylander . . . . 21
5.2.1 Beam B90SB17-2-45 . . . . 27
5.2.2 Beam B91SC4-2-69 . . . . 27
5.2.3 Beam B91SD3-4-66 . . . . 28
6 NLFEA of the beams with Atena 2D 29 6.1 Indata . . . . 29
6.2 Modelling limitations, assumptions and simplifications . . . . 31
6.3 Material models . . . . 31
6.3.1 Concrete . . . . 31
6.3.2 Reinforcement . . . . 32
6.4 Boundary conditions and loading . . . . 33
6.5 Mesh and elements . . . . 34
6.6 Iterative solver . . . . 35
7 Results of the NLFEA 37 7.1 Beams by Holmgren and Nylander . . . . 37
7.1.1 Deep beam 1 . . . . 37
7.1.2 Deep beam 7 . . . . 40
7.1.3 Deep beam 9 . . . . 42
7.1.4 Comparison of safety formats . . . . 44
7.2 Beams by Hallgren . . . . 45
7.2.1 Beam B90SB17-2-45 . . . . 45
7.2.2 Beam B91SC4-2-69 . . . . 47
7.2.3 Beam B91SD3-4-66 . . . . 49
7.2.4 Comparison of safety formats . . . . 50
7.3 Model uncertainty factor . . . . 52
8 Discussion 53 9 Conclusion 57 9.1 Future research . . . . 57
Bibliography 59
Appendix A Calculations for the ECOV method 61
Appendix B Analytical calculations according to Eurocode 2 65
Appendix C Illustrations of beams from ATENA 2D 81
Symbols
Uppercase letters
Ac Area of the concrete cross section [m2]
As Reinforcement area [m2]
As,1 The area of one reinforcement bar [m2]
Asl Area of the tensile reinforcement [m2]
Asw Area of the shear reinforcement [m2]
CRd,c Coefficient [–]
C Force in strut [kN]
C1 Force in strut [kN]
C2 Force in strut [kN]
D Global displacement vector [mm]
Ec Modulus of elasticity for concrete [GPa]
Es Modulus of elasticity for steel [GPa]
Esh Modulus of elasticity for hardening of steel [GPa]
Fc Compression force [kN]
Fs Tension force [kN]
Gf Fracture energy [Nm/m2]
I Moment of inertia [m4]
I Internal forces [N]
Ia Internal force [N]
[K] Global stiffness matrix [N/m]
K0 Tangential stiffness [N/m]
L Length [m]
Mi Bending moment in element node [Nm]
NEd Axial force in the cross section due to loading or prestressing [N]
Ni Axial load in element node [N]
P External load [N]
Pd Ultimate design load [kN]
Pexp Ultimate load obtained in experiment [kN]
Psim Ultimate load obtained from simulations [kN]
Pu,k Ultimate load from analysis with characteristic values [kN]
Pu,m Ultimate load from analysis with mean values [kN]
R Global load vector [N]
R1 Reaction force at support 1 [kN]
R2 Reaction force at support 2 [kN]
RA Reaction force at support A [kN]
Ra Force residual [N]
Rd Design resistance [kN]
Rk Structural resistance for characteristic values [kN]
T1 Force in tie [kN]
VR Coefficient of variation of the resistance [–]
VRd,c Design value for the shear force in members without shear rein- forcement
[kN]
VRd,s Design value for the shear force in members with shear reinforce- ment
[kN]
VRd,maxMaximum design value for the shear force [kN]
Vg Coefficient of variation for geometry [–]
Vf Coefficient of variation for material uncertainties [–]
Vi Shear force in element node [N]
Vθ Coefficient of variation for model uncertainties [–]
Lowercase letters
a Shear span [mm]
a2 Width of the compression strut [m]
as Distance between node point and edge of the beam [mm]
bf Width of the flange [m]
bw Smallest width of the cross section of the structure [mm]
ca Displacement correction [m]
cb Displacement correction [m]
d Effective height [mm]
d Local displacement vector [mm]
fc Compressive strength [MPa]
fcc Compressive cylinder strength [MPa]
fcd Design compressive strength [MPa]
fck Characteristic compressive strength [MPa]
fct Tensile strength [MPa]
fcube The compressive cubic strength [MPa]
ft Tensile strength [MPa]
fy Yield strength [MPa]
fy,8mm Yield strength for 8 mm reinforcement bars [MPa]
fy,12mm Yield strength for 12 mm reinforcement bars [MPa]
fy,16mm Yield strength for 16 mm reinforcement bars [MPa]
fyk Characteristic yield strength [MPa]
fyd Design yield strength [MPa]
fywd Design yield strength of the shear reinforcement [MPa]
fctd Design tensile strength for concrete [MPa]
h Height of the beam [m]
k Coefficient [–]
k1 Coefficient [–]
k2 Coefficient [–]
[k] Local stiffness matrix [N/mm]
la Length of support [m]
m1 Material coordinate [–]
m2 Material coordinate [–]
n Number of reinforcement bars [–]
r Local load vector [N]
s Spacing between stirrups [mm]
CONTENTS CONTENTS
u1 Displacement [m]
ua Displacement [m]
ub Displacement [m]
v Coefficient [–]
v1 Strength reduction factor [–]
vmin Lower limit of shear capacity [MPa]
w Length [m]
wc Crack opening displacement [m]
wd Plastic softening compression [m]
x Global coordinate [–]
xdi Design value for a parameter [–]
xki Characteristic value for a parameter [–]
y Global coordinate [–]
z Internal lever arm [mm]
Greek symbols uppercase
∆P Load increment [N]
Greek symbols lowercase
α The inclination of the shear reinforcement [◦]
α Safety factor [–]
αcw The coefficient for the state of the stress in the compression chords [–]
αl Coefficient [–]
αR Sensitivity factor for the reliability of reticence [–]
β Reliability index [–]
β Reduction factor [–]
γC Partial factor for concrete [–]
γR Global resistance factor for material and geometrical uncertainties [–]
γRD Global resistance factor for model uncertainty [–]
γS Partial factor for steel [–]
γi Partial factor that are material dependent for the specific parameter [–]
0 Strain [–]
1 Strain [–]
2 Strain [–]
c Strain for compression [–]
c0 Strain at compressive cylinder strength [–]
d Strain [–]
t Strain for tension [–]
crack Non-linear strain from crack opening [–]
elastic Elastic strain from uncracked concrete [–]
total Total strain [–]
θ Model uncertainty [–]
θ The inclination the compression struts [◦]
ν Possion’s ratio [–]
σcc,1 Compression in strut [MPa]
σcc,2 Compression in strut [MPa]
σcc,A Compression at support [MPa]
σcp Compressive stress in concrete form axial loading or prestressing [MPa]
σt Stress in tension [MPa]
τ Shear stress [MPa]
φ Diameter of reinforcement bars [mm]
Abbreviations
1D One dimension
2D Two dimension
3D Three dimension
CPU Central Processing Unit CV Characteristic values
EC2 Eurocode 2
ECOV Estimation of Coefficient Of Variation
FE Finite Element
FEA Finite Element Analysis FEM Finite Element Method MC10 Model code 2010
MC90 Model code 1990
MV Mean values
NLFEA Non-Linear Finite Element Analysis PF Partial Factor
RC Reinforced Concrete
1. Introduction
Non-linear finite element analysis is a powerful tool that can be used in the design of reinforced concrete structures. However, when looking at safety formats, there does not exist a single method that should be used. Instead there exists several different methods where each one of these different methods have their own level of uncertainty in the finite element models.
1.1 Background
Non-linear finite element analyses is today widely used when assessing reinforced concrete structures. There are mainly three different methods that can be used for the calculation of the safety format for the non-linear finite element analysis. The method that are the most used today are the Partial factor method (in short form known as PF), since it can be found in most design codes. The other two methods that can be used are the Global resistance method with estimation of coefficient of variation of resistance (in short form known as ECOV) and the Full probabilistic analysis. The ECOV method is a simplified probabilistic method that can be used instead of a full probabilistic analysis ( ˆCervenka, 2013). All of three methods are describe in the fib Model Code 2010 and the ECOV method has been proposed to be included in the new version of Eurocode 2 (CEN, 2017), that is under development.
In the last years some evaluations of the different methods have been performed. In these evaluations different types of reinforced concrete structures have been analysed with the different methods and the results have been compared with eachother. The structures that have been analysed are previously tested in laboratories which makes it easier to verify the results given from the different simulations. From these different investigations a discussion has been raised about the model uncertainty. The different analyses have given different suggestion on how to consider the model uncertainty. It is clear that the model uncertainty is one of the biggest factors that influence the safety analyses (Schlune et al., 2012). For the ECOV method the model uncertainty is given as a factor with a value of 1.06 or higher according to the fib Model Code 2010 (FIB, 2013). Test done by Rijkswaterstaat (the Dutch Ministry of Infrastructure and Water Management) showed that model uncertainty factor on average should be around 1.15 which is a bit over 1.06 (Belletti et al., 2017). The model uncertainty factor is dependent on which type of failure mode the reinforced concrete structure achieves. For instance will a bending failure get a lower value for the model uncertainty factor while a shear failure will get a higher value ( ˆCervenka, 2013).
1.2 Aim
The aim of this master’s thesis is to investigate the model uncertainty factor in the ECOV method and try to propose a reasonable value for it. There will also be a comparison of the ECOV method with the PF method and with analytical calculations according to the current version of Eurocode 2. The investigation will be done by modelling different types of beams in the finite element programme ATENA 2D. The beams that will be used to perform the analyse have previously been tested in laboratories.
1.3 Limitations
In this thesis a few limitations were made, which mainly regards the NLFEA. All analyses were performed in 2D and not in 3D which excludes effects from actions acting perpen- dicular to the XY-plane (i.e. in the Z-direction) such as twisting or buckling. For more detailed specifications about limitations, assumptions and simplifications regarding the NLFEA see chapter 6. The thesis is limited to only analyse simply supported beams.
The number of simulated beams were limited to six divided between two different beam types, deep beams and slender beams. The study focuses on beams which failed by shear during previous experimental tests.
1.4 Outline
This master thesis is divided into eight chapters and three appendixes. Each chapter represents different phases during the work process. The different chapters are shortly described below.
Chapter 2 presents different types of safety formats that is available.
Chapter 3 introduces the concept of finite element method (FEM) and the implemen- tation of it in a finite element analysis (FEA) together with a description of a NLFEA.
Relevant theoretical information is also presented.
Chapter 4 describes the theory of shear in the different tested beam types.
Chapter 5 describes the beams that will be analysed. Their mechanical and geometrical properties are presented along with illustrations of the placement for the reinforcement bars.
Chapter 6 contains information of how the nonlinear finite element analyses in ATENA 2D were performed. The input values, used material models, boundary conditions and mesh type are presented. A brief theoretical description is also provided for each section.
Chapter 7 presents the results from the simulations of the beams in ATENA 2D. The results for the different safety formats and the model uncertainty factor are also presented.
Chapter 8 presents a discussion about the results and their validity.
1.4. OUTLINE
Appendix A presents the calculations performed for the ECOV method.
Appendix B presents the analytical calculations performed according to Eurocode 2.
Appendix C contains the illustration of the beams from the different simulations in ATENA 2D.
2. Safety formats
2.1 Partial safety factor method (PF)
The PF-method is a method that is used to obtain design values for different material dependent mechanical properties such as strength by reducing them with a coefficient (FIB, 2013). Design values are derived by dividing the materials characteristic value of the strength parameters by a material dependent partial factor, which usually is greater than 1, according to the following formula for a material i :
xdi = xki
γi (2.1)
Where:
xdi is the design value for a parameter.
xki is the characteristic value for a parameter.
γi is a material dependent partial factor for the specific parameter.
This method assigns the same reduced value of the mechanical properties to every point in the whole structure. These values are often extremely low and very conservative which may result in deviations regarding structural response (Hendriks et al., 2016). For example may an analysis, which utilizes the PF-method, show a failure in bending while it would break of a shear failure in reality (FIB, 2013).
2.2 Global resistance method with estimation of coefficient of variation of resistance (ECOV)
This method estimates a coefficient of variation from two separate non-linear finite ele- ment analyses with different input values regarding the materials mechanical properties.
One simulation uses the characteristic values and the other one uses the mean values ( ˆCervenka, 2008). The coefficient of variation, VR, is then produced according to follow- ing formula:
VR= 1
1.65lnPu,m
Pu,k (2.2)
Where:
Pu,m is the ultimate load for a structure obtained by a numerical simulation with mean material strength properties.
Pu,k is the ultimate load for a structure obtained by a numerical simulation with charac- teristic material strength properties.
The coefficient of variation is then used to produce a global resistance factor, γR, according to the following formula:
γR = eαR·β·VR (2.3)
Where:
αR is a sensitivity factor for the reliability of resistance and is usually set to 0.8.
β is a reliability index dependent on the safety class of the structure and a reference period. Values ranges between 3.3-5.2
The design resistance, Pd, is then calculated according to the following formula:
Pd= Pu,m
γRD · γR (2.4)
Where:
γRD is a model uncertainty coefficient which has a fixed value of 1.06.
The ECOV method is based on the assumption that the random distribution of the re- inforced concrete members resistance can be described by a two-parameter lognormal distribution which is bound downwards at the origin (Hendriks et al., 2016). This as- sumption originates from a theory formulated thru probabilistic studies. In the ECOV method the two random parameters are described by the mean resistance,Pu,m, and the coefficient of variation, VR (FIB, 2013).
2.3 Other methods
In addition to above mentioned methods there exists two other possible methods to use in non-linear FE analyses. These two methods are the Global resistance factor method and the Full probabilistic method.
2.3.1 Global resistance factor method
This method is similar to the earlier mentioned ECOV-method with the difference that this method only uses relative mean values for the mechanical properties in a non-linear FE analysis (CEN, 2004b). The relative mean values are specified as fym = 1.1 · fyk for the reinforcements yield strength and fcm = 0.85 · fck for the compressive strength of the concrete (CEN, 2004b). The results from the NLFEA is then used to obtain a design load according to:
Pd= Pu,m
γRD · γR (2.5)
Where:
Pu,m is the ultimate load for a structure obtained by a numerical simulation with relative mean material strength properties.
γRD is a model uncertainty coefficient which has a fixed value of 1.06.
γR is a global resistance factor which has a fixed value of 1.2.
2.4. MODEL UNCERTAINTY FACTOR
2.3.2 Full probabilistic method
The probabilistic method is based on statistics and probability theory. Simply described, this method uses a reliability index that is obtained from several randomized simulations to evaluate the safety. By doing this the reliability index considers all uncertainties due to random variations regarding materials mechanical properties, proportions and other effects (FIB, 2013). The reliability index is used in the following formula to obtain a design resistance:
Rd= R(α · β)
γRD (2.6)
Where:
R(α · β) is the resistance corresponding to reliability index, β.
α is a sensitivity factor for the reliability of resistance.
β is the reliability index.
γRD is a model uncertainty coefficient which has a fixed value of 1.06.
2.4 Model uncertainty factor
When the ECOV method was presented for the first time it did only have one uncertainty factor, which was the global safety factor of resistance γR ( ˆCervenka, 2008). It was later showed that γRonly considered the uncertainties in the properties of the material and the dimensions of the geometry. Therefore, the global safety factor of resistance should be formulated as γR= γR∗ · γRd as presented in the fib Model Code 2010 and also in the new proposition to the next version of the Eurocode 2 (CEN, 2017). Where γR∗ is the global safety factor of resistance that accounts for the uncertainties of the material properties as well as the geometrical dimensions and γRd is the global safety factor of resistance that take the model uncertainty into account. In the fib Model Code the suggested value for γRd is 1.06 for well validated numerical models (FIB, 2013). This value is also suggested to be used by Eurocode 2-2 (CEN, 2004b). These values are considered to be a bit low since they are based on the analyses of simple beams (i.e. beams with rectangular cross sections and no irregularities) only ( ˆCervenka, 2013). However, in the proposition to the next version of Eurocode 2 the following values for γRd are suggested (CEN, 2017):
• 1.06 for well validated numerical models.
• 1.10 for low-level validated numerical models.
• 1.35 for other models.
In the ECOV method, the global safety factor is assumed to be log-normal distributed with a coefficient of variation that are based on the mean and characteristic data values.
There exist a suggested alternative way to calculate a global safety factor of resistance.
In this alternative method the global safety factor is also assumed to be log-normal distributed. The difference is how the coefficient of variation is calculated. For the al- ternative method the coefficient of variation is calculated as follows (Schlune et al., 2012):
VR= q
Vg2+ Vθ2+ Vf2 (2.7)
Where:
Vg is the coefficient of variation for the geometry.
Vg is the coefficient of variation for the model uncertainties.
Vg is the coefficient of variation for the material uncertainties.
The problem with this method is that it has not considered the effects of the model val- idation for the model uncertainty ( ˆCervenka, 2013).
To be able to decide a model uncertainty factor a large statistical analysis has to be done.
The model uncertainty is often defined according to the following expression:
θ = Pexp
Pu,m (2.8)
Where:
θ is model uncertainty.
Pexp is the ultimate load for a structure obtained by an experimental test.
Pu,m is the ultimate load for a structure obtained by a numerical simulation with mean material strength properties.
The way of describing model uncertainty that is presented above is a standardised way of presenting it. If a model uncertainty for a structure is going to be compared with other structures and their model uncertainties, it is preferable to divide them in different groups after the type of failure that have occurred in the structure (Belletti et al., 2017).
3. Finite element method
3.1 Finite element analysis
Finite element analysis (FEA) is an analysis that utilizes the finite element method (FEM) to find approximate numerical solutions to mathematical problems described by differen- tial equations or integral expressions, also known as field problems (Cook et al., 2002).
This type of mathematical problems can for example be to calculate stress distribution in a beam or temperature distribution in a pipe. The FEM has several advantages in comparison with most other numerical methods especially regarding versatility. FEM can be used on any type of field problem and on any kind of geometrical shape and there are no restrictions regarding boundary or loading conditions. One FEM-model can consist of various components with different structural behaviour (bar, beam, plate, cable and fric- tion elements). One model can also consist of several materials with different mechanical properties, even within an element (Cook et al., 2002).
In a FEA, a FE-model is created by discretizing a mathematical model which is a sim- plified and idealized representation of a real structure. A FE-model consists of several smaller components called elements which are assembled in a specific pattern, called a mesh, to represent the real structure. All elements are connected to one another at points called nodes ( see figure 3.1). Each node has multiple degrees of freedom which simpli- fied can be compared with possible movements. These movements can be horizontal and vertical displacements, rotations and translations (Cook et al., 2002).
Figure 3.1: Schematic picture of a structure and its components in a FEA.
In FEA there is a general procedure on how to create the FE-model. First the structures geometry and reinforcement are defined, thereafter each structural member are assigned with their corresponding material and its mechanical properties. Boundary conditions and loads are then defined for the entire model. Finally, a mesh with a suitable element size is generated (Malm, 2016).
When performing FE-analyses there is a general solving procedure that simplified can be described by the seven following steps (Pacoste, 2017):
1. Assemble local stiffness matrix [k] for each element.
2. Transform and assemble all local stiffness matrices [k] to a global stiffness matrix [K]:
X
elements
[k] = [K] (3.1)
3. Determine the global load vector R from the local load vectorr, i.e. set all nodes to zero except the ones being exposed to external forces:
X
elements
{r} = {R} (3.2)
4. Reduce the system due to the boundary conditions by removing all fixated nodes from [K] and R.
5. Find global displacement vector D by solving the equation of equilibrium:
[K] · {D} = {R} → {D} = [K]−1· {R} (3.3) 6. Calculate section-forces (and/or stresses):
{˜r} = [NiViMiNiViMi]T (3.4)
{D} → {d} (3.5)
{˜r} = [k] · {d} (3.6)
7. Calculate reaction forces {RF }:
{RF } = [K] · {D} − {R} (3.7)
All capital letters symbolises global vector/matrices and all lower letters symbolizes local vector/matrices.
The quality of the analysis can be affected by three different errors. Discretization errors occur due to a poor mesh (i.e. a large element size) and can be reduced by using smaller elements. Modelling errors arise from simplification of the geometry and is reduced by adding more details. Numerical errors are frequently present since the analyses are
3.2. NON-LINEAR FINITE ELEMENT ANALYSIS
3.2 Non-linear finite element analysis
FEA where the FE-models are assumed to be linear (i.e. materials has an elastic be- haviour) usually obtain approximations that are considered sufficient for a lot of different structural problems. Nevertheless, it is quite common with problems that are non-linear, for instance in structural analyses where there exists several factors that impose a non- linear behaviour. These factors can be material yielding or creep, local buckling that may occur and gaps that may open or close (Cook et al., 2002).
There are three types of nonlinearity in structural mechanics (Cook et al., 2002):
• Material nonlinearity − in which material properties (e.g. elasticity, plasticity and creep) are functions of the state of stress or strain.
• Contact nonlinearity − in which the contact forces changes due to a change in the contact area between parts. This occurs when a gap between two adjacent parts opens or closes. There might also be a sliding contact which involves frictional forces.
• Geometric nonlinearity − in which deformation is large enough that the de- formed structural geometry need to be considered when formulating the equilibrium equations. Loads may also change direction when increasing.
Although NLFEA is a powerful tool, it has a few drawbacks. A non-linear behaviour is complicated to describe in a mathematical or numerical manner which in turn will make the analyses more difficult to perform. The analyses will also be more time-consuming due to the higher complexity and it will require a greater level of effort from the analysts to perform them. Consequently, this will increase the computational cost of a NLFEA (Cook et al., 2002). The higher complexity of a NLFEA can partly be derived from the fact that the stiffness and possibly the load as well are functions of displacements or deformations. In a mathematical manner, this can be illustrated with a structural equation:
[K] · {D} = {R} (3.8)
Where the stiffness matrix [K] and possibly the load vector R are functions of the de- formation vector D. Hence, a potential solution to the equation for D is obtained thru iterations. Since the magnitude of the load and the deformation changes for each iteration (in comparison with linear FEA where the load and the deformation increase linearly) it is necessary to divide the load application into several smaller load steps called increments.
In order to find a solution to the equation above the iterative procedure is repeated for each one of these increments. Thus, that is probably a significant contributing factor to the increase of computational costs in NLFEA (Cook et al., 2002). An iteration in a NLFEA almost has the same computational costs as a complete linear FEA. In addition
3.2.1 Iterative procedure
In a structure that has a linear behaviour there is always a force equilibrium. This means that if a structure is exposed to an external load, all internal forces and moments must be equal to the external load according to (Malm, 2016):
P − I = 0 (3.9)
Where:
P is the external load.
I is the internal forces.
When regarding the same case but from a non-linear perspective there can’t be a complete force equilibrium. This results in a modified approach where a force residual is calculated for each increment, a, according to:
P − Ia = Ra (3.10)
There are a few calculations preceding the calculation of the force residual and determine the non-linear response of a structure. If the structure is initially exposed to a small load increment ∆P , a displacement correction factor ca is determined by extrapolation using the tangential stiffness K0. The displacement correction factor ca is in turn used to update the initial displacement u0 to ua (see figure 3.2). After this, the internal forces Ia can be calculated for the corresponding displacement ua and in turn can the force residual be calculated (Malm, 2016).
Figure 3.2: Iteration of an increment in a NLFEA. Recreated from (Malm, 2016).
Since the force residual can’t be zero it is compared with a predefined tolerance value that basically will determine if the solution is good enough or if another iteration has to be made. A tolerance value are usually somewhere around 0.1-1% dependent on the magnitude of the load. If the residual is lower than the tolerance value the structure is considered to be in equilibrium and the calculations can proceed with the next load increment (Malm, 2016).
3.2. NON-LINEAR FINITE ELEMENT ANALYSIS
If the force residual is larger than the tolerance another iteration has to be made. The new iteration begins with updating the stiffness K0 based on uato the stiffness Ka which is used together with Ra to determine a new correction factor cb. The new correction factor will in turn give a new displacement and finally a new force residual that most likely are closer to equilibrium (see figure 3.3). This new force residual are compared with the tolerance value and the solution is either accepted or rejected. This iterative process is continued until all load increments have acceptable solutions (Malm, 2016).
Figure 3.3: Illustration of iteration procedure for an increment in a NLFEA. Recreated from (Malm, 2016).
3.2.2 Non-linear material behaviour
A material behaviour that is defined as linear is often described by another structural term namely elastic. Properties of elastic materials is for instance that the material will recover from deformation completely if the structure is unloaded. In opposite to elasticity, there is plasticity. A plastic material behaviour includes creep and yielding which results in permanent deformations when loaded. Many materials as reinforcement steel, concrete and wood are both elastic and plastic in different stages of loading. A linear FEA only considers the elastic part of a materials behaviour while a NLFEA considers both the elastic and the plastic part of a materials behaviour.
3.2.2.1 Non-linear behaviour of concrete
There are two fundamental methods that can be used to model cracking of concrete, the discrete crack method and the smeared crack method. In the discrete crack method, the elements are described with a linear elastic behaviour and the crack opening is described by an interface between the elements. In a model, the cracks will be visible and its fairly easy to determine the crack width since it is defined as the distance between to nodes on
appear in advance to be able to define the interfaces between the elements. This may be rather complicated in larger structures. The smeared crack method is also considered as the most common and is also the method that is used in this report (Malm, 2016).
Figure 3.4: Discrete and Smeared crack method. Recreated from (Malm, 2016).
3.2.2.1.1 Smeared crack models
As mentioned earlier the effect of a crack is smeared over a whole element and the total strain in one element consists of two parts according to:
total = elastic+ crack (3.11)
Where:
elastic is the elastic strain from the uncracked concrete.
crack is the non-linear strain from the crack opening.
In the smeared crack method there exists two different approaches that describes the development of cracks, the fixed crack model and the rotated crack model (see figure 3.5).
Figure 3.5: Fixed and rotated crack models. Recreated from (Malm, 2016)
In both crack models, the development of a crack initiates when the maximum principal stress is equal to the tensile strength of the concrete. The initiation of a crack occur in the same direction as the maximum principal direction. It is how the crack develops that
3.2. NON-LINEAR FINITE ELEMENT ANALYSIS
the subsequent loading, regardless of how the stress state changes. This results in the presence of shear stresses at the crack surface due to the different directions between the crack and the principal stress. In the rotated crack model, the crack direction in the element will rotate along with the principal stress direction which also rotates due to the stress state changes that follows from the subsequent loading. This will yield a crack surface without any shear stresses (Malm, 2016).
3.2.2.2 Non-linear behaviour of reinforcement
Reinforcement can be modelled in mainly two different ways, as discrete or as smeared (also cyclic reinforcement is available, but it is similar to the discrete) ( ˆCervenka et al., 2016). The smeared reinforcement creates special element layers where the amount of reinforcement is distributed (smeared) across an entire layer. The smeared layers are placed at the same vertical coordinat as the reinforcement bars are placed in reality and the thickness of the smeared layers are determined by the cross-sectional area of the bars since they should be equal ( ˆCervenka et al., 2016). The discrete reinforcement is represented as regular reinforcement bars modelled by 1D truss elements, which is basically a line for each bar. ( ˆCervenka et al., 2016).
4. Shear beam theory
4.1 Deep beams
According to the Eurocode a deep beam is defined as a beam which has a span length that is smaller than the height of the beam multiplied by three. Since a deep beam differs a bit from an ordinary beam geometrically it has a different kind of behaviour compared to an ordinary concrete beam.
In deep beams there are differences regarding the stress distribution in comparison with an ordinary beam. The stress distribution in an uncracked ordinary beam behaves linear in contrary to a deep beam where the stress distribution is non-linear. Generally, as the span/height relationship decreases the stress distributions linearity also decreases. The internal lever arm in the beam are also affected by the span/height relationship. When the deep beam has started to crack a redistribution of the stress in the beam occurs and the internal lever arm increases (Ansell et al., 2014).
It is common for deep beams to fail due to compression of concrete in the support region of the beam. This occurs since the beam can carry large loads in relation to its span length, reinforcement and section width, therefore a large force must pass through the support region that usually have small dimensions. In a deep beam the force in the rein- forcement are always constant, this differs from a ordinary beam where the force varies along the moment curve. Since the reinforcement force are constant through the whole beam the internal lever arm will vary along the moment curve instead. This leads to that the beam will behave as an “arch with tendon”. In an ordinary beam a part of the load can be carried by arch action but in deep beams the whole load can be carried by arch actions (Ansell et al., 2014).
According to Eurocode 2, deep beams should not be calculated as an “arch with ten- don”. Instead Eurocode 2 recommends that deep beam should be treated as areas of discontinuity, so called D-zones. These D-zones are areas where theory of linear strains cannot be applied. To deal with such zones a truss model should be used. A truss model is where external forces are carried by a system of compressed concrete struts and ties of reinforcement. An example of a truss model (or a strut-and-tie-model which is another name for it) for a simply supported beam are shown in figure 4.1. For a truss model to work it has to fulfill some conditions, as equilibrium has to be fulfilled, the truss model should imitate the “natural” behaviour of failure but it should not exclude any possible failure mechanisms. A difference between a truss model and an “arch with tendon” is that the truss model only considers the resultant of the different force while the “arch with tendon” consider the distribution of the external loads (Ansell et al., 2014).
Figure 4.1: Illustration of a truss model for a deep beam. Recreated from (Ansell et al., 2014).
Then both a truss model and an “arch with tendon” are dependent on the internal lever arm and the inclination angel of the compressive struts. Theses both parameters are related to each other according to the following:
z = a · tan (θ) (4.1)
Where:
z is the internal lever arm.
θ is the inclination angel of the compressive struts.
a is the distance between the support and the load point.
Both the internal lever arm and the inclination angel of the compressive struts can be chosen quite freely with in the different methods to calculate load capacity for deep beams. For an “arch with tendon” it is more advantageous to choose the internal lever arm freely and for the truss model it is the inclination angel or the compressive struts that are more advantageous to choose freely (Ansell et al., 2014).
4.2 Slender beams
Loading of a beam give rise of cracks, the first cracks are often vertical ones. When the load increases cracks start to develop near the support. Older cracks that already have appeared starts to bend towards the loading point. When the beam has reached shear failure, a crack will go from the tensile reinforcement to the compression zone at the top of the beam. This crack will be oriented diagonally, and it is denominated as a flexural shear crack. There also exists another type of shear crack that is called a web shear crack.
The difference between these two shear cracks are that web shear cracks develops from diagonal cracks near the neutral axis, while a flexure shear crack develops from bending cracks (Ansell et al., 2014).
A flexural shear failure is influenced by several different factors such as shear slenderness, concrete strength, reinforcement content, beam height and reinforcement bond. The in-
4.2. SLENDER BEAMS
The reason for this is that an arch effect develops in shorter beams that makes it possi- ble for the beam to carry larger loads. That the concrete strength would influence the shear failure is pretty obvious for a concrete beam without any shear reinforcement. The influence from the reinforcement content can be described by three different factors. The first one is that the height of the compression zone is increased, the second one is that it will be more difficult for cracks to open and the last one is that the resistance of the reinforcement against transvers displacement are increased. When the height of the beam increases the load carrying capacity of the beam will decrease. The load carrying capacity of the beam will be higher for a low reinforcement bond than for a high reinforcement bond if the reinforcement is fully end-anchored. This is however not practical since a special anchoring device is necessary to achieve full end-anchoring which reinforcement with a high bond does not need (Ansell et al., 2014).
To calculate the load carrying capacity for flexural shear failure the following formula from Eurocode 2 will be used:
VRd,c = (CRd,c· k ·p3
100ρ1· fck+ k1· σcp) · bw· d (4.2) Where:
CRd,c= 0.18γ
c
γc= 1.5 is a partial coefficient for concrete.
k = 1 + q200
d ≤ 2.0 ρ1 = bAsl
w·d ≤ 0.02
Asl is the area of the tensile reinforcement.
fck is the characteristic compressive strength of the concrete.
k1 is a coefficient.
σcp = NAEd
c
NEd is the normal force from loading or pre-stressing of the beam.
Ac is the area of the concrete cross section.
bw is the smallest width of the cross section.
d is the effective height.
The formula above can be calculated with a minimum value as follows:
VRd,c = (vmin+ k1· σcp) · bw· d (4.3) Where:
vmin = 0.0035pk3· fck
If the load point is placed close to the support (a ≤ 2d) part of the load will be carried by arch action (as mentioned above). This will give a higher loading capacity for the beam.
Since loads that are carried partly of the arch action are not contributing to the flexure shear failure the ultimate load can be reduced. The reduction factor (β) are decide of the ratio 2da (Ansell et al., 2014).
