Master thesis: VARIATIONAL MECHANICS AND NUMERICAL METHODS FOR STRUCTURAL ANALYSIS

The past 25th June I had to present my (unexpected) master thesis in the UPC.
I want to thank my mentors Prof. Jaume Roset and Prof Vojko Kilar for their support and tutoring during the whole process of writing that led to a mark of 8 over 10 at the UPC Physics Department.
The following is an excerpt of the content of the document that will hopefully soon be published in the Physics department's webpage (http://fcia.masters.upc.edu/TFM/llista):

VARIATIONAL MECHANICS AND NUMERICAL METHODS FOR STRUCTURAL ANALYSIS

Abstract

This work focuses on the particular application of the variational principles of Lagrange and Hamilton for structural analysis. Different numerical methods are compared in their computation of the elastic energy through time.
According to variational mechanics, the difference between the stored elastic energy and the applied work should be null on each time step, so by computing this difference we can account for the level of accuracy of each combination of numerical methods. Moreover, in some situations when numerical instabilities are difficult to perceive due to high complexities, this procedure allows for the control and straightforward visualization of them, being an excellent source of hindsight on the behaviour of the analysed system.
The purpose of this dissertation is to present a scheme where the current numerical methods can be benchmarked in a qualitative as well as in a quantitative manner. It is shown how different combinations of methods, even for a simple model, can give very different results, particularly in the field of dynamics, where often also instabilites arise.
The first half of the thesis is a thorough explanation of these concepts and their application in terms of structural analysis. In the second part, a review on the numerical methods in general and of those implemented for our experiments is provided, followed by the experimental results and their interpretation. The model of choice, for simplicity and availability of analytical results is one cantilever column. Bending elastic energy of the column is monitored under transient regimes of different shapes, computing the total action of the system as its integral through time.

 Introduction

Targets and interest of our research

Variational mechanics date back as far as the XVIII century, when Leibniz, Euler, Maupertuis and eventually Lagrange devised the calculus of variations and the principle of least action. This methodology of treating physical phenomena is based on the notion that everything in Nature tends to a state of minimal energy(1).
In this work we will be focusing on its particular application in structural analysis, where one deals with “engineering scales” whose dimensions span between 100 times bigger or smaller than those of a human being. This is in contrast with other areas of applied physics like astronomy or molecular dynamics but will be shown how those variational principles still apply and even become powerful tools for the comprehension of the behaviour of our built environment.
Numerical methods, on the other hand, have proliferated since the 1950s alongside with the ever increasing power of computers as a means to simulate physical phenomena. This ceaseless growth in number and terminology has given place to a cumbersome mix of mathematical, physics and computer science often difficult to grasp.
Choosing one simple cantilever beam as our test model, we will utilize and compare different combinations of these methods to compute its elastic energy under transient loading regimes of different shapes. This will render useful in future research in nonlinear analysis of more complex structural systems.

Variational mechanics

According to the principles of variational mechanics (2) , the difference between the measured energy and the applied work should be minimal, so by accounting this difference in each time step of our simulations we should be able to infer the degree of accuracy provided by each combination and discuss the reasons that lead to differences in result using the energy as the natural norm for analysing the error(3).

Numerical methods for structural analysis

In a previous work by the authors (4)(5), it was shown how the vast amount of existing numerical methods can be grouped into three main sets according to the kind of physical phenomena they model and the type of differential equations they discretize: matter integration techniques (Partial Differential Equations), constraint integration techniques (Algebraic Differential Equations) and time integration techniques (Ordinary Differential Equations).
According to this, we will particularize in the following matter integration implementations: Finite Element (FEM), Finite Differences (FDM) and Mass Spring Systems (MSS). For the constraint integration we will be comparing Penalty Method (PM) and Lagrange Multipliers (LM). And for the time integration techniques we will employ Newmark Beta (NB), Houbolt's (HBT), and the Linear Acceleration Method (LAM).
Other combinations are also possible, as the proposed scheme is easily extensible, but for our current purposes it should suffice.

Discussion and future work

A numerical comparison of methods commonly employed in structural mechanics was presented.
It was made on the basis of energy principles and eventually the total action of a system under transient loading has been computed for each possible combination of methods.
It was shown how variational principles and an energetic norm can be employed in the benchmarking and assessment of the accuracy and stability of different implementations.
The scheme provided, tested on a simple example, is easily extensible to more complex systems with more elements. The advantage of this approach is that it allows for the monitoring of the global behaviour by means of one simple scalar, whose value is to be compared against that of an analytical computed from external forces or accelerations.
Also, a conceptual framework for the classification and treatment of numerical methods, grouping them into time, matter and constraint integrators, was used for the systematic analysis of the results.

Future work aims at the application of the same methodology in nonlinear analysis and more complex structures.
The combination with stochastic techniques for the integration of the action and the search of minimal energy states is one of the final targets of the current research.

References

1) C. Lanczos, 1970, “The variational principles of mechanics”, University of Toronto Press, Canada
2) W. Wunderlich, W Pilkey, 2002. “Mechanics of structures. Variational and computational methods”, CRC Press, pp.: 852-877
3) G. Bugeda, 1991, “Estimacion y correccion del error en el analisis estructural por el MEF”, Monografia 9, CIMNE, Barcelona
4) R. Andujar, J. Roset, V. Kilar, 2011. “Beyond Finite Element Method: An overview on physics simulation tools for structural engineers“, TTEM 3 / 2011. BiH.
5) R. Andujar, J. Roset, V. Kilar, 2011. “Interdisciplinary approach to numerical methods for structural dynamics”, WASJ Vol 14 Num.8, 2011. Iran.
6) J.H. Argyris, 1960, “Energy theorems and structural analysis”, Butterworths Scientific publications, London, UK.
7) R. Aguiar, 2005, “Analisis estatico de estructuras”, Colegio de ingenieros civiles de Pichincha, Quito, Ecuador.
8) H.P. Gavin, 2012, “Strain energy in linear elastic solids”, Duke University, Department of Civil and Enironmental Engineering, USA.
9) A. Nealen, M. Müller, R. Keiser, E. Boxerman, M. Carlson, “Physically Based Deformable Models in Computer Graphics” in Computer Graphics Forum, Vol. 25, issue 4, 2005
10) T.M. Wasfy, A.K. Noor, “Computational strategies for flexible multibody systems” in Appl. Mech. Rev. vol 56, no 6, Nov 2003
11) G.R. Liu, 2003, “Mesh free Methods: Moving Beyond Finite Element Methods”, CRC Press, USA.
12) J. S. Przemieniecki, 1968,  “Theory of matrix structural analysis”, McGraw-Hill. Inc, USA
13) G.D. Smit, 1978, “Numerical Solution of Partial Differential Equations by Finite Difference Methods”, 2nd ed. Oxford Applied Mathematics and Computing Science Series, UK.
14) M. Müller, 2008, “Real Time Physics course notes”, Siggraph USA.
15) K.J Bathe, 1995, “Finite element procedures in engineering analysis”, Prentice Hall, USA.

7. Key words
Finite Element, Finite Differences, Variational mechanics, Euler-Bernoulli beam

Nonlinearities in structural mechanics

Two years have already gone in the preparation of this thesis and now it seems in order to pass to the next part: the study of stochasticity.
Altogether with the preparation of our next article I am already summarizing and setting things up for the account for uncertainties in the structural design.
For such purpose, and after so much research into elasticity, dynamics and numerical methods, I have done the following scheme where most of the nonlinearities (which not surprisingly are often associated to randomness) can be fit into.
The scheme is based on the ubiquitous "Governing Equations" that serve as basis for every mechanician.

• System's internal balance equations (Statics or dynamics)
• Nonlinear damping
• Constitutive equations (Material)
• Plasticity
• Viscosity
• Creeping
• Hiperelasticity
• Kinematic equations (Geometry)
• Large deformations with small strains
• Large strains
• Buckling
• Instability
• Boundary conditions
• Forces
• Pressure
• Loads
• Wind
• Waves
• Friction
• Displacements
• Contact
• Link rupture

In the following months I will try to expand this list at the same time as I will try to understand how these nonlinearities are tackled in a deterministic manner.

The main point of my thesis is that stochastic methods are way more suitable for this purpose than deterministic ones...we will see...:)

Soon I will be getting back to my loved Ljubljana for a while, far from noisy Barcelona...

The Rayleigh-Ritz method explained for humans

Preliminaries
At last, I feel ready to explain the key concept needed to follow the steps in the Finite Element Method: the notion of STRAIN ENERGY!
The problem with this concept is that is as simple as it looks, ergo all the attached complications must come from somewhere else...
For the animation below I have chosen a classic spring to illustrate the concept (yes I know, is also a classic example...). Could have been anything else capable of storing energy, as a phone battery or a dam, but this seems more adequate given we are talking about structures.
The graph (click on it to see it moving) represents the relation between the applied force to stretch the spring and the actual displacement of the spring. So far, this is pure Hooke's law from secondary school (F=Kx)...
Also, the area that gets shadowed when the spring gets "tensioned" is, in mathematical terms, an integral of the drafted curve.

Displacement e vs Force F on a spring

of elastic material. Shaded area represents

stored strain energy. When pulled with

a force F its tip deforms an amount e.

When released it pulls back with the

same force F.

Given this particular, although very common case of a material (an elastic one), there is a straight tilted part in the curve that is called "elastic domain". If one keeps pulling the spring further then we still wouldn't break it, but instead we cannot expect it to recover its initial shape. The spring has plastified and we are in the horizontal part of the graph...but that is another story.
So, this "triangle" underneath the tilted part has an area, an integral, a scalar value in squared centimeters. The value of such integral is, not in centimeters but in N per meter or Joules (we are comparing force vs displacement), so is an Energy!

Applied force on the tip F vs potential energy

Displacement e vs potential energy

 And the shape of its curve, against the same displacement e or force F (click the animations) happens to be that of a parable (wasn't it obvious?)...
Well, and so far for the basics. Let's now think of a solid beam or column and put it in the place of our spring. When loaded, it would deform, would store elastic energy in its material ready to be released as soon as the applied force leaves. Normally things tend not to greedily accumulate energy but just the opposite. In Nature, a state of balance is that where things store a minimum amount of potential energy.
This basic principle in general applies for almost everything in life, and always does in physics. A drop of water in a cloud that has absorbed thermal energy and stored some of it as potential energy would release this as soon as it can, becoming rain when there is a lot of them, a capacitor in a circuit or even the adapter of a laptop store a differential of potential that make them useful, etc. This is the famous principle of minimum total potential energy (not to be mistaken with the principle of minimum energy, which also applies but in other fields).

Relating strain energies and forces

The way engineers have found a beam or column or piece of metal can store energy is by means of five types of loading: axial, bending, torsion, direct shear and traverse shear. Actually this is the way they have found they can stress things. As for every stress there should be associated a strain and hence some potential energy, this is tantamount as to say that they store energy. Though disputable, this is a convention and so we will accept it in the absence of a better one.

Resume of the strain energies associated to each type of applied force in engineering. Traverse shear has been omitted as is generally neglected.

 The good thing about this convention is that the job of assimilating stored potential energies for each type has already been done for us and put up in the shape of a formula. The table above gives the relation for each one of them in relation to the applied force and in relation to the caused strain.

The compressed beam under axial force F stores energy U of a value according to the formula.

The same beam, under flexural bending moment.

Similarly for shear

...and now for torsion

 The Rayleigh-Ritz method

And now we have all the concepts and are ready to work.
The Rayleigh–Ritz method (after Walther Ritz and Lord Rayleigh-Wikipedia), is considered a variational method. This is so because it is based in the calculus of variations. When applied to structural mechanics the varying value is the value of the strain potential energy, that is to be equilibrated against the work done by the external forces.
What we will try to do is to find what is the configuration, the state of the beam or group of beams, where the total energy of the system "external forces + internal stresses" gives a lowest value for the value of the energy.

It has eight fairly simple steps (this is also arguable, but this is how I understand it best) which I will try to picture as clearly as possible:

Step 1: Define the internal stored energy function (always an integral from the table above or a combination of them) and the external (a function of the caused displacement times the applied force).

The internal strain energy is the sum of all the possible combinations of internal strains from the table above

The external applied energy is the result of multiplying the external forces by the displacement

Step 2: Arbitrarily and happily (well, almost...) choose one function where the dimensions of the beam and the variable in the internal stored energy function get related. Multiply it by a dummy parameter (we will call it u overlined). The function chosen v is an approximation, so we overline it also to denote this.

The "dummy function" we will introduce into the energy formula can be any from above. The more elaborated, the higher the precision.

Two possibilities of approximation functions (without the U parameter). It is known that beams deflect quadratically.

Step 3: Solve, by means of integral calculus, the integral in step 1 with the chosen function in step 2. This will provide a relation between the dummy parameter and the internal energy. For the example below, where normally only bending moment is considered:

Laterally loaded beam

 
Internal strain energy due to the weight on the tip is given here as a function of the angle of rotation of the beam, but also is generally rounded as the second derivative of the displacement (v'').

Obtention of the second derivative of the "dummy" function

A very important remark: it is assumed from the beginning that the result will not be accurate!!

Step 4: Do the same for the external applied energy.

The "dummy" function and its second derivative substituted into the potential strain (Ui) and the external (Ue) energies.

Step 5: Add the two energies in order to get the total energy of the system, V (column internal + external forces)

Step 6: Applying the theorem of minimum total potential energy, the derivative of the total energy with respect to the parameter must be 0.

                                  

Step 7: Get the value of the parameter from the equation above.

Step 8: Substitute its value in the step 2 to get the real function of the displacement.

                                      

And that is it. For whatever value of x along the beam now we have a value of the deflection. Unfortunately the approximation chosen v(x)=u1(x/l)^2 is not very accurate, so a calculated displacement of the tip would diverge at least 25% from the analytical...but that is another story...
The most difficult steps are actually 1 and 2, where one should memorize the functions in order to master the method. The rest is easy algebra and very simple first order derivatives.
The Rayleigh-Ritz method eventually evolved into the nowadays ubiquitous Finite Element (did you recognize the shape function there?).

Well, I really hope this helps anyone who reads it as much as it helped me to write it...

Enjoy!

A brief summary on the stiffness method for structural engineering

Well, another year has gone by, and it has been a while since the last post. Life in Barcelona is being quite busy, so I am not dealing very well with my updates...This doesn't mean I am not researching, but apparently the topics I am getting to are deeper and deeper, so conclusions are getting harder to write...

The following is a summary of what I found here: www.civil.iitb.ac.in/~sghosh/CE317/history-stiffness-method.pdf

It is authored by the famous Zienkiewicz, so it must be a warranty of rigor. Out of it I have written my own summary, on pure chronological order, that goes as follows:

• 1862: Glebsch, Alfred wrote an algorithm for 3D trusses where no moments were involved. The models one would get out of it is as putting lots of springs connected by pins.
• 1883: Saint Venant, Andre wrote an article with comments on Glebsch's methodology.
• 1880: Manderla, Heinrich proposed an iterative method for the solution of the set of equations involved in pin-jointed trusses
• 1892: Mohr, Otto introduced the notion of "secondary stresses". He was referring to the moments. They were called secondary then as in steel and iron trusses the main ones were axial.
• 1910: Richardson provided the study of the stresses in a dam by means onf the Airy stress function.
• 1914: Bendixsen, Axel published a method similar to the slope deflection for structures with internal hyperstaticity.
• 1915: Wilson&Money proposed the slope deflection method, apparently unaware of Bendixsen's.
• 1922: Ĉaliŝev gave an iterative solution of successive approximations for frames without side-sway.
• 1930: Cross, Hardy published his famous method for redistributing unbalanced moments according to the internal stiffnesses of the composing bars. This provided the engineer with the necessary joint rotations.
• 1932: Grinter wrote a method for multi-story frames with side-sway called "method of successive corrections"
• 1935: Southwell tackled the problem of the constraints in a similar manner giving his famous method of "systematic relaxation of constraints", commonly known as dynamic relaxation.
• 1941: Courant, Richard explained his "variational methods for the solution of problems of equilibrium and vibrations". As he was a mathematician, the practical application of his work results still of some difficulty.
• 1943: McHenry in U.S. and Hrenikoff in Canada proposed a lattice analogy for the solution of stress problems so that continuum could be represented as a set of beams.
• 1944: Kron, Gabriel provided a complete algorithm in matrix form for 3D frames. In the 1930s it became fashionable to show the whole process of solving the linear equations in this systematic manner.
• 1954: Argyris wrote his "Energy theorems ans structural analysis", a comprehensive presentation of both the force and the displacement method for complicated airplane structures.
• 1967: Zienkiewicz presented his summary of the "Finite Element Method in structural Mechanics".