LOCKING IN FINITE ELEMENT ANALYSIS prepared by S

LOCKING IN FINITE ELEMENT ANALYSIS prepared by S. S. Nantha Kumar M. VIGNESWARI N. SHANMUGA PRIYA Thiagarajar College of Engg.

INTRODUCTION All physical problems can be expressed mathematically in terms of ordinary or partial differential equations. All differential equations cannot be solved analytically. Many DE’s do not possess closed form solutions. Finite Element Method is an effective numerical method for solving differential equations. The differential equation in structural mechanics is the Differential equation of equilibrium along with displacement and force boundary conditions.

How Does Finite Element Analysis Work? The mathematical representation of the physical system in FEM is by using nodes and elements. The process of creating nodes and elements is said as discretization or meshing. The number of degrees of freedom in each node depends on the element type. As the FEM model has finite DOF as against infinite DOF in actual structure, FEM solutions are approximate and FEM model is stiffer

AIM OF THE PROJECT : In spite of providing approximate solutions to problems which cannot be solved analytically, there are some situations for which convergence is poor in FEM. The aim of our project is to bring out such an error called “locking” and how it is being overcome by the commercial FEA codes. The FEA softwares used are ANSYS and ABAQUS.

LOCKING : Locking means the effect of a reduced rate of convergence in dependence of a critical parameter. eg. , slenderness of the plate in case of transverse shear locking. Because of locking , wrong displacements , false stresses, spurious natural frequencies are encountered. Locking may be �Volumetric locking �Membrane locking �Shear locking

Contd… Shear locking occurs on use of fully integrated linear elements that are subjected to pure bending. The formulation of these elements “promotes” shear strains that do not really exist. Volumetric locking occurs in fully integrated elements when the material behavior is nearly or fully incompressible( poisson’s ratio approaches or equals 0. 5) Membrane locking occurs only in Shell and Curved beam elements.

Shear locking: Timoshenko beam theory works for all l/d ratio whereas Bernoullie beam theory works for only large values of l/d ratio. In thick beams , shear deformation becomes predominant and Bernoulli’s theory wont give accurate result. Shear locking occurs in elements which account for shear deformation , when these elements are used in modeling bending predominant beams.

Modelling • A Cantilever beam subjected to loading at the free end is analysed with the following properties : Young’s Modulus = 200 GPa Poisson’s Ratio = 0. 3 • Ansys Element SOLID 45 and SOLID 185 is considered.

MESHING : The beam is meshed such that there are 4 elements in height , 10 elements along the length and 1 element along the width.

RESULTS USING SOLID 45: L/d ratio = 25 L/d ratio = 5

The analytical solution for the tip deflection in case of a cantilever problem is calculated from • The Answer obtained analytically for a cantilever beam of L/D ratio 25 is 0. 312 cm The free end deflection obtained from SOLID 45 is 0. 308 cm whereas from SOLID 185 element the tip deflection is found to be 0. 1514 cm. Thus it is apparent that the results from SOLID 185 cannot be relied upon. The difference between the results of Solid 45 and Solid 185 is due to shear locking and it is not an approximation error.

ANALYTICAL DISCUSSION : • Consider Stiffness matrix, K , for the beam element and the expression is given by o. The bending stiffness can be exactly evaluated by one point gauss quadrature, on the other hand the shear stiffness matrix (shear strain energy) contains second order terms and two point gauss quadrature is required to exactly integrate it.

GAUSS POINT INTEGRATION: • The bending stiffness term is given as where E is the young’s modulus The Shear Stiffness term is given as where G is the rigidity modulus and α is the shear correction factor.

One point quadrature : One point quadrature “under integrates” the shear stiffness (energy) , thus resulting in closer solution. The stiffness matrix for the beam is given by, Applying the boundary conditions and defining ψ = and β= , we get the equilibrium equation for the beam as

Contd. . Solving we get, W= and θ= v In the limiting condition i. e. , tending to thin beam limit β >> ψ, the equation becomes, v W=l and expression for θ remains unchanged. Thus the beam deformation is solely due to bending as given by the equation and this is correct.

Two point quadrature: Two point quadrature integrates exactly the shear stiffness causing an overly stiff element. This leads to erroneous results as they contain only the coefficient corresponding to shear deformation. The stiffness matrix can be calculated as, The equilibrium equation for the beam after imposing boundary conditions is given by, • Solving for w , we get, W = P +

Contd. . Solving for θ , we get θ = In the thin beam limit β >> ψ , the expressions become, W = θ Ø = It can be observed that these expressions lead to erroneous results as they contain only the co-efficient β corresponding to shear deformation.

One point integration leads to a more flexible mesh and in this case only one constraint is imposed on the element. In two point integration , the stiffness is completely dominated by the shear stiffness and the mesh gets ‘locked’. Fully integrated lower order elements exhibit “overstiffness” in bending problems. This formulation includes shear strains in bending which do not physically exist, called parasitic shear.

Tip deflection values in cm s. no L/D ratio 1. 25 2. ANALYTICAL SOLUTION SOLID 45 SOLID 185 (full integration) % Error SOLID 185 (reduced integratio n) 0. 3128 0. 308 0. 1514 50. 8% 0. 3328 20 0. 16028 0. 1578 0. 0953 39. 6% 0. 1705 3. 15 0. 0677 0. 0667 0. 0489 26. 68% 0. 0720 4. 10 0. 02014 0. 0196 0. 0145 26. 28% 0. 0214 5. 5 0. 002569 0. 002523 0. 00235 6. 8% 0. 002757

Interpretation of the results : It is clear from the table, as the beam thickness increases (as shear becomes predominant) the error is minimized. The error in the tip deflection values of bending predominant beams is due to Shear locking.

How ansys Overcomes locking : In FEA codes , some element technologies have been incorporated Extra shape function (ESF) Reduced Integration SOLID 45 has extra shape function inbuilt in it, thus it gives exact solution irrespective of the L/d ratios. SOLID 185 has reduced integration option which when enabled produces closer results.

Extra shape function Here u 1, u 2 , u 3 are node less variables

ABAQUS MODELLING The 8 noded brick element in abaqus is the continuum 3 -dimensional element – C 3 D 8 This element is studied further for the locking problem. Cantilever beam is modeled for various L/D ratios & the results from both ansys and abaqus packages are tabulated.

RESULTS USING C 3 D 8 : L/d Ratio = 25 L/d Ratio = 5

DEFLECTION RESULTS: L/D ratio ANALYTICAL SOLUTION Solid 45 (mm) C 3 D 8(mm C 3 D 8 R ) (mm) C 3 D 8 I(m m) 25 0. 3128 3. 08 0. 9189 3. 32 3. 054 20 0. 16028 1. 578 0. 6312 1. 701 1. 565 15 0. 0677 0. 667 0. 362 0. 7183 0. 6612 10 0. 02014 0. 1967 0. 1454 0. 2137 0. 1968 5 0. 002569 0. 0252 0. 0235 0. 02752 0. 02523

ABAQUS SOLUTION Abaqus has got the following element technologies to overcome locking: Reduced integration(C 3 D 8 R) incompatible modes(C 3 D 8 I) C 3 D 8 element has default setting as reduced integration. Enabling incompatible modes option, the answer obtained is still more accurate.

INCOMPATIBLE MODES: 8 -node brick elements are frequently used with incompatible modes to improve the bending behavior in case of pure displacement problems. These elements additional degrees of freedom to improve the accuracy of the element. The internal nodes are condensed out of the element equations, so the final element matrices have their original DOF s but the element retains its accuracy.

REDUCED INTEGRATION ELEMENTS: Fully integrated linear hexahedral elements(C 3 D 8)use two integration points in each co-ordinate system. The reduced integration elements use fewer integration points and will decrease shear locking because some terms in the gauss integration are eliminated.

A numerical approach to locking: The tendency of an element to over-stiff behavior is measured with help of a constrain ratio r. this number is defined as r = Ndof / Ncon For an element with i nodes, the constraint ratio is given by

Contd. , It can be seen that the constrain ratio r gets smaller as more quadrature points are used. As the number of quadrature points decrease, the value of r increases & hence element becomes flexible. This explains reduced integration as remedy to locking issue.

MEMBRANE LOCKING: Occurs only in curved beam and shell elements. In curved elements , both flexural & membrane deformation takes place Such elements perform poorly in case where inextensional bending becomes predominant i. e. Membrane strain becomes vanishingly small when compared to bending strain.

Theoretical Proof • The Strain – Displacement Relations are Where ε is the membrane strain and χ is the bending strain • The circumferential and radial displacements are given by • The strain field interpolations

INEXTENSIONAL BENDING If a thin and deep arch such that L/T>>1 and R/H is small is modeled by a curved shell element , the physical response is inextensional bending such that the membrane strains tends to vanish. It tends to the following constraints Each in turn implies the condition These are the spurious constraints. Each in turn implies the conditions

Contd. , A shell was modeled using ansys elements SHELL 93 & SHELL 63 having the following properties: young’s modulus = 30 Gpa Poisson’s ratio = 0 shell thickness = 3 cm radius of curvature = 300 cm length of the shell = 600 cm pressure applied = 6250 N/m 2

Shell Elements :

ANSYS MODEL : shell 63 angle 40 degrees shell 63 angle 160 degrees

STAAD PRO MODEL :

DEFLECTION AT EDGE NODE : S. NO ANGLE STAAD (DEGREES) PRO RESULT SHELL 63 -1. 3526 SHELL 93 % ERROR -0. 7837 41. 8 1. 160 -1. 311 2. 120 -1. 809 -1. 7769 -1. 4312 19. 4 3. 80 -2. 4012 -2. 23 -2. 0979 5. 9 4. 40 -1. 385 -1. 3582 -1. 3509 0. 0053

Interpretation : As the angle becomes larger (as the shell becomes deeper), the difference between the results of edge deflection obtained from shell 63 and shell 93 is greater. But Shell 63 results are almost closer to the Staad. Pro results. This shows that for deeps shells Shell 63 should be used instead of Shell 93 has spurious membrane strains introduced in it as the shell becomes deep and Shell 63 donot encounter such an issue.

REMEDY : Reduced integration techniques are not always successful in overcoming locking behavior and result in overstiffness. Shell 63 is a linear element whereas Shell 93 is a curved element with midside nodes. Faceted modeling(planar elements)is currently the most reliable and successful method. Curved modeling(curved elements)has only limited success.

CONCLUSION In Ansys , Element technologies like Reduced Integration and Extra Shape Function has been adopted to overcome Shear locking. Similarly in Abaqus , Reduced Integration and Incompatible Modes can be used to overcome the issue. For Membrane Locking , the remedy is the judicious selection of appropriate elements i. e. Linear elements instead of curved elements in deep arches and shells.

References : Finite Element Analysis by C. S. Krishnamoorthy , IIT , Madras. Strength of Materials- Elementary theory and problems by Stephen Timoshenko Eric Quili Sun, “ Shear locking and hourglassing in MSC Nastran, Abaqus and Ansys” www. cmmacs. ernet. in/cmmacs/pdf/ www. imechanica. org www. mscsoftware. com/events/vpd 2006/na/presentations/. . . /27. p df www. eng-tips. org www. springerlink. com www. openpdf. com www. sciencedirect. com

THANK YOU