Failure models for metals – Implementation of the T-Failure criterion

19 02 2010


Why do materials actually fail? Is there a single definition of failure? Is there an approach to failure that will not fail itself from microscopic (atomic) level up to macroworld structures? In the last decade, large effort is dedicated to formulate material failure in terms of atomic dislocations, inclusions, microvoids, crack initiation and propagation, thermodynamics and macroscopic quantities such as stresses and strains. Unified approaches have emerged from these efforts, which compare well to experimental evidence. Lately, some of these approaches have been ported to the finite element method, particularly for crack propagation studies. In this chapter, a small excursion into material failure models is presented, the T- failure criterion is incorporated in Abaqus/Explicit and numerical verification is conducted. Also, a novel method for calculating threshold values for the T-failure criterion used for the simulations is presented.

In modern research and industrial applications, the estimation of an ultimate safety factor under extreme loading conditions, plays an important role in the overall safety assessment of a structure. In structural problems, where the response should be determined beyond the initiation of nonlinear material behaviour, fracture may be of profound importance for the determination of the integrity of the structure. Such problems include deep drawing, impact and crushing. Finite element methodologies, provide analysts with sophisticated tools for performing both material and geometry nonlinear calculations, which are yielding fairly accurate results. On the other hand, due to the lack of globally accepted fracture criteria, the determination of the structure’s damage, due to material failure using finite element methods, is still under intensive research. One could distinguish material failure in two broader categories depending on the scale in which the material is examined:

  • Failure is defined in terms of crack propagation and initiation. In this case, the microscopic behaviour of the material is examined, and respective theories are translated to finite element decoupling or even special element formulations that can reproduce crack propagation. In this case, fine meshes are required to capture the evolution of cracks. Such methodologies are useful for gaining insight in the cracking of specimens and simple structures under well defined global load distributions.
  • Failure is defined in terms of element removal. In this work, failure of an element is considered to occur when the element is no longer capable of carrying further loading, or equivalently, cannot store any further deformational energy. In this case, rather macroscopic failure criteria should be applied. When the conditions of such criteria are met, the element is deleted from the analysis. Such a failure mechanism is incorporated in most explicit time integration dynamic codes available today.

The second approach is most appealing for determining an upper limit to a particular loading on a structure. This work, provides a methodology for incorporating failure criteria in Abaqus/explicit, for metals following the von Mises yield conditions, as well as an application of the T-failure criterion [2] to [5]. Abaqus/Explicit is chosen because it is the only available explicit finite element code at NTUA that allows the implementation of user subroutines that describe material failure, in the sense defined earlier. It also requires the user to fully describe plasticity calculations. This has the drawback of requiring extensive benchmarking of the user subroutine, but on the other hand provides the user with the full stress and strain tensors at each time step. This makes the implementation of any failure criterion that exists or will occur much easier.

Read the rest of this entry »


An approach to the response Analysis of Shafts – Numerical Examples

18 02 2010

3Numerical examples

In the following examples, finite element simulations are used to demonstrate this concept. These simulations concern shaft loading with torsional moment and axial compression. This is a common loading condition for shafts, for example in rotating machinery used in power plants. For performing the simulations, the commercial finite element code Abaqus/Explicit 6 has been employed, using double precision arithmetic. Further, the POD processing is performed with custom software based on the LAPACK library 6. In all cases, POD analysis is performed over at least 300 snapshots (taking care to always observe the lower bound of the Nyquist–Shannon sampling theorem) including all the degrees of freedom of all nodes in the models. The results presented hereafter include the extracted mode shapes and plots of variable field distributions, amplitude vs. time fluctuation for each mode and Singular value percentage for each mode. It is noted here that in all cases multi-field POD is performed. This means that the input snapshot vectors include three displacement variables for each node. The resulting POMs include fields that correspond to each of the input variables. The amplitude vs. time curves for each mode represents the variation with time of the mode participation in the time-space domain of the simulation. In vibration problems this is usually an oscillating curve and a Fourier analysis is performed on this time variation so as to calculate the excited frequencies for each mode 6. Finally, the singular values relative percentage presents an overall estimated participation factor for the mode. Considering the fact that in the results there are three fields to each POM, it is interesting to calculate the norm of each field so as to determine its participation to the POM. In the results presented, the sum of the squared field norms equals unity and the field norms are presented next to the singular value percentage so as to easily determine which field is dominating the POM.

Read the rest of this entry »

An approach to the response Analysis of Shafts – Discussion

16 02 2010

4 Discussion

The study of shaft behaviour under dynamic loading and rotation is of profound importance in predicting resonance, in system control and system monitoring. The purpose of the work presented in this paper is to introduce the method of Proper Orthogonal Decomposition (POD) as a tool that can be effectively used in characterising the dynamics involved in the above tasks and extracting useful information for the real-time behaviour of the structure. The method is used in the form of meta-processing of finite element simulation results. It has been shown that the POD method is used to reproduce modes and corresponding frequencies that are systematically correlated under variations of initial conditions in dynamical problems of free vibration. Even when combined loading is applied, the POD method correctly discriminates and classifies these modes and frequencies. It has also been shown that the frequencies are affected by dynamical effects and pre-strain, a behaviour that is expected but is often difficult to calculate.The Proper Orthogonal Decomposition presents a considerable advantage: it is indifferent to the system that generates the input to the method. Nevertheless, it succeeds in extracting dominant modes and classifying them from the time-space response of the structure. Throughout this text, POMs have not been considered to coincide necessarily with natural modes of vibration. Rather, POMs are appropriate combinations of modes and therefore feasible configurations of a body. POMs, being orthogonal, they form a basis of the space where the configurations of the body in the particular process lie. Moreover, they are classified in an eigenvalue sense. These two properties are very important. In combination, they identify dominant POD modes in the response. This information can be interpreted in two ways, especially  in cases where a unilateral behaviour is desired. The first interpretation is that a strong dominant mode depicts a process that is consistent and “robust” to that mode. The second interpretation is that singular value percentage dispersion over more than one mode signifies a process that includes strong interference with the dominant mode.

Read the rest of this entry »

Is optimisation so bad?

27 07 2009

Tools are not bad, the way we use them may be. Let’s first take a look at the question: What is optimisation?

If we seek for a general definition of optimisation, this could be “to provide the best possible answer to a question according to a criterion.” What are the ingredients of this definition?

  1. Question
  2. Answer
  3. Criterion
  4. Optimisation method

The answer (2) depends on (1), (3) and (4) above. It is possible for (2) to be different (but not imperative) if you change (1), (3) and (4). Therefore optimisation per se is not a good or bad practice but merely a way of answering questions. For it to provide a reasonable answer, you need to use an appropriate criterion and method to a “well posed” question. “Well posed” is rather a motivation to make sure the question you ask is the one you really care about, i.e. it is sufficient to describe the problem you have in hand.

In a more mathematical manner, we could describe optimisation as follows:

“within a given range, find the extremum of a function, where the extremum satisfies some requirements”.  (a)

This implies that

  1. There is a function that describes our problem
  2. There is an extremum within the given range

For an analytical function, the extremum can be accurately calculated analytically.

For a numerical function an approximation can be given using numerical methods.

Such a numerical optimisation needs additional input, aka the accuracy to which the extremum is to be sought. Depending on the method used, the accuracy can be treated in a different manner.

Now, most problems can be translated to the statement (a). However, it is important to know what kind of problem one needs to solve and not confuse the definitions. Let’s take an example where equally interesting questions can yield different results that under some circumstances may be confused with each other.

Let’s suppose you’re standing on the edge of a cliff over a gorge and you want to get to the bottom of that gorge. This is a rather vague need and we might stipulate the following problems:

  1. Find the deepest spot in the gorge
  2. Find a sufficiently deep spot in the gorge
  3. Find the fastest/slowest/safest way to the deepest spot
  4. Find the fastest/slowest/safest way to a sufficiently deep spot
  5. If the spot you want to visit is known (you have seen it in a book) and maybe you even know where it is, find the fastest/slowest/safest way to that spot.
  6. Find a spot which is deep enough, but not that deep so that you can climb back up.
  7. Find a spot which is deep enough but on a bump where you can stand and take pictures, which bump is not that steep so you totter.
  8. If a friend is where you want to be, and he can see you with his binoculars but you cannot see him, you might want him to come and get you rather than you try to find him.

In some situations above a “sufficiently deep spot” is sought. A sufficiently good solution is sometimes the best solution to a problem, as it provides a good balance between the result and the effort to reach that result. For example let’s say you’re using the Navier-Stokes equations to choose which shape has minimum drag and maximum lift between a circle, a square, a line and an airfoil and perturbations of these. Any airfoil is a sufficiently good solution.

The CAE case

Computer Aided Engineering means using a numerical system to simulate a real system and take decisions based on the outcome of tests performed on the numerical system. One can only hope the numerical system simulates the actual system sufficiently well.

Assuming the numerical system can be a good enough approximation, one has to compose the right system for his problem. Considering the case of the airfoil, a designer could stop when he finds that among the attainable shapes, the airfoil is the one that provides more lift for the least drag. But he could go further and optimise the shape of the airfoil for a particular speed range and angle (or other service requirements). The result might be a quite peculiar cross-section. Let’s say that we construct one such cross-section and test it, and let’s assume that the results of the tests are identical to the simulation results. In this ideal case, the numerical or virtual system acts like a gauge to the real or physical system. But in order to achieve this result, a high precision manufacturing process was needed, for which the time length and cost of it might not be possible to justify.

The designer, answered the question

“Which is the airfoil that performs best under the given service conditions”?

The assumption here is that if the airfoil must be constructed, there is potentially infinite time and funding for trial and error until the result conforms with the simulation, given that the numerical system accurately simulates reality.

If a series of airfoils are to be manufactured though, a less time consuming and costly method is required. This additional restriction results in departure from the ideal case and landing into the reality of

System + measurement = variation

In this situation the “optimum” might change, which is normal, because manufacturing is a different problem to designing.

Now the question is

“Which is the airfoil that performs best under the given service conditions and can be manufactured with the particular process”?

Obviously this is a different problem. The system is no more the airfoil alone but the manufacturing is introduced as well. But this is still an optimisation problem.

So identifying the problem you actually want to solve is a first step to getting it solved in a satisfactory manner. Obvious statement but often ignored. Take another example: consider the case of a metal forming process, i.e. forming a part from a coil, and try to look for the optimal thickness of the coil so as to minimise weight. The possible coil thicknesses are 0.5, 0.6 and 0.7 mm. Let’s say a designer performs a simulation for each case and no cracks are identified for any of these. Let’s assume that an optimisation algorithm is also used which identifies 0.500 as the optimal thickness, with 0.499 giving cracks.  The designer would pick 0.5mm as the optimum thickness. So far so good, for what we asked we got an answer, there’s no statement here that we actually want to produce this part. Now, if we go to production, we might find out that the manufacturer gives us a variation of 0.01 for the coil thickness, which will result to scrap. Assuming that the manufacturer gives the same variation for all coils, the production engineer would pick 0.6 as the optimum thickness that eliminates cracks and minimises weight. Is this the best solution? It is good enough for the question asked. If the coil supplier is asked to provide statistical data in terms of a median and standard deviation for his coil we might calculate a scrap rate of a% for the thickness of 0.5. Then we should see if this scrap rate is acceptable and try to optimise the cost of production if we want to reduce it. That’s yet a different problem.

So, all in all, optimisation is a method for seeking the best solution for a given problem. As with all tools, one has to know how to use it. It is not bad or good. However, sometimes it is misused.

Response Analysis of Shafts – The POD method

8 06 2009

2 The POD method

The POD is a powerful method for system identification aiming at obtaining low-dimensional approximate descriptions for multidimensional systems. It provides a basis for the modal decomposition of a system of functions, usually data obtained from experiments, measurements or numerical simulations. The basis functions retrieved are called Proper Orthogonal Modes (POMs). From a different point of view, the method provides an efficient way of capturing the dominant components of a multidimensional system and represents it to the desired precision by using the relevant set of modes, effectively reducing the order of the system. The POD is interpreted in quite a few ways and appears in the literature in three equivalent forms: the Karhunen-Loeve decomposition (KLD), the Principal Component Analysis (PCA) and the Singular Value Decomposition (SVD). The latter is preferred in this text for its simplicity. The method is well known and its mathematical background and numerical implementation are explained in detail in many classic textbooks and papers. In practical applications, the input to the method is the distribution of a variable in time and space in matrix form, either measured or calculated from a simulation. For example, in structural response analyses, the distribution of nodal displacements over a body at consecutive time steps is used. The POD method provides then a decomposition of this time-space variation into Proper Orthogonal Modes (POMs) in space, fluctuation in time and participation factors. Moreover, since the method operates on the actual response, any complex dynamical effects are included and represented in the results as soon as the space and time discretisation is adequate. The advantage of the method is therefore, that it is indifferent to the system itself: it does not introduce any restrictions or assumptions on the mechanical or material response. My attempt here is to illustrating these unique properties with simple examples of continuum bodies, such as solid cylinders, which are nevertheless quite common and critical components in rotating machinery assemblies. The simplicity of the structures does not necessarily lead to equally simple responses.

In order to formulate the decomposition from the SVD aspect, it is assumed that the A is an m-by-n matrix defined over the field Ω, which is generally either the field of real or complex numbers. Then there exists a factorization of the form:

A=UV* (1)

where U is an m-by-m unitary matrix over Ω, the matrix is a m-by-n diagonal matrix of the positive singular values sorted in descending order and V* is the conjugate transpose of V, an n-by-n unitary matrix over Ω. U and V are orthogonal matrices and the columns u1,…,um of U yield an orthonormal basis of Ω m and the columns v1,…,vn of V yield an orthonormal basis of Ω n. The existence and uniqueness theorems of this decomposition have been extensively analyzed and can be found in the literature.

In the finite element context, the SVD can be used to reduce finite element matrices by expressing them using a much smaller set of approximation basis functions. In this way the solution of nonlinear finite element equations can be significantly accelerated.  But using the SVD purely for facilitating numerical computations is not its only application and indeed not the application of interest here. The SVD is a powerful tool that can be used to process finite element results just as any other measurement. From this aspect,  given a discretized deformable body under a set of external forces, its distinct configurations at times ti can be arranged in a matrix form A and analyzed using the SVD. The physical meaning of the three matrices produced by the decomposition of such a matrix has been extensively discussed. V contains the “empirical eigenvectors” or POMs of A, U contains the time coefficients of the corresponding POM and the singular values in can be thought of as the intensities or magnification factors of the participation of each mode in the whole process. For this latter case, the correlation of POMs with normal modes is largely argued. It is considered that only the first or dominant mode is approximated correctly while the deviation of the other modes depends on several parameters, most dominant being the sampling points and the structure of the mass matrix.

I will try to take another view on the POMs through the following simple observation: let matrix A be a row ordered matrix, where each row ai for i[1,m] is an n-sized vector containing the displacement values calculated at each node or sampling point of a structure. By definition, vector ai Ω, where Ω is the space of all possible configurations of a body. If K is the generalized stiffness matrix of a deformable body and ωi its eigenvectors,


On the other hand, using (1), POMs can be expressed as:

vi= (3)

(3) shows that vi can be expressed as linear combinations of ai. Since ai Ω, also vi Ω, this observation means that POMs are linear combinations of the eigenvectors ωi. A first conclusion from this observation is that POMs are feasible configurations of the deformable body. The space generated by the POMs, the latter being an orthogonal basis, is a subspace of Ω. The orientation of the basis is such that the POMs are those combinations of eigenvectors whose participation in the process described by A can be classified on the basis of the singular values. This observation is important in cases where the response of the body is so complex (due to loading, material and structural arrangement) that it is difficult to calculate normal modes or even impossible. POMs guarantee a decomposition of the response to modal forms, albeit linear. It remains to investigate the usefulness of this decomposition: “The attractiveness of the POD lies in the fact that it is a linear procedure. The mathematical theory behind it is the spectral theory of compact, self-adjoint operators. This robustness makes it a safe haven in the intimidating world of non-linearity; although this may not do the physical violence of linearization methods[14]”.

In most mechanical systems, each component that performs a particular function and is under a particular loading condition, there are only a few – possibly one or two – modes of ideal motion, deformation or both. These can be represented in a convenient coordinate system. This functional behaviour can be desirable or undesirable. A simulation of this process will provide a set of configurations that, when analysed with the POD method will possibly return a dominant mode as well as other modes. The subset of dominant modes corresponds to the functional behaviour of the component or structure. The singular values demonstrate the relative participation of each mode in the process. These features together, illustrate the consistency or “robustness” of the process with respect to the particular mode or modes. If singular values are comparable or in other words there is dissipation of the contribution between the modes, then the process is not robust with respect to the desired modes. Therefore the issue here is not whether POMs approximate normal modes, but rather whether the particular process is consistent in the sense that the simulated response approximates the desired response. Ideally, for example in transmission systems, one would expect a single mode that represents rigid body rotation or equally in a rotating coordinate system, no deformations. An analysis of such a simulation shows that the POMs and singular values quantify the departure from the ideal condition. This notion can be extended to situations where normal modes are difficult or impossible to calculate. Such a case may be a rotating shaft driving a ship propeller, where the effective mass of the water on the blades varies with the depth, and the propeller is prestrained under gravity, inertia, pressure and turbulence effects. Even in such complex situations, the desired behaviour of the system is known by requirements and comparing the departure from this behaviour under parameter changes can be very helpful in designing and controlling such systems. On the other hand, it may be desirable to mitigate the dominant mode and minimize the difference between the largest and smallest singular values. In either case, parametric studies can be used to investigate the convergence or divergence from the state in question and can be used as an effective tool in design and optimisation. Having seen that POMs are feasible configurations of a body, numerical computations are used to demonstrate that the different aspects in the loading of rotating machinery are actually represented and quantified by the modes.

To be continued…

An approach to the response Analysis of Shafts – Introduction

8 06 2009

1 Introduction

There are several situations where investigating the dynamical behaviour of shafts and attachments is important: resonance prediction is very important in the dimensioning and design of shafts but equally important is the prediction of shaft response in system control and system monitoring. The finite element method is frequently employed in numerical investigations aiming in analysing the influence of parameters such as section properties, materials and loading conditions. A few theories are presented in the literature that cope with particular problems. Such applications are presented hereinafter but certainly there should exist much more. Gubran studies the dynamic performance and cross-section deformation of shafts made of steel and aluminium, composites (CFRP and GFRP) and hybrids of metals and composites. A layered degenerated shell finite element with transverse shear deformation is employed. Results obtained show that improvement in dynamic performance and reduction of cross-section deformation of hybrid shafts over metallic and composite shafts is possible. Khulief et al identify drill-string vibration as one of the root causes of a deteriorated drilling performance. In order to demonstrate the complex vibrational mechanisms in a drilling system, as well as control its operation and improve its performance, they employ a dynamic model of the drill-string including both drill-pipe and drill-collars. The model accounts for gyroscopic effects, torsional/bending inertia coupling and the effect of the gravity. Modal transformations are applied to obtain a reduced order modal form of the dynamic equations.

Lee et al present a method for obtaining the out of balance response orbit of a gear-coupled two-shaft rotor-bearing system, based on finite element simulations. Bumps in the out of balance response are observed at the first torsional natural frequency because of the coupling between the lateral and torsional dynamics due to gear meshing. Neogy et al present work related to the stability analysis of finite element models of asymmetric rotors (non-axisymmetric shaft and non-isotropic bearings) in a rotating frame, which rotates synchronously about the undeformed centre-line of the rotor. Chatelet et al outline the necessity to predict at the design stage the dynamic behaviour of rotating parts of turbomachinery in order to be able to avoid resonant conditions at operating speeds. In their study, the global non-rotating mode shapes of flexible bladed disc-shaft assemblies are used in a modal analysis method for calculating the dynamic characteristics of the corresponding rotating system. The non-rotating mode shapes are computed using a finite element cyclic symmetry approach. Rotational effects, such as centrifugal stiffening and gyroscopic effects are taken into account. Turhan et al study the stability of parametrically excited torsional vibrations of shafts connected to mechanisms with position-dependent inertia. The shafts are considered to be torsionally elastic, distributed parameter systems and discretised through a finite element scheme. The mechanisms are modelled by a linearised Eksergian equation of motion. A general method of analysis is described and applied to examples with slider-crank and Scotch-yoke mechanisms. Li et al carry out finite element simulations using a speed-increase gearbox. Profile error curves are determined by the gear precision and the time-variable stiffness curves of the gears were acquired using 3D contact. After the stiffness excitation and error excitation are determined and the dynamic model of the speed-increase gearbox is established, the normal frequencies and vibration response of the whole gearbox and transmission shafts are analysed. Sekhar outlines the importance of dynamics and diagnostics of cracked rotors. In his study a model-based method is proposed for the on-line identification of two cracks in a rotor. The fault-induced change of the rotor system is taken into account by equivalent loads in the mathematical model. The rotor has been modelled using finite elements, while the cracks are considered as local flexibility changes. The cracks have been identified for their depths and locations on the shaft. The nature and symptoms of the crack are ascertained using FFT. Nandi presents a simple method of reduction for a finite element model of non-axisymmetric rotors on non-isotropic spring support in a rotating frame. In this frame the stiffness matrix, mass matrix and Coriolis matrix for the non-axisymmetric rotor (rotor with rectangular cross-section, cracked rotor, etc.) is independent of time but the support forces become periodic. Therefore, a large set of linear ordinary differential equations with periodic coefficients at support degrees of freedom is formed, which requires substantial computational effort to solve. To effectively handle this large system, a reduction method is introduced that keeps the essential information almost intact. Jun introduces Timoshenko beam theory for modelling shaft behaviour. Complex variables are used to represent the displacement, slope, moment and shear force. The complex transfer matrix between the variables at both ends of the shaft element is derived. The influence coefficients are analytically derived for the general flexible rotor having two resilient bearings at both ends. Modelling and derivation of the influence coefficients is based on the transfer matrix method. Simulated influence coefficients are compared to the results using the finite element method based on Timoshenko beam theory. Hu et al present a finite element-based formulation for modelling the dynamic behaviour of a rotating flexible shaft supported by a flexible structure. The coupling effect between the rigid-body rotation and the flexible deformation of the shaft is considered and represented by non-linear coupling terms in the mass matrix and forcing vectors in the global system of equations. The rigid-body rotation is treated as one of the degrees of freedom (d.o.f.) of the entire system. The interaction between the rotating shaft and the flexible support is modelled by either linear or non-linear springs distributed around the circumference of the shaft. The coupling between the flexibility of the shaft and the flexibility of the support structure are considered. The equations of motion are solved in the time domain using a modified Newmark scheme. Analyses are performed to validate the new development for different combinations of load condition, spring type, and rigid-body rotation. Mohiuddin et al present the finite element formulation of the dynamic model of a rotor-bearing system. The model accounts for the gyroscopic effects as well as the inertia coupling between bending and torsional deformations. A reduced order model is obtained using modal truncation. The modal transformation invokes the complex mode shapes of a general rotor system with gyroscopic effects and anisotropic bearings. The reduced modal form of the dynamic model is numerically simulated and the dynamic responses due to different excitations are obtained.

It is obvious that there is considerable research work devoted to the subject of rotating shafts and attachments and this is only part of it. Many of the afore-mentioned approaches involve considerable assumptions on the response of shafts and subsystems as well as for the involved physical phenomena. In this article I try to establish a way of dealing with such problems, either at the design stage or the monitoring stage by use of the Proper Orthogonal Decomposition (POD) as a tool that can provide insight to the behaviour of systems by appropriately processing data, either from simulations, experiments or measurements. Towards this, the ideal configuration of the system is established and deviations from this either in simulations or reality are quantified. The advantage is that no special algorithms or measurement devices are required by the method itself. Rather the mothod extracts such deviations from the available data. The quality and completeness of the result depends on the ohysical content of the data.

Next: the POD method