FEM analysis of buckling

FEM analysis of buckling

Jerzy Pamin e-mail: JPamin@L5.pk.edu.pl

With thanks to:

M. Słoński, M. Radwańska, A. Wosatko ANSYS, Inc. http://www.ansys.com ROBOT http://www.autodesk.com

Lecture contents

Buckling phenomenon

Algorithm of FEM buckling analysis

FEM in buckling simulations

Nonlinear anaysis of RC shells

Buckling phenomenon [1,2]

Assumptions of linear buckling analysis:

I one-parameter loading, varying proportionally to load parameter λ P = λP^∗

I loading is conservative, i.e. does not change direction during structure deformation

I structure (bar, panel, shell) is ideal, with no geometrical, material or load imperfection which would disturb ideal pre-buckling state

Buckling phenomenon cont’d

Buckling occurs when increasing load reaches critical value Pcr = λcrP^∗, where P^∗denotes so-called configurational load for which λ = 1.

Characteristic feature of buckling as one of loss of stability phenomena is the significant change of deformation mode of structural system which experiences compressive stresses as a whole or in some part.

Source: E. Ramm, Buckling of Shells, Springer-Verlag, Berlin 1982

Examples of buckling phenomenon

Static criterion of the buckling as one of loss of stability types consists in examination of close pre- and post-buckling states. The phenomenon is presented for

I simply supported bar,

I deep cantilever beam,

I unidirectionally compressed panel, simply supported along circumference,

I cyllindrical shell under normal pressure, clamped along lower edge.

Buckling of a bar

Before buckling

The bar:

I has straight axis,

I is only compressed (is not bent).

After buckling

The bar:

I has curved axis,

I exhibits compression and bending.

Buckling of deep cantilever beam

Before buckling

I The beam is bent in plane by vertical force applied at the free end

X

Y

Beam displacements in pre-buckling state

After buckling

I lateral buckling (warp, twist) occurs due to coupled bending and torsional deformation

Buckling of deep cantilever beam cont’d

Z X

Buckling modes (forms)

Buckling of panel compressed in one-direction

Before buckling

Ideal membrane state:

I Panel with ideal medium plane,

I Constant compressive loading along one direction in the medium plane.

After buckling

Bending occurs:

I non-zero displacements perpendicular to medium plane,

I non-zero curvatures and bending moments.

Buckling of compressed panel (ANSYS, [3])

First and second buckling mode

Third and fourth buckling mode

Buckling of cylindrical shell under external radial pressure

Before buckling

In the shell:

I axisymmetric conditions,

I in large part of the long shell pure membrane state,

I bending in vicinity of clamped edge (flexure) state.

After buckling:

Significant disturbance of axisymmetry:

I waves along circumference,

I number of half-waves is different for subsequent critical multipliers of the applied loading.

Buckling of cylindrical shell cont’d (ANSYS, [3])

Subsequent buckling modes

General buckling analysis [1,2]

Energetic buckling criterion

Energetic buckling criterion consists in the analysis of an increment of potential energy Π during transition from pre- to post-buckling state.

One considers two adjacent equilibrium states:

I pre-buckling state I for which

δΠ^{(I )}= 0

I post-buckling state II for which

δΠ^{(II )}= δΠ^{(I )}+ δ∆Π = 0

I energetic criterion of critical state δ∆Π = 0.

Algorithm of FEM buckling analysis

Matrix equation describing the loss of stability viz. buckling [K₀+ λK_σ(s^∗)]v = 0

or

{K0+ λ[Kσ(s^∗) + Ku1(g^∗)]}v = 0

where:

I linear stiffness matrix of the system K₀

I initial stress matrix Kσ(s^∗) and initial displacement matrix Ku1(g^∗)

I critical loading multiplier to be determined λcr

I respective post-buckling form represented by eigenvector v = ∆d

Pre-buckling statics

Stage I of algorithm:

1. Compute the global stiffnesss matrix K₀

2. Compute nodal forces representing initial loading configuration P^∗, i.e. for loading multiplier λ = 1 (one-parametr loading assumed P = λP^∗)

3. Take boundary conditions into account

4. Solve equation set K0· d^∗= P^∗, to obtain nodal displacements in pre-buckling state: d^∗= K⁻¹₀ · P^∗

5. From displacements of the system d^∗ extract element dofs d^e∗and compute in each element:

I displacement gradients g^e∗and

I generalized stresses s^e∗.

Buckling analysis

Stage II of algorithm:

1. Generate:

- initial stress (geometrical) matrices for each element K^e_σ(s^∗e) and the whole structure K_σ(s^∗)

- optionally initial displacement matrix K_u1(g^∗)

2. Formulate non-standard (generalized) eigenproblem representing linearized buckling problem: [K0+ λ(Kσ+ Ku1)]v = 0

or initial buckling problem: [K₀+ λK_σ]v = 0

3. Solve the eigenproblem to determine the pairs (λ1, v1), . . ., (λN, vN) where:

I N – number of dofs

I λi – eigenvalue - critical loading multiplier

I vi= ∆di – eigenvector - post-buckling deformation mode

Buckling of ideal panel/plate [4] – input data

I dimensions: Lx = Ly = 1.16 m, h = 0.012 m

I material data: E = 2.05 · 10⁸kN/m², ν = 0.3

I configurational loading along perimeter which represents in-plane bending: |px ,max ,min^∗ | = 1.0 kN/m

I two options of boundary conditions along circumference:

a) simply supported (hinged, right) b) fully supported (clamped, left)

Panel buckling

Assumptions:

I ideally flat medium plane,

I loading acts exactly in the medium plane,

I the one-parameter loading is governed by λ factor.

Buckling analysis of ideal panel under pure in-plane bending

Loading which causes pure in-plane bending prior to buckling Computations:

I numerical (FEM packages ANKA and ROBOT) – approximate solution

I analytical – exact solution

Panel buckling for in-plane bending

Computation of critical load:

Loading and deformation in pre-buckling state

Analytical solution for:

I simply supported panel p^ib,analitcr = ^{25.6·π2·Dm}

L²_x = 6077 kN/m

I clamped panel p^ib,analitcr =^{39.0·π2·Dm}

L²_x = 9259 kN/m Numerical solution (ANKA, mesh 8 × 8 ES) for:

I simply supported panel p^ib,MEScr = 6028 kN/m

I clamped panel p^ib,MEScr = 11304 kN/m

Numerical solution (ROBOT, mesh 12 × 12 ES) for:

I simply supported panel p^be,MEScr = 6241 kN/m

be,MES

In-plane bending in pre-buckling state

Distribution of membrane force n_x

for the simply supported (left) and clamped (right) panel

In-plane bending, buckling modes

First two buckling modes for simply supported panel (ROBOT)

In-plane bending, buckling modes

First two buckling modes for clamped panel (ROBOT)

Buckling of I-beam – input data

I dimensions: L_x = L_y = 1.16 m, h_s = 0.012 m, h_p= 0.018 m

I material data: E = 2.05 · 10⁸kN/m², ν = 0.3

I configurational loading along beam sections:

|p_{x ,min,max}^∗ | = 1.0 kN/m

I two options of buckling analysis (ROBOT):

option 1: local buckling of the web

option 2: buckling of beam segment (web+flanges)

Option 1: web buckling

Local web buckling:

I for isolated web panel, in reality connected to flanges and ribs, different boundary conditions can be imposed along the connection lines

I in limiting cases one can assume:

a) hinged support along whole circumference b) clamped support along whole circumference

I hence the actual situation is inbetween

I former computations can be used to consider the web buckling

Option 2: I-beam buckling

Buckling analysis for the beam:

I discrete model in ROBOT for I-beam composed of web (12 × 12 elements) and two flanges (4 × 12), loaded by bending in the plane of the web

I numerical results (ROBOT):

I pcr^be,MES= 9068 kN/m

I comparison of critical forces computed with FEM (ROBOT):

I for isolated web:

- simply supported (ss) - clamped (cl)

I whole beam segment (be) p^ss,MES < p^be,MES < p^{cl ,MES}

6241 kN/m < 9068 kN/m < 11666 kN/m

Bending of I-beam in pre-buckling state

Membrane force distribution nx in pre-buckling state

Buckling modes for I-beam

Two buckling modes for I-beam segment subjected to bending (ROBOT)

Geometrically and physically nonlinear analysis [5]

K Q Q P R Q

t Q

T T

tN

( )⋅ − ⋅ = ( , )

+ =

∗ +

∆ ∆

∆ ∆ ∆

λ λ

λ τ

1

∆Q, Q

A^-1

∆q, q

B

∆εεεε, εεεε

P L E S

K^e, R=P-E

A

k^e, f^e

∆σσσσ, σ σ σ σ, D^e-p

Scheme of computation strategy at the levels of:

I structure I finite element I layer I point Effects considered:

I stress evolution in cross-section I elastic-cracking concrete I elastic-plastic reinforcement I large displacements and their

gradients Aims:

I computation of displacement evolution

I determination of damage mechanism

I estimation of load-carrying

RC shell model

I Degenerated 8-noded shell element (Mindlin-Reissner theory)

I Layered RC shell model (5 concrete layers, 4 steel layers representing two reinforcement grids)

I Continuum elastic-cracking model for concrete layer (concrete softening, reduction of shear stiffness)

I Elastic-plastic model for steel layer

Numerical analysis of cooling tower shell [5]

Numerical analysis of cooling tower shell

Diagrams λ − wK obtained using two FEM packages using force or displacement control for loading g + λ(w + s)

Cooling tower loads:

I self-weight g

I wind w

I internal suction s

I temperature variations

I subsidence

Numerical analysis of RC shell

Results for shell with technological opening

Numerical analysis of RC shell

Directions of principal stresses in external layer Smeared cracks visualization

Damaged cooling tower shell [6]

Analized cases for load combination g + λ(w + s)

I designed shell

I built shell with zones of weak concrete (fcm=11 MPa)

I shell with two circumferential openings (25m and 14m in length)

I repaired and strengthened shell (5cm reinforced shotcrete in height

Linear buckling of cooling tower shell (DIANA)

Loading Designed Constructed Damaged Repaired

λg 25.52 22.22 18.59 19.11

λ(w + s) 13.15 14.72 6.24 20.71

λ(g + w + s) 11.29 11.39 5.49 20.84

Critical load multipliers λ1

Nonlinear analysis results

The construction error did not have a significant influence on the short-term load carrying capacity of the cooling tower, but it affected its durability due to local concrete overload and reinforcement corrosion.

Cracking zone prediction (DIANA)

Cracking zones in inner and outer concrete layer for λ = 3.2 (DIANA)

Model of shell with holes (DIANA)

For smeared cracking computations diverged for λ ≈ 1.0 - it is necessary

References

[1] M. Radwańska. Ustroje powierzchniowe, podstawy teoretyczne oraz rozwiązania analityczne i numeryczne. Wydawnictwo PK, Kraków, 2009.

[2] Z. Waszczyszyn, C. Cichoń, M. Radwańska. Stability of Structures by Finite Elements Methods. Elsevier, 1994.

[3] M. Bera. Analiza utraty stateczności wybranych tarcz i powłok sprężystych metodą elementów skończonych. Praca dyplomowa, Politechnika Krakowska, Kraków, 2006.

[4] M. Radwańska, E. Pabisek. Zastosowanie systemu metody elementów skończonych ANKA do analizy statyki i wyboczenia ustrojów powierzchniowych. Pomoc dydaktyczna PK, Kraków 1996.

[5] Z. Waszczyszyn, E. Pabisek, J. Pamin, M. Radwańska. Nonlinear analysis of a RC cooling tower with geometrical imperfections and a technological cut-out. Engineering Structures, 2, 480-489, 2000.

[6] A. Moroński. Analiza zarysowania i utraty stateczności uszkodzonej powłoki żelbetowej chłodni kominowej. Praca dyplomowa, Politechnika Krakowska, Kraków, 1996.

[7] A. Moroński, J. Pamin, M. Płachecki, Z. Waszczyszyn. Fracture and loss of stability of a partly-damaged cooling tower shell. Proc. 2nd Int. DIANA Conf. on Finite Elements in Engineering and Science, Eds M.A.N. Hendriks et al, 107-110, Kluwer Academic Publishers, Dordrecht, 1997.