# 3D Finite Element Stress Analysis

### 3D Finite Element Stress Analysis

### Computational Continuum Mechanics

Born in Brunswick, Germany, Gauss is regarded as a great mathematician for work on number theory, geometry, probability theory, geodesy, planetary astronomy, the theory of functions, and potential theory including electromagnetism.

Young Gauss excelled at classics, languages, and mathematics with his ability to do complex calculations. Astounding his teachers he was introduced to the Duke of Brunswick in 1791, who became his patron. While at school, in 1792, Gauss made his first significant finding, the construction of a seventeen sided polygon using only a compass and ruler. The proof to this work relied on analysis of the factorisation of polynomial equations which later led to the ideas of Galois Theory. The Duke's financial support enabled Gauss to study mathematics at Gottingen University from 1795 to 1798. Here Gauss produced a proof of the fundamental theorem of algebra. His thesis of 1797 demonstrated that every polynomial equation with real or complex coefficients has as many roots as its degree. During Gauss's life he gave three further proofs of this thesis.

Possibly, Gauss's most remarkable work came in 1801 when he produced two publications. The first was a systematic textbook on algebra number theory, Disquisitiones Arithmeticae. The second was his rediscovery of the asteroid Ceres proving his ability to calculate the asteroid's orbit accurately enough to determine its reappearance. Subsequently, Gauss worked on astronomy and in 1809 published a book on planetary motion and the computation of orbits.

Gauss married in 1805. He had three children. After the birth of his second son his wife died followed by the death of his child. Shortly after, Gauss remarried.

Gauss undertook a project surveying the state of Hanover between 1818 and 1832, during which he invented a heliotrope to measure points at great distances and improve the accuracy of observations. He also formulated the concept of the curvature of a surface. Gauss established an interest in terrestrial magnetism participating in 1830 in a survey of Earth's magnetic field, which led him to invent the magnetometer. With the help of physicist Wilhelm Weber he developed important mathematical research in electromagnetism and gravitation now known as potential theory. In 1823 Gauss was awarded the Danish Academy of Sciences' prize for his investigations into angle-preserving maps.

Important pieces of work by Gauss were never published, often being left for others to recreate and claim credit. After his death in Gottingen the wealth of Gauss's unpublished papers extended his influence through the remainder of the century.

### 1b. George Green (1793 - 1841)1,2

Green was an English born mathematician, whose work heralded the start of modern mathematical physics in Britain, pioneering the application of mathematics to physical problems. He was the only son of a Nottinghamshire baker. Showing great enthusiasm and mathematical talent he was sent to the town's leading school in 1801. Green learnt everything the school could offer and was, subsequently, set to work in his father's mill. Subscribing to the Nottingham Library in 1823 by 1828 Green produced his first and most important work, ‘An Essay on the Application of Mathematical Analysis to the Theories of Electricity and Magnetism', which introduced the concept of potential and what is now known as Green's theorem. Published by private subscription it attracted little interest. However, it gained the attention of Sir Edward Bromhead who recognised the work's originality.

Green's life was not conducive to studying. His first child's birth in 1824 was followed by his mother's death in 1825. In 1827 a second child was born followed, in 1829, by his father's demise. Green had a further five children. Green returned to milling after the disappointing response to his paper, but in 1830 he met Bromhead who disseminate Green's work amongst Cambridge mathematicians.

Encouraged by Bromhead, Green wrote three subsequent papers. The first, ‘Mathematical Investigations Concerning the Laws of Equilibrium of Fluids Analogous to the Electric Fluid' was published in 1833; his second, ‘The Determination of Exterior and Interior Attractions of Ellipsoids of Variable Densities', in 1835; and his third, ‘Researches on the Vibrations of Pendulums in Fluid Media', in 1836.

Green enrolled at Caius College, Cambridge in 1833, taking the Mathematical Tripos in 1837, and emerging as fourth Wrangler. He was awarded the Perse Fellowship at Caius. Green produced six further papers on hydrodynamics, reflection and refraction of light, and reflection and refraction of sound. Each was published in the Transactions of the Cambridge Philosophical Society, of which he became a member.

Ill-health forced him to leave Cambridge and return to Nottingham where he died in May 1841.

### 2. Gauss-Green Derivation[3]

Starting with the three dimensional form of Gauss' divergence theorem: V div q dV= SqTn dS (2.1)

Where q is an arbitrary vector and n is the surface normal pointing outwards from the defined boundary [ 1], defined as:

q= qxqyqz n= nxnynz

Therefore, Gauss' divergence theorem can be written in longhead as:

V ∂qx∂x+ ∂qy∂y+∂qz∂z dV= SqxqyqznxnynzdS (2.2)

Equation (2.2) can be expressed in terms of stresses where q is equivalent to σ.

We also know that [σ]=D[ε] and [ε]=LNd[4].

Thus giving, σ= DLNd where L= ∂∂x000∂∂y000∂∂z∂∂y∂∂x00∂∂z∂∂y∂∂z0∂∂x

Therefore, div can be defined as LT= ∂∂x00∂∂y0∂∂z0∂∂y0∂∂x∂∂z000∂∂z0∂∂y∂∂x=div

Substituting into equation (1): Ω LTDLNd dΩ= ΓDLNdTn dΓ (2.3)

Now applying the Gauss-Green theorem and multiplying DLNd by a scalar cT

Ω LTcTDLNd dΩ= ΓcTDLNdTn dΓ (2.4)

Now using the identity divϕq= ϕdiv q+ ∇ϕTq in terms of stress and strain:

LTcTDLNd= cTLTDLNd+ LcTDLNd (2.5)

Substituting equation (2.5) into equation (2.4):

Ω cTLTDLNd dΩ+ΩLcTDLNd dΩ = ΓcTDLNdTn dΓ (2.6)

Finally, rearranging to give the Gauss-Green theorem

Ω cTLTDLNd dΩ = -ΩLcTDLNd dΩ +ΓcTσn dΓ (2.7)

The transpose has been neglected from this final form as σ= σT

### 3. The Weak Form[5]

The weak form is a permissive variation of the strong form governing partial differential equation of equilibrium, subject to chosen boundary conditions. A solution of the strong form can also satisfy the weak form, however, a solution to the weak form will not satisfy the strong form. Using the weak form the finite element equations three dimensional elasticity can be derived.

The strong form [equation (3.1)] constitutes a set of differential equations of equilibrium for all points within the subject domain. Therefore, satisfying the strong form asserts that all the differential equations at all points in the domain are also satisfied. However, this statement is not always true[6].

LTDLNd + Fb = 0 (3.1)[7]

The weak form [equation (3.3)] can be derived by firstly pre-multiplying the strong form by a weight function, C T and subsequently integrated over the volume [equation (3.2)]. This arbitrary vector, C T, is contrived as a result of the possibility that not all the points in the domain are in equilibrium and, therefore, their equations of equilibrium do not equal zero. The weak form, therefore, finds a solution which satisfies the strong form on average over the volume.

Ω cTLTDLNd dΩ +ΓcT Fb dΩ= cT0=0 (3.2)5

Integration by parts, for example by using the Gauss-Green theorem, is then performed and finally substituted back into equation (3.1), providing the weak form which holds for any three dimensional or two dimensional problems:

ΓcTσn dΓ-ΩLcTDLNd dΩ + ΩcTLTDLNd dΩ=0 (3.3)5

The weak form has been adopted for the finite element method because the order of differentiation of the function is lower than it is for the strong form. This facilitates the approximation process because a function which can be easily integrated once may have a second order derivative which tends to infinity, therefore, creating difficulties during integration. The weak form can often be a more realistic representation than the original strong form partial differential equation because a large number of differential equations do not admit a smooth solution, which is implied when using the strong form. Furthermore, the weak form applies without changes for continuous as well as discontinuous problems, implying the important conclusion that the finite element method also holds for continuous as well as discontinuous problems, whereas modifications are required of the strong form in the presence of discontinuities.

### 4. Material Elastic Anisotropy[8]

A material will exhibit different values of a property in different crystallographic directions when anisotropy occurs. There are five distinct classes of elastically anisotropy[9].

i. Triclinic materials exist with no rotations of symmetry. The [D] matrix is shown as

D11D12D13D14D15D16D21D22D23D24D25D26D31D32D33D34D35D36D41D42D43D44D45D46D51D52D53D54D55D56D61D62D63D64D65D66 (4.1)

ii. Monoclinic materials have symmetry in the z axis only. This results in a rotation of π. The [D] matrix is shown as

D11D12D1300D16D21D22D2300D26D31D32D3300D36000D44D450000D54D550D61D62D6300D66 (4.2)

iii. Orthotropic materials have two pairs of coordinate directions. Therefore, rotations of π and (nπ2+θy) will be symmetrical. The [D] matrix is shown as

D11D12D13000D21D22D23000D31D32D33000000D44000000D55000000D66 (4.3)

iv. Transversely Isotropic materials have a finite number of axes of rotation and elastic symmetries. Therefore, a rotation will not have any effect on the material's properties. The [D] matrix is shown as

D11D12D13000D21D22D23000D31D32D3300000012(D11- D12)000000D55000000D66 (4.4)

Samantha Osbaldiston11/12/09

v. After the Transversley Ispotropoic class of anisotropy the material will become fully Isotropic where it will possess an infinite number of axes of rotations.

Samantha Osbaldiston11/12/09

### 5. Modelling Near - Incompressible (and Fully Incompressible) Media

Incompressibility is one of most common constraints encountered in finite element analysis. An incompressible material or a material near incompressibility, such as rubber-like materials or materials which flow, can sustain large deformations without palpable volume changes. A material becomes incompressible as the corresponding Poisson's ratio approaches a value of 0.5. Near to this value stresses become strongly dependent on Poisson's ratio and can even double between the values of 0.48 and 0.5. Unless a problem is one of plane stress only, analysis of incompressible materials cannot be tackled because a Poisson's ratio of 0.5 will cause the denominator of the shear modulus and the [E] matrix to become zero. Consequently, as the denominators tend to zero they act as a penalty number warning of the approaching incompressibility constraint. Also, some displacement based finite element approximations are unstable when dealing with incompressible or near compressible materials. Ultimately, as additional nodes are brought to the mesh, it will bring extra nodal degrees of freedom. However, a mesh will ‘lock' because the addition of nodes and elements does not increase the number of ‘free' degrees of freedom. This problem can be overcome in a number of ways, which will be explained below.

Firstly, ‘locking' can be surmounted if the element bending stiffness matrix associated with the problem is singular. To achieve singularity it is necessary to perform underintegration.[10]

Secondly, a mixed method approach can be used based on the deformation and pressure present within the media, and by choosing suitable finite element approximations for each variable. Here, the mean stress is separated from the volume's total stress field and utilised as an independent variable. The mixed form (5.1)[11] can then be found by applying the independent stress and displacements to the weak form of the equilibrium equations.

ACCT-Vup=f1f2 (5.1)

Underintegration and the mixed form methods are the more favourable solutions to the problems encountered whilst a material's properties are incompressible or near incompressible. However, the use of transformations is an alternative method to impress relations amongst the structure's degree of freedoms. This method, however, is found to be less attractive because it adds substantial processing. Lagrange's method of undetermined multipliers can be also be used. This method finds the maximum or minimum of a function whose variables have an association but are not independent. Many Lagrange multipliers are added to the list of unknowns as a result of using this system and it is therefore a less attractive alternative of the solution to the problem than other methods.

### 6. The Patch Test

The patch is used as a tool to determine how an element will perform under different conditions by assessing the convergence of any finite element approximation imitation.[12] This test continues to be considered the most important check for practical finite element codes. Analysing the results of the test can also assess the degree of robustness of the finite element code, whereby an element which occasionally fails is termed non-robust. One additional use of the patch test was suggested by Babuska et al, whereby the test can establish the efficiency of gradient (stress) recovery processes which are important when estimating the errors, such as round-off.

When undertaking a patch test an arbitrary ‘patch' (cluster) of elements is selected from a larger mesh spanning an entire geometry. The mesh under consideration is typically fabricated of small elements, and consequently, the stress in all the neighbouring elements should be similar. It is mandatory that at least one node selected in the ‘patch' is within the volume rather than along the edge or on the surface of the geometry. The ‘patch elements' must also have adequate supports to prevent rigid body motion. To avoid spurious success an arbitrary state of stress is selected. Work equivalent loads, consistent with the assumed stress, are then applied to the boundary nodes of the ‘patch'. To assure convergence it is also necessary that the approximation fulfil both consistency and stability requirements[13].

* The consistency requirement ensures that as the size of the elements decreases, the model will still represent the exact differential equation and the boundary conditions.

* The stability requirement needs the solution of the mathematical model to be unique and avoid spurious mechanisms which may pollute the solution for all sizes of elements.

For a model to be successful in the patch test the results must exactly represent the correct constant stress in the loaded direction, within numerical error. All other stresses must be zero at all stress calculation points. In a practical finite element analysis problem, an element which has passed the patch test will approach results which agree with the mathematical model through convergence as the mesh is repeatedly refined. Therefore, before convergence the finite element model will disagree with its mathematical model because of an error known as discretization error, which tends to zero as the mesh is refined.[14]

Completing the patch test successfully ascertains that a valid element used in a mesh, rather than in isolation, will be able to display a state of constant strain and rigid body motion without strain. Success will also denote that when a state of constant strain occurs the element will be compatible with all the adjoining elements. Alternatively, the displacements of an incompatible valid element are only compatible when it is in a constant strain state. The element's incompatibility will transpire when there is a strain gradient. However, with repeated mesh refinement, the strain gradient across the element will become negligible and the element will act as a valid compatible element.

Although passing a patch test demonstrates that an element works, it does not take into consideration the rate of convergence during the test. Hence, a successful patch test does not establish how well the element will work; simply, it concludes that an element does work.

[1] O'Neil, P.V., Advanced Engineering Mathematics, third edition, pages 560, 566, 1006, 1008, 1251, and 1455

[2] Encyclopedia Britannica, http://www.britannica.com, [date visited 11/11/09]

[3] Reference made to Ottosen N., Petersson H., Introduction to the finite element method, 1992, pages 74-75

[4] Notation from Crouch R., 3D Finite Element Analysis Course Notes, 2009

[5] Reference made to Ottosen N., Petersson H., Introduction to the finite element method, 1992, pages 56-64 and 292-306

[6] Reference made to Zeinkiewicz, O.C., Taylor, R.L., Zhu, J.Z., The Finite Element Method: Its Basis & Fundamentals, pages 57-59

[7] Crouch R., 3D Finite Element Analysis Course Notes, 2009

[8] Rand, O. Rovenskii, V. Y. Analytical Methods in Anisotropic Elasticity: with symbolic computational tool, 2005, pages. 53-62

[9] William D., Callister, Jr., Materials Science and Engineering an Introduction, pages 64-66

[10] Cook, R.D., Concepts and Applications of Finite Element Analysis, Second Edition, pages 427 - 429

[11] Zeinkiewicz, O.C., Taylor, R.L., Zhu, J.Z., The Finite Element Method: Its Basis & Fundamentals, page 385

[12] Irons B., Ahmad S., Techniques of the Finite Elements, 1980, Page 151

[13] Zeinkiewicz, O.C., Taylor, R.L., Zhu, J.Z., The Finite Element Method: Its Basis & Fundamentals, pages 329-334

[14] Cook, R.D., Concepts and applications of Finite Element Analysis, Second Edition, 1981, pages 86-67