Dynamic Blast Load Analysis using RFEM: Software evaluation

Dynamic Blast Load Analysis using RFEM

Software evaluation

Oskar Dädeby

Civil Engineering, master's level (120 credits) 2021

Luleå University of Technology

Department of Civil, Environmental and Natural Resources Engineering

Abstract

The purpose of this Master thesis is to evaluate the RFEM software and determine if it could be used for dynamic analyses using blast loads from explosions. Determining the blast resistance for a structure is a growing market and would therefore be beneficial for Sweco Eskilstuna if RFEM could be used for this type of work. The verification involved comparing the RFEM software to a real experiment which consisted of a set of blast tested reinforced concrete beams. By using the structural properties from the experiment project with the experiment setup the same structure could be replicated in RFEM. RFEM would then simulate a dynamic analysis loaded with the same dynamic load measured from the experi- ment project in two different dynamic load cases caused by two differently loaded explosions. The struc- tural response from the experiment could then be compared to the response simulated by the RFEM software, which consisted of displacement- and acceleration time diagrams. By analysing the displacement and acceleration of both the experiment and the RFEM software the accuracy was determined, and how well RFEM preformed the analysis for this specific situation. The comparison of the displacement and acceleration between the experiment and RFEM was considered acceptable if the maximum displace- ment was consistent with the experiments result and within the same time frame. The acceleration was considered acceptable if the initial acceleration was consistent with the experiment result. These criteria needed to be met for the verification that RFEM could simulate a dynamic analysis. If the software managed to complete a dynamic analysis for two dynamic load cases, then the software could be evaluated which consisted of determining if the post blast effects could be determined and if the modelling method was reliable.

The acceleration from RFEM were in good agreement with the experiment test at the initial part of the blast, reaching a close comparison for both load cases after 3 ms. Then the RFEM acceleration had a chaotic behaviour reaching no similarities for the duration of the blast. The displacement managed to get a close comparison of the maximum displacement with a margin of 0,5 mm for both load cases within a 1 ms time margin. RFEM managed in conclusion to simulate a blast load analysis, the displacement and acceleration gave acceptable results according to the criteria.

With the method chosen a fast simulation was achieved and with the same model complying with two different load cases for the same model gave indication that the first result was not a coincidence. The steps taken in the modelling method was straight forward, but two contributing parameters were deter- mined to devalue the reliability. First parameter was the material model chosen for the concrete, which was chosen to a plastic material model. The two optional material model’s linear elastic and non-linear elastic both caused failed simulations. Also, the better model for the material model would have been a diagram model which insured that the concrete lost is capacity in tension with maximum capacity, but this was not available in a dynamic analysis with multiple load increments. Which is the reason why a plastic material model was chosen for the concrete. The second reason was the movement of the beam in the supports. This data was not recorded in the experiment but was determined to be a contributing part of the test. This however gave big differences of the result depending on how much the beam could move. In the end the best possible result was chosen to comply with the first load case where the same RFEM model was used in the second test. The second load case showed just as good results as the first load case, but with the big variation in results depending on the movement of the beam in the supports made this part unclear.

For the evaluation the question if the RFEM could provide a post blast analysis needed to be addressed, where the answer is no. The failure mode was chosen to comply with the choice of modelling method which required the analysis of the plastic strain in the reinforcement bars. This information was not available using the add-on module DYNAM-PRO and could therefore not provide the answer if the model structure resisted the blast.

For future work of this master thesis is to build a model that would give a more detailed post blast analysis, where this thesis was made to test the software. For this more work would be necessary by the creators Dlubal to further improve the add-on-module, which involves more extractable results and more detailed tools when using a dynamic load case, where some important functionality is only usable in a static load case. Other than that, RFEM managed to complete the dynamic analysis, and with further improving of the modelling method a more detailed analysis can be made and then be usable in real projects in the future.

Keywords: Blast load, Dynamic, RFEM, time history, experiment, time step

Preface

When designing a structure for dynamic blast loads finite elements programs are recommended to use.

Sweco Eskilstuna use a finite element software called RFEM created by a German company Dlubal which might have had the potential to execute these kinds of dynamic analyses. RFEM had an add-on module called DYNAM-PRO which preforms dynamic analysis, but it was still untested by anyone in the Eskils- tuna office.

This Master thesis was the final part in my master’s degree in Structural and civil engineering at the division of building design and production and- division of civil engineering and natural resources at Luleå Technical university (LTU). This thesis started in the beginning of August 2020 and was completed in 30-05-2021. This thesis was written in collaboration with Sweco Eskilstuna where I had had my internship from June to December 2019.

This thesis was written during the eye of the Corona pandemic which certainly provided its challenges.

I want to give a special thanks to Frank Axhag and Carolina Eklöf at Sweco Eskilstuna who has guided me through this entire process. I also want to thank my supervisor at LTU Gabriel Sas, other teachers and people as well as my parents who have helped me with my work.

Luleå, 30-05- 2021 Oskar Dädeby

Table of content

ABSTRACT ... V PREFACE ... VII

1 INTRODUCTION ... 6

1.1 Background ... 6

1.2 Aim and Objective ... 6

1.3 Limitations ... 7

2 THEORY OF BLAST LOADS ... 9

2.1 Blast loads mechanics ... 9

2.1.1 Definition of explosion ... 9

2.1.2 Pressure-time profile and Reflected pressure ... 9

2.1.3 Impulse load and Pressure load ... 10

2.2 Dynamic mechanics ... 11

2.2.1 Dynamic formulation ... 11

2.2.2 SDOF elastic response to impulse load ... 12

2.2.3 SDOF plastic response to impulse load ... 13

2.3 Time history analysis... 14

2.3.1 Explicit ... 14

2.3.2 Implicit Newmark method ... 15

2.3.3 Non-linearities in DYNAM-PRO ... 16

2.3.4 Damping ... 17

3 METHOD FOR MODELLING ... 19

3.1 Experiment project ... 19

3.1.1 Introduction to the experiment project ... 19

3.1.2 Test objects ... 19

3.1.3 Blast layout ... 20

3.1.4 Results from tests ... 21

3.2 RFEM Simulation introduction ... 22

3.2.1 Variation of modelling ... 22

3.2.2 Variation of material model ... 23

3.2.3 Material description ... 25

3.2.4 Support boundary conditions ... 25

3.2.5 Load ... 26

3.2.6 FE-mesh settings ... 27

3.3 RFEM Add-on module DYNAM-PRO ... 28

3.3.1 Introduction to the add-on module ... 28

3.3.2 Time history analysis ... 28

3.3.3 Time diagram ... 28

3.3.4 Dynamic load case ... 28

4 RESULTS ... 30

4.1 Post blast results ... 30

4.2 Test results ... 30

4.2.1 Variation of modelling ... 30

4.2.2 Variation of material model ... 31

4.2.3 Variation of support boundary conditions ... 32

4.2.4 Variation of FE-mesh settings ... 33

4.2.5 Variation of time history model ... 34

4.3 Experiment result and Final model ... 36

5 ANALYSIS AND DISCUSSIONS... 39

5.1 Simulation ... 39

5.1.1 Modelling ... 39

5.1.2 Material model ... 40

5.1.3 Support boundary conditions ... 40

5.1.4 FE-mesh ... 41

5.1.5 Time history ... 41

5.2 Simulated Results ... 41

6 CONCLUSIONS AND FUTURE WORK ... 43

7 REFERENCES ... 45

APPENDICES A Step by Step Guide for Modelling A-1 B RFEM Concrete and Steel Parameters A-13 C Recorded Time Diagrams from Experiment Test D2 & D3 A-15 D Pressure Data A-17 E Photos from experiment Test A-19

Notations and abbreviations

Notations

Symbol Description M, m Mass (matrix) (𝑘𝑔) D Damping matrix (-) K, k Stiffness (matrix) (𝑚⁴) 𝑢̈, 𝑎 Acceleration (𝑚 𝑠⁄ ) ² 𝑢̇, 𝑣 Velocity (𝑚 𝑠⁄ ) u, d Displacement (𝑚) p(t), F Excitation force (𝑁) R Resistance force (𝑁) 𝑊_𝑦 Outer work (𝐽) 𝑊_𝑖 Inner work (𝐽)

𝜔 Angular frequency (𝑅𝑎𝑑 𝑠⁄ ) h Time step (𝑠)

T Wave- length time (𝑠) t Time (𝑠)

𝜉 Damping ratio (-)

𝛼 Damping mass coefficient (-) 𝛽 Damping stiffness coefficient (-) Abbreviations

FEM Finite element method SDOF Single degree of freedom NVC Natural vibration case DLC Dynamic load case SLC Static load case

List of Figures

Figure 1: Pressure-time profile of the explosion wave (Draganic & Sigmund, 2012) (Johansson, 2012) 10

Figure 2: Variation of blast wave pressure with distance (Draganic & Sigmund, 2012) ... 10

Figure 3: damage curves based on Impulse load or pressure load, in Swedish (Johansson, 2012) ... 11

Figure 4: Load variation dependent of duration (a) Impulse load (b) Pressure load (Johansson, 2012) .. 11

Figure 5: single degree of freedom system for dynamic definition, SDOF (Johansson, 2012) ... 12

Figure 6 elastic response: (a) SDOF (b) Force- displacement relationship (c) equalization of energy .... 13

Figure 7 plastic response (a) SDOF (b) Force- displacement relationship (c) equalization of energy (Johansson, 2012) ... 14

Figure 8 Relationship between Rayleigh coefficients (Dlubal Software GmbH 2020, 2020) ... 18

Figure 9 Distribution of reinforcement (Magnusson & Hallgren, 2003) ... 20

Figure 10 Shock tube setup (Magnusson & Hallgren, 2003) ... 21

Figure 11 placement of gauges (Magnusson & Hallgren, 2003) ... 21

Figure 12 Model structure with (gray surface) concrete surface, (blue bars) extruded reinforcement, (green pyramids) nodal supports, (red dots) nodes. ... 23

Figure 13 Plastic material model of steel ... 24

Figure 14 Plastic material model of concrete ... 25

Figure 15 Distribution of load ... 26

Figure 16 Reflective pressure of beam D3 ... 27

Figure 17 Reflective pressure of beam D2 ... 27

Figure 18 Displacement of beam without reinforcement ... 30

Figure 19 fracture critera for SLC ... 31

Figure 20 Displacement of beam with an elastic material model ... 31

Figure 21 Normal force of members during time ... 32

Figure 22 Displacement of fully rigged beam ... 32

Figure 23 Displacement of beam with variating stiffness factor in supports ... 33

Figure 24 Maximum movement of supports ... 33

Figure 25 Displacement with variating FE mesh size ... 34

Figure 26 Linear time history model ... 34

Figure 27 Displacement of Explicit vs Implicit ... 35

Figure 28 Displacement comparison between RFEM and beam D3 ... 37

Figure 29 Acceleration comparison between RFEM and beam D3 ... 37

Figure 30 Displacement comparison between RFEM and beam D2 ... 38

Figure 31 Acceleration comparison between RFEM and beam D2 ... 38

Figure 32 RFEM image of parallel surface and steel bars ... 39

Figure 33 RFEM error message using diagram type materal model ... 40

List of Tables

Table 2 Steel properties ... 20

1 Introduction

1.1 Background

An accidental load is an action upon a structure caused by human error or unforeseen circumstances (EN 1991-1-7, 2006). To protect structures against accidental loads EN 1997-1-7 describe three consequence classes describing the importance of the building protection. The lowest class is CC1 which require no specific consideration and the highest class is CC3 which requires special cases and might need advanced methods like dynamic analyses, non-linear models and interaction between the load and the structure.

This standard does not however deal with accidental actions of type blast loads which might be induced by accidental explosions, warfare, sabotage or terrorist activities.

Certain structures of the highest consequence class need protection against blast loads, otherwise an entire society can be negatively affected by the destruction of an important piece of infrastructure. Sweco Eskils- tuna had the intention to improve their knowledge about blast loads, specifically in combination with using their standard finite element software RFEM. RFEM works like many other finite element soft- ware’s used for structural engineering but in addition it also has an add-on module named DYNAM- PRO which performs a dynamic analysis. Since blast loads are categorized as a dynamic load case this add-on module, DYNAM-PRO was tested to assess if it could be used for analyzing a structure affected by a blast load. If RFEM could perform this type of dynamic analysis the need for a separate program would not be necessary for future work, where certain software’s can be restrictive and expensive.

To see if the software could perform this type of analysis RFEM had to be verified. For this a method of comparison was chosen to compare the software results to a real experiment. The project chosen was the scientific study published by the Swedish Defense Research Agency performed by Magnusson, J and Hallgren, M (Magnusson & Hallgren, 2003). They had performed a blast test on a set on reinforced concrete beams and provided results in form of deflection- time diagrams and acceleration- time diagrams.

RFEM can produce these types of diagrams which can be used to compare experimental results with simulated results from RFEM. If the results from RFEM compares well with the experimental results, the software could be evaluated to potentially be taken into consideration for blast loads analysis when designing real structures.

1.2 Aim and Objective

The aim of this master thesis is to assess if the response of the experimental test can be replicated in RFEM and thereby determine if the software can be used for blast load analysis. This is determined by comparing the acceleration and displacement during a time period. The comparison of the displacement and accel- eration between the experiment and RFEM was considered acceptable if the maximum displacement was consistent with the experiments result and within the same time frame. The acceleration was con- sidered acceptable if the initial acceleration was consistent with the experiment result. Two load cases were analysed for this and if both follow the above criteria then the question if RFEM and the experiment had a good comparison could be answered. When the results from the experiment and RFEM has been compared the software could be evaluated. For the evaluation the following questions needed to be answered.

• Can information be extracted from RFEM to give an indication if the structure resisted the blast?

• Did the method of modelling create problems which devalued the reliability of the simulated re- sult?

To reach the goal the blast load mechanics and the dynamic mechanics were assessed followed by under- standing how the experiment was performed. Then the structure model was created in RFEM and pa- rameters within the RFEM model were varied to replicate the response from the experiments. The varied parameters consisted of initially trying for example a type of structure model, or a type of material model and then compare the result with the experiment project. By analysing the response, the model could be

improved and the correct path for each step could be taken. When all steps provided a unified acceptable result, a discussion was made to get an accurate representation between the experiment and RFEM me- chanically. Lastly the final simulations were run with two DLCs using the same structure.

The hypothesis was that the results would have some similarities, either the correct deflection mode or values when comparing the time- diagrams. Some differences were expected but the work was strived towards to have the results from the RFEM analysis, to coincide with the experimental results, as close as possible while having good arguments and assumptions for the steps taken.

1.3 Limitations

The focus of the thesis was to analyse the similarities between the response of RFEM and the experimental but not necessarily determining the post blast effects of the structure, after the duration of the blast load.

A professor at LTU (Zhang, 2020) said that many simulated blast tests measured the cracking in the concrete. In this thesis the similarities of the displacement and acceleration were basis of the analysis.

However, the plastic strain in the reinforcement would have been of interest to analyse, to see if the rebar failed or not and give a simple indication if the structure resisted the blast. Unfortunately, no data for steel strain were available using DYNAM- PRO, so no failure modes in RFEM could be analysed in this thesis.

The limitations of the scope of the assignment was linked with the experiment data available, which limited the options of analysed structures to the beams in the experiment project. This also limited how accurate the most accurate model could be as some information needed was not provided. As will be seen in this thesis shell (2 dimensions) and beam (1 dimensions) elements were used instead of a solid element (3 dimensions). A solid model was used initially but took too much time to optimise correctly but was initially the preferred method.

2 Theory of blast loads

2.1 Blast loads mechanics 2.1.1 Definition of explosion

All matter has potential energy which is defined by energy which has the potential to convert to some- thing else. After ignition of an explosive material a rapid release of potential energy occurs which creates a hot dense gas which expands the material (Draganic & Sigmund, 2012). This expansion causes high over-pressure which is the source of a powerful blast wave which can expand over a large volume at incredible high speed.

The initial ignition of the explosive material causes a chain reaction which spread through the explosive material (Physics, 2018). The first chemical reaction releases energy which triggers the next chemical reaction, the more rapidly this goes the stronger the reaction will become. The reasoning why an explo- sive material can release its energy rapidly has to do with its chemical structure. A log of wood for example doesn’t explode because it has no inbounded oxygen, which creates the needs to outsource oxygen to create is chemical reactions. This reaction happens slowly which creates a non-destructive outcome, even though it has a lot of potential energy to be released. For an explosive material like nitro-glycerine which doesn’t need an outside reactant the chain reaction happens rapidly which creates the sudden surge of energy and makes the experience destructive.

Dependent on which explosive material that’s in question the reaction acts differently which makes some explosives more destructive than others. To categorise the different explosives a correction factor for the explosive mass has been implemented based on the original explosive material TNT. PETN for example has a correction factor of 1,27 making 1 kg of PETN equivalent to 1,27 kg of TNT.

2.1.2 Pressure-time profile and reflected pressure

The acting pressure from an explosion will have an arrival time tA which is the time it takes for the blast wave to reach the subject (Draganic & Sigmund, 2012). Initially the pressure will be at its maximum and thereafter diminishes until it reaches the ambient(atmospherically) pressure. This is called the Positive phase, see figure (1) which happens during a time period, t0. This is followed by the negative phase which reaches a minimum point and turns back to its original state at the ambient pressure during a time period t^-0. This is because of the law of conservation of energy which carries the momentum to swing the structure in the opposite direction of the blast until it returns to its original state. This diagram is called the pressure-time profile which describes the duration of a blast on a structure. This model shows that a blast load has a dynamic behaviour and will linger until the energy is fully absorbed by the structure.

Figure 1: Pressure-time profile of the explosion wave (Draganic & Sigmund, 2012) (Johansson, 2012)

A detonation will increase its pressure fast, spanning just a couple of milliseconds before a nearby target is reached (Löfquist, 2016). A free air explosion will have a spherical shape if no constraints is added expanding its volume until it reaches ambient pressure. The pressure wave reduces its strength and speed as distance is increased, see figure (2). Therefore, a nearby blast will be more destructive than a blast further away.

Figure 2: Variation of blast wave pressure with distance (Draganic & Sigmund, 2012)

As of a normal open-air detonation the wave spreads without any wave amplifications (Draganic &

Sigmund, 2012). But if the detonations happen near the ground the pressure wave reflects away from the ground in various ways. If reflections are acting in a detonation the initial pressure will be significantly increased, this is defined as reflective pressure which is the wave amplification of the wave reflective of nearby surfaces. This can also be a contributing factor where the pressure wave reflects from nearby structure to either reach point around the front of a structure or multiple reflective waves that interacts, creating a slight increase at some point after the initial pressure.

2.1.3 Impulse load and Pressure load

A structure can be damaged by the initial pressure, but the impulse can also be damaging in the right circumstances (Johansson, 2012). Impulse is the integral of the pressure-time profile, see figure (1) which can be damaging if the blast has a long duration. This is mostly common in tunnel systems where the reflective behaviour of a blast wave can create a long duration of a blast. Figure (3) describes the damage zones based on the Pressure F (or p) and Impulse I. The green field represents the undamaged zone and blue the damaged zone. If the blast has a short duration the pressure required to reach the blue zone is high, while if the blast has a low intensity it needs a long duration to reach the blue zone.

Figure 3: damage curves based on Impulse load or pressure load, in Swedish (Johansson, 2012)

These loads are defined as ideal impulse load and pressure load. Impulse load is indexed with Ik and is said to be high intensity and short duration and pressure load is indexed with Fk which is low intensity and long duration, see figure (4).

Figure 4: Load variation dependent of duration (a) Impulse load (b) Pressure load (Johansson, 2012)

2.2 Dynamic mechanics 2.2.1 Dynamic formulation

The definition of dynamics is a property which depends on time, in contrast to static which is constant (Johansson, 2012). A traditional model to explain structural dynamics is d’Alembert’s principle, see figure (5), which is a SDOF linear elastic system exposed to an excitation force. The resisting forces of this system is based on the mass, the stiffness and energy loss mechanics called damping.

Figure 5: single degree of freedom system for dynamic definition, SDOF (Johansson, 2012) By equalising the resisting forces with the acting forces, equation (2.1) can be described. Here the iner- tia force is a product of mass matrix and acceleration, the damping force is the product of damping con- stant matrix and velocity and the elastic force is the product of spring stiffness matrix and displacement.

𝑀𝑢̈ + 𝐷𝑢̇ + 𝐾𝑢 = 𝑝(𝑡) (2.1)

The effectiveness of a structures probability to resist a blast load has to do with its energy absorbent properties. By determining the inner and outer work caused by the blast load a breaking point can be determined which describes where the structure starts to overpower the load at a certain displacement of the structure. The outer work can be described as the initial kinetic energy dependent on the impulse load Ik and mass m in equation (2.2).

𝑊_𝑦= 𝐼_𝑘²

2𝑚 (2.2)

As said the breaking point can be found by equalising the inner and outer work in equation (2.3)

𝑊_𝑦= 𝑊_𝑖 (2.3)

The energy absorbent properties of a structures corelates to the acting force and displacement, which is the two main components of work (Johansson, 2012). In many cases when a structure is loaded with a high intensity impulse load it is hard to fully resist a structure loaded by an equivalent static force. The load is higher than its static counterpart, but it only exists for a very short period, therefore it’s important to use a structures deformation capability to counteract the increase in load. If the structure has good deformation capabilities it can result in big energy absorbent capabilities, while a structure with high resisting force with bad deformation capacity could be weak against a blast.

2.2.2 SDOF elastic response to impulse load

If the structure has a linear elastic behaviour the inner resisting force can be describes in equation (2.4)

𝑅(𝑢) = 𝐾𝑢 (2.4)

Where K is constant stiffness and u is displacement. The elastic inner work can then be described in equation (2.5) which is then combined with equation (2.4)

𝑊_𝑖 =𝑅(𝑢_𝑒𝑙)𝑢_𝑒𝑙

2 =𝐾𝑢_𝑒𝑙²

2 (2.5)

When the inner and outer energy is equal the inner work starts to take over and the elastic displacement stops. Where the elastic displacement can be described in equation (2.6)

𝑢_𝑒𝑙 = 𝐼_𝑘

𝑚𝜔 (2.6)

Where 𝜔 is the angular frequency which is defined in equation (2.7)

𝜔 = √𝐾

𝑚 (2.7)

This brief explanation of inner and outer energy in a SDOF elastic system describe the theoretical way for how the structure deforms during a dynamic load. As can be seen in figure (6) the inner work in- creases with linearly elastic displacement. The elastic response doesn’t have a maximum resisting force.

Therefore, the structure can indefinitely deform with increasing force and return to the original state when unloaded. In (c) the equalization of work can be seen where after deformation uel is reached the deformation stops.

Figure 6 elastic response: (a) SDOF (b) Force- displacement relationship (c) equalization of energy 2.2.3 SDOF plastic response to impulse load

If the structure has a plastic behaviour the inner resisting force can be describes in equation (2.8)

𝑅(𝑢) = 𝑅 (2.8)

Since R is constant the inner work can be described in equation (2.9)

𝑊_𝑖 = 𝑅(𝑢_𝑝𝑙)𝑢_𝑝𝑙= 𝑅𝑢_𝑝𝑙 (2.9)

Using equation (2.8) in combination with equation (2.9) the plastic deformation needed for equilibrium to be reached can be described in equation (2.10)

𝑢_𝑝𝑙= 𝐼_𝑘²

2𝑚𝑅 (2.10)

The plastic response in contrast to the linear elastic response does have a maximum resisting force avail- able, see figure (7). This is a more realistic view of the behavior of a structure where the deformation upl

will be reached and not providing any more work after that. This is the balance of resisting force and plastic deformation which describe the energy absorbent capabilities. As said the structure does not always resist the initial excitation as seen in (c) where the force is much larger than the resisting force. But depending on the energy absorbent properties the work can be equalized as the structure deforms. This require that the structure can deform enough to create the work needed, creating a more ductile behavior when designing against blast loads.

Figure 7 plastic response (a) SDOF (b) Force- displacement relationship (c) equalization of energy (Johansson, 2012)

2.3 Time history analysis

There are many ways of solving a system affected by dynamic formulations, two of this is the implicit Newmark method and the explicit method, which roughly means the “implied” and “clear” method (Clough & Penzien, 1993).

2.3.1 Explicit

The explicit method starts by reformulating equation (2.1) to have its initial properties at timestep t = 0 described in equation (2.11)

𝑀𝑢̈₀+ 𝐷𝑢̇₀+ 𝐾𝑢₀= 𝑝(0) (2.11)

The initial acceleration is then singled out from equation (2.11) in equation (2.12) 𝑢̈₀= 1

𝑀[𝑝₀− 𝐷𝑢̇₀− 𝐾𝑢₀] (2.12)

The method gone through here is the second central difference method. This is an approximated method which initiates with approximate the velocity in the middle of each time step before and after t = 0 described in equation (2.13)

𝑢̇_{−1 2⁄} =̇𝑢₀− 𝑢₋₁ ℎ 𝑢̇1

⁄2=̇𝑢₁− 𝑢₀

ℎ (2.13)

Where h [s] represents the time step between t0 and t1, which is the increments of the analysis. This is then transformed to an expression of the acceleration described in equation (2.14)

𝑢̈₀=̇𝑢̇1

⁄2− 𝑢̇_{−1 2⁄}

ℎ =̇ 1

ℎ²(𝑢₁− 𝑢₀) − 1

ℎ²(𝑢₀− 𝑢₋₁) (2.14)

Equation (2.14) is simplified to equation (2.15) 𝑢̈₀= 1

ℎ²(𝑢₁− 2𝑢₀+ 𝑢₋₁) (2.15)

Combining equation (2.12) and (2.15) gives the following expression in equation (2.16) 𝑢₁− 2𝑢₀+ 𝑢₋₁=ℎ²

𝑀(𝑝₀− 𝐷𝑢₀̇ − 𝐾𝑢₀) (2.16)

The displacement is singled out from equation (2.16) described in equation (2.17)

𝑢₁=ℎ²

𝑀(𝑝₀− 𝐷𝑢̇₀− 𝐾𝑢₀) + 2𝑢₀− 𝑢₋₁ (2.17) Missing variable from this expression is the displacement at time t = -1, which can be determined from the following velocity expression described in equation (2.18)

𝑢̇₀ =𝑢₁− 𝑢₋₁

2ℎ (2.18)

The displacement at time t= -1 is singled out in from equation (2.18) described in equation (2.19)

𝑢₋₁= 𝑢₁− 2ℎ𝑢̇₀ (2.19)

Substituting equation (2.19) into equation (2.17) gives the following expression described in equation (2.20)

𝑢₁= 𝑢₀+ ℎ𝑢₀+ ℎ²

2𝑚(𝑝₀− 𝑐𝑢̇₀− 𝑘𝑢₀) (2.20)

The velocity at t = 1 is also required which is achieved by assuming that the average velocity at time t = 1 and t = 0 is equal to the velocity within the time step described in equation (2.21)

1

2(𝑢̇₀+ 𝑢̇₁) =𝑢₁−𝑢₀

ℎ (2.21)

Where the velocity at t = 1 can be singled out from equation (2.21) described in equation (2.22) 𝑢̇₁=2(𝑢₁− 𝑢₀)

ℎ − 𝑢̇₀ (2.22)

The second central difference method uses equations (2.20) & (2.22) to approximate from one step to the next throughout the required timespan. This is a simple explicit step-by-step method, which has the requirement for an appropriate time step for it to be stable. It’s said to be conditionally stable and can fail if the time step is not made small enough. Rule of thumb is to use the following expression of the time step which represents around three-time steps per vibration period, if the dynamic force is harmonic.

ℎ 𝑇≤1

𝜋 (2.23)

2.3.2 Implicit Newmark method

The implicit method uses a method of integration which calculates values based on the initial and final conditions for each time step (Clough & Penzien, 1993). The essential concept is represented by the velocity integration in equation (2.24) and displacement integration in equation (2.25)

𝑢̇₁= 𝑢̇₀+ ∫ 𝑢̈(𝑟)𝑑𝑟

ℎ

0

(2.24)

𝑢₁= 𝑢₀+ ∫ 𝑢̇(𝑟)𝑑𝑟

ℎ

0

(2.25)

Theses integrals express the final velocity and displacement dependent of the acceleration and velocity.

This method requires an assumption of the acceleration progression during the timesteps. A step-by-step formulation was produced by Newmark also called the Newmark implicit method, since the assumption of acceleration makes this an implicit method over an explicit method. In the Newmark formulation, the basic integration equations for the final velocity and displacement are expressed as follows:

𝑢̇₁= 𝑢̇₀+ (1 − 𝛾)ℎ𝑢̈₀+ 𝛾ℎ𝑢₁ (2.26) 𝑢₁= 𝑢₀+ ℎ𝑢̇₀+ (1

2− 𝛽) ℎ²𝑢̈₀+ 𝛽ℎ²𝑢̈₁ (2.27)

The 𝛾 in equation (2.26) and 𝛽 in equation (2.27) provides a linearly varying weighting for the final acceleration and velocity.

Studies showed that the gamma factor controlled the amount of artificial damping induced. If 𝛾 has a value of ½ no artificial damping is implied in the system, which is the recommendation for SDOF analysis using this method. If the beta value is set to ¼ it shows that the formulation reduces directly to the expression shown for the final velocity and displacement in equation (2.26) and (2.27). Therefore, this variation is called Newmark 𝛽 = ¼ method which induces a constant average acceleration.

Another variation is 𝛽 = 1/6 with gamma still on 𝛾= ½, the expression of the final velocity and displace- ment is described in equation (2.28) and (2.29)

𝑢̇₁= 𝑢̇₀+ℎ

2(𝑢̈₀+ 𝑢̈₁) (2.28)

𝑢₁= 𝑢₀+ 𝑢̇₀ℎ +ℎ²

3 𝑢₀+ℎ²

6 𝑢̈₁ (2.29)

The derivation can also imply that the acceleration is linear during the time step. Which makes the Newmark b = 1/6 method known as the linear acceleration method. Both these methods are widely used while the b = ¼ is said to be unconditionally stable while the b = 1/6 has a time step requirement described in equation (2.30)

ℎ 𝑇≤√3

𝜋 (2.30)

2.3.3 Non-linearities in DYNAM-PRO

The two methods explicit in section 2.3.1 and implicit in section 2.3.2 is both linear solvers, but a dy- namic blast analysis require a non-linear model because of the constant changing non-linearities. Both solvers can perform nonlinear solvers which preform a large strain deformation analysis. For a non-linear formulation it is said to only keep the physical properties constant for short increments of time. By refor- mulating the terms to incremental formulations, the physical properties can have a non-linear outcome.

As said in section 2.3.1 the explicit method in conditionally stable depending on the choice of time step.

Non-linear analysis in RFEM requirement the time increments Δt(h) to be less than the stable time increment Δtstable. An estimation of the stable time increments is described in equation (2.32) with finite element length L^e and the dilatational wave speed for a linear elastic material mode cd.

∆𝑡 < Δ𝑡_{𝑠𝑡𝑎𝑏𝑙𝑒}=𝐿_𝑒

𝑐_𝑑 (2.31)

Where the dilatational wave speed depends on the Young’s modulus E and material density ρ described in equation (2.32)

𝑐_𝑑= √𝐸

𝜌 (2.32)

This is an estimation which might create the need for an even smaller time step. The time incrementa- tion is fixed with can be initially stable but chances in the analysis creates the need for a smaller time step. If Δtstable can as said be initially stable with the choice of time step Δt, but if properties change along the way the stable time can be decreased to be lower than the time step. Therefore, it is said to be con- ditionally stable. The implicit Newmark analysis is however unconditionally stable as said in section

2.3.2, where the time step has no upper limit of the size of the time step. A reasonable time step is still required for accuracy reasons, but it has no intendencies to fail for having to small increment in time.

2.3.4 Damping

Internal damping is a hindrance and brake for a structure to sway (Löfquist, 2016). Structures has mode- shapes in which the structure resonance, the different mode shapes describe different periodic variations in which the structure can sway and are usually created with periodically unfavourable excitation forces.

For example, half a wavelength can be called mode shape one, with the acting frequency called the eigenvalue.

In big building a dampener can be added to stop the building for resonating, where instead of the exci- tation force setting the structure in motion some energy sets a large pendulum in motion instead which creates a uneven distribution of force stopping the building for resonating (Engineering, 2016). Structures has some internal damping which is the source of its energy loss mechanics.

One way of describing internal damping is with Rayleigh damping which is based on one or more shape- modes and a guessed damping ratio. Rayleigh uses a α & β which describes a relation between the angular frequencies and the damping ratio.

α described in equation (2.31) corresponds to the mass of the structure and is unfavorable with higher damping ratio which can be shown in figure (8)

𝜉_𝑥 = 𝛼

2𝜔_𝑥 (2.31)

β described in equation (2.32) corresponds to the stiffness of the structure and is favorable with higher damping ratio which can be shown in figure (8)

𝜉_𝑥 =𝛽𝜔_𝑥

2 (2.32)

Combination of equation (2.31) and (2.32) gives the following expression described in equation (2.33) 𝜉 = 𝛼

2𝜔+𝛽𝜔

2 (2.33)

By analyzing two damping ratios with two shape- modes with corresponding angular frequencies the following matrix can be created described in equation (2.34)

{𝜉₁ 𝜉₂} =1

2[1 𝜔⁄ ₁ 1 𝜔⁄ ₁ 1 𝜔⁄ ₂ 1 𝜔⁄ ₂] {𝛼

𝛽} (2.34)

Taking the determinant and solve for α & β in equation (2.34) gives an expression of internal dampening in equation (2.35). The formulation describes the energy losing properties the structure has between the two shape-modes, which creates a system of how the internal damping locks like in the material.

{𝛼

𝛽} = 2 𝜔₁𝜔₂

𝜔₂²− 𝜔₁²[ 𝜔₂ −𝜔₁

−1 𝜔⁄ ₂ 1 𝜔⁄ ] {₁ 𝜉₁

𝜉₂} (2.35)

Figure 8 Relationship between Rayleigh coefficients (Dlubal Software GmbH 2020, 2020)

3 Method for investigation

3.1 Experimental project

3.1.1 Introduction to the experiment project

In 2003 Johan Magnusson & Mikael Hallgren issued by the Swedish defense research agency performed a test of 16 high-strength concrete beams of which 12 where induced to an air blast test (Magnusson &

Hallgren, 2003). The beams were divided into sets of three which involved different ways of reinforcing and casting. First set was normal reinforced, second set where normal reinforced with added steel fibres, third set had high-strength concrete casted in the bottom layer and medium-strength concrete in the top layer with added steel fibres and the last set was an higher strength concrete in the top layer then the previous set with added steel fibres.

Only the first set will be examined in this thesis since this is the easiest to replicate in a software. The first beam tested resulted in a failed test, which caused that data to be corrupt. From the other two named D2 and D3 the D3 test yielded the best results and will be used for the comparison. To prevent the issue of only getting an acceptable comparison between RFEM and the experiment for overworking the model the same model will be used with the inputs from D2, to see if both experiments managed to get a compatible result.

3.1.2 Test objects

All the beams had the same dimensions of 1720 mm x290 mm x160 mm (Length, width, thickness). The concrete was developed at the Norwegian University of Science and Technology which had an aggregate of crushed mylonite which categorised it as a high strength concrete. In table (1) the characteristic com- pressive strength can be shown for the concrete used for beam D2 and D3.

Table 1 Concrete properties Con-

crete grade

Age

(Days) fc,cube150

MPa fc,cube100

MPa fc,cyl100

MPa fspl,cube150

MPa fspl,cyl100

MPA ρc

kg/m³

140 73 107,2 +-

5,4(3) 106,6 +-

3,5(3) - - - 2428

140 97 115,2 +-

1,4(3) 115,2 +-

0,6(3) - - - 2458

140 103 - - - - 8,1 +-

0,4(3) 2487

fc,cube150 = characteristic compressive strength from test of 150mm concrete cube

fc,cube100 = characteristic compressive strength from test of 100mm concrete cube

fc,cyl100 = characteristic compressive strength from test of Ø100x200 mm concrete cylinder

fspl,cube150 = The splitting tensile strength determined on 150 mm concrete cube

fspl,cyl100 = The splitting tensile strength determined on Ø100x200 mm concrete cylinder

ρc = Concrete density

The beams were reinforced with 6Ø12 mm tensile and 2Ø10 mm compressive rebars of steel grade B500BT along with stirrups of 6Ø s200-N of the same quality. The bars where tested and had a yield strength fsy of around 550 MPa, see table (2).

Table 2 Steel properties Ønom

mm Øreal

mm fsy

MPa fsu

MPa Es

GPa 𝜀sy

‰

𝜀su

%

12 11,65 +-

0,05 544 +- 18 656 +-12 210 +-8 2,7 +-0,1 19,2 +-3,2 ønom = nominal steel bar diameter

øreal = measured steel bar diameter fsy = steel yield strength

fsu = steel ultimate steel strength Es = Young’s modulus of steel 𝜀sy = the yield strain

𝜀su = ultimate steel strain

Only the tensile reinforcement in the bottom layer will be included in the model. When designing a concrete beam, the structure needs to be resistant to moment and shear. The moment will be at its max in the middle of the beam. After a certain period of load the concrete will crack because of the low- tension capacity of concrete in the middle of the beam, and the tensile reinforcement will take over.

Then the steel will elongate as the beam deforms until it reaches the steels ultimate strength and start to be close to flexural failure. The shear forces acting will be resistant because of the shear reinforcement (stirrups) which will cause diagonal cracks largest near the supports, by placing vertical reinforcement bars the steel can intersect with these cracks and therefore make it more resistant to shear. The compressive steel will affect the structure but since they are in compression and has a small overall area, they are therefore deemed to not be necessary. Since only the middle-span is investigated only the tensional rein- forcement will be added to the model, as will be explained in chapter 3.1.3. The placing of the bars can be seen in figure (9)

Figure 9 Distribution of reinforcement (Magnusson & Hallgren, 2003) 3.1.3 Blast layout

The beams were anchored with two bolts and nuts at each support 1,5m apart vertically with an explosive charge 10m away inside a shock tube at the beams centre, see figure (10). The charge was a plastic

explosive which consisted of 86% PETN and 14% fuel oil. The charge weight was 2,0 kg for the beam D1, 2,5 kg for the beam D2 and 3,0 kg for the beam D3. In Appendix E there is a picture of the used shock tube.

Figure 10 Shock tube setup (Magnusson & Hallgren, 2003)

To measure the pressure there were two pressure gauges at each side of the beam. As discussed in sec- tion 2.1.2 the blast will reflect at the surface from the inside of the shock tube, see figure (10). This will cause a non-linear outlook of the pressure- time profile. To measure the acceleration three accelerome- ters where one was at the midspan and the other two at each quarter span of the beam length. And to measure the deflection two deflection gauges where used in the midspan of the beam, see figure (11) for the placing of the gauges.

Figure 11 placement of gauges (Magnusson & Hallgren, 2003) 3.1.4 Results from tests

The results from the test will be compared to the results from RFEM. This result will be the comparison to reach the goal if RFEM can produce a similar outcome as the experiment. The result consists of time diagrams which include Reflective pressure, acceleration and deflection, see Appendix B for the time diagrams. To compare the results the analysis boundary was set to the period up to around where the

beam reached its max displacement point, after that the potential damage caused by an impulse load would have occurred.

3.2 RFEM Simulation introduction

The three critical factors for analysing the modelling is:

• Accurate results

• Computing time

• No failed simulation

When using a finite element software like RFEM the simulation times can be relatively long depending on the type and size of the calculation. The software needs a detailed model to accurately represent the result and to not cause a failed simulation. A failed simulation can be caused by the pathway of the calculation to be destroyed, which means that the current position of the calculation has nowhere to go and therefore stops the simulation. Using the results from the experimental project which include the acceleration-time diagram and the deflection- time diagram the models’ improvements can be analysed.

These where used to analyse early stage test to see if the model responded correctly.

The methodology for reaching the result was a practical approach. Different kinds of modelling tech- niques and settings were tried. In this following chapter all different parts of creating the model are shown.

Some test gave results which lead to the realisation that something was wrong, by doing this the best method for all steps could be identified. Other methods gave frequent error and was overall hard to work with which also narrowed down the method for the final simulation. The simulation was said to be completed when all steps had been acknowledged which then produced a good result. This process was made using the results from experiment D3. Then the same model was used in experiment D2 so the results can be analysed from the accuracy of two different load cases.

The steps taken was the following:

• Variation of modelling

• Variation of Material model

• Material description

• Support boundary conditions

• Load

• FE-mesh settings

• DYNAM-PRO

In Appendix A there is a step-by-step guide how the modelling was executed.

3.2.1 Variation of modelling

The concrete structure can either consist of a surface (shell) elements or solid elements (Dlubal Software, 2020). A surface is bounded by n numbers of boundary lines connected by n numbers of nodes. Resulting from material and thickness properties gives the surface a stiffness and mass. Solids are 3D elements created by “null” surfaces but have the same properties as a surface. The difference is the dimensions of the FE mesh and where solids has the properties to connect to other solids. The solid also needs to calculate the entire structure for the stiffness and mass.

Reinforced concrete design for solids is not implemented in RFEM according to the RFEM software manual (Dlubal Software, 2020). Since there is not possible way to insert reinforcement automatically the steel bars need to be placed manually. When designing concrete structures in RFEM the concrete is added and using a reinforcement add-on the appropriate reinforcement can be determined, here the reinforcement is already chosen. When placing the reinforcement, it is important to get an integration (connection) between the concrete and steel, otherwise the materials will be calculated individually. Dif- ferent methods were tried but the most successful method was using a surface with steel beam elements which were placed inside the surface at the correct position in the yz-plane that was then extruded down

to the correct placement in the x- direction. Since surface elements are 2D element only the result in that specific strip is saved and can be analysed. But the surface has a cross-section and stiffness as a 3D object, it’s just not visible in the software. Any integrated objects will be included in the surface cross- section which makes the integration acceptable. Beam elements are predefined geometries taken from the RFEM geometry library, see figure (12) for the final model.

The process of modelling was to place four corner nodes, connected with lines and placing a concrete surface bounded by the lines. Then placing nodes and reinforcement bars in the surface so that they will be integrated with the surface, see figure (12).

Figure 12 Model structure with (gray surface) concrete surface, (blue bars) extruded reinforcement, (green pyramids) nodal supports, (red dots) nodes.

The varying in parameters for the structure consisted of getting a correct integration. For this a test was run with and without reinforcement were a static and dynamic test was run to compare the differences in a DLC using the add-on module DYNAM-PRO and a SLC using original RFEM calculation.

3.2.2 Variation of material model

The material model describes the stress- strain relationship for the two isotropic (same properties in all directions) materials. In the software there is 3 different material models to choose from. The default setting is the Isotropic Linear Elastic model. The linear elastic model had a linear increase on stress with

plasticization of the steel and no upper limits of compression and tension capacity. The second model is the isotropic plastic model. This has an upper limit of compression and tension capacity and creates a plastic deformation which creates the response of the deformation not returning to its original state. The last one is a non-linear elastic model. This has upper limits for compression and tension capacity, but it does not insure a plastic response.

For the steel beams the Isotropic Plastic 1D material model was used and Isotropic Plastic 2D/3D for the concrete. This will ensure that there is a permanent deformation. When the yielding point is reached the Young’s modulus will convert to a strain hardening modulus. At this point RFEM will recalculate with the new modulus.

As shown in figure (13) there is a shaded dotted arrow, this is the way back for the strain when stress is decreased. The yield strain shown in table (2) was around 550 MPa with a Strain hardening modulus of 0,21 MPa.

Figure 13 Plastic material model of steel

Since concrete and steel has different properties when it comes to the difference in compression and tension a separated strain hypothesis is set for the concrete. For the concrete the Drucker-Prager is set, this will give two different properties dependent if a part of the beam is in compression or tension, see figure (14). The compression strength shown in table (1) was set to 115 MPa with a tension strength of 5,2 MPa, where the strain hardening modulus was 0,45 MPa.