ELSEVIER

Contents lists available at ScienceDirect

Results in Physics

journal homepage: www.journals.elsevier.com/results-in-physics

Fatigue life prediction using multiaxial energy calculations with the mean stress effect to predict failure of linear and nonlinear elastic solids

Marko Nagode, Domen Seruga *

University of Ljubljana, Faculty of Mechanical Engineering, ASkerceva 6, SI-1000 Ljubljana, Slovenia

ARTICLE INFO ABSTRACT

An approach is presented that enables the calculation of elastic strain energy in linear and nonlinear elastic solids during arbitrary thermomechanical load cycles. The approach uses the simple fact that the variation of both strain and complementary energies always forms a rectangular shape in stress-strain space, hence integration is no longer required to calculate the energy. Furthermore, the approach considers the mean stress effect so that predictions of fatigue damage are more realistically representative of real-life experimental observations. By doing so, a parameter has been proposed to adjust the mean stress effect. This parameter a is based on the well-known Smith-Watson-Topper energy criterion, but allows consideration of other arbitrary mean stress effects, e.g. the Bergmann type criterion.

The approach has then been incorporated into a numerical method which can be applied to uniaxial and multiaxial, proportional and non-proportional loadings to predict fatigue damage. The end result of the method is the cyclic evolution of accumulated damage. Numerical examples show how the method presented in this paper could be applied to a nonlinear elastic material.

© 2016 The Author(s). Published by Elsevier B.V. This is an open access article under the CC BY-NC-ND

license (http://creativecommons.org/licenses/by-nc-nd/4XI/).

CrossMark

Article history: Received 8 April 2016 Accepted 10 June 2016 Available online 23 June 2016

Keywords: Fatigue Temperature Energy

Multiaxial stress Mean stress

Introduction

Mechanical components are usually subjected to variable loads during operation. These cause stresses, strains and temperature rises in the component as a reaction to such loading. Depending on the size of the load and the exposure (operating) time, fluctuating stress-strain fields in a component can eventually lead to a crack where the damage is greatest - the critical location, which will continue to grow under continued loading and eventually result in the failure of the component [1-5]. Although various mechanisms lead to the deterioration of mechanical products once they are put into operation, fatigue is still one of the main sources of failure for products that operate over longer amounts of time, e.g. months, years or hundreds of thousands of load cycles [6,7]. Fatigue mechanisms prosper due to the changeable load and environmental conditions and can ultimately lead to a complete stoppage of functionality of these products [7-11].

Predicting fatigue damage of a product is therefore not only important directly prior to manufacture but is an integral step in the early stages of product development. However, the identification of critical locations and the quality of the prediction, e.g. the predicted number of cycles to failure, will only be as good as the following: level of complexity of the temperature dependent

* Corresponding author. Fax: +386 1 2518 567. E-mail address: domen.seruga@fs.uni-lj.si (D. Seruga).

stress-strain calculation; reproducibility of the damage accumulation modelling; level of detail included in the fatigue damage prediction; and accuracy of the input data of the material properties [9,11]. The final experimental verification of the prediction is always valuable before the component enters the manufacturing stage but prior to this stage, computer aided predictions are necessary as a means of reducing financial outlay and shortening development times [6,12].

The majority of fatigue damage predictions are still based on uniaxial approaches or transformations of multiaxial stress-strain states into equivalent (uniaxial) cases either assuming a failure theory (e.g. signed von Mises stress) or applying the critical plane approach [5,12-17]. They usually give satisfactory predictions, especially if they incorporate various influences on the fatigue damage prediction such as e.g. the mean stress correction [5,1721]. However, under more complex conditions some of the commonly used techniques may no longer be capable of producing accurate predictions [16,17,19,21]. Alternatively, the invariance of the energy (which is independent of the coordinate system of observation) and its dissipation during cyclic loading have proven to be a suitable tool for predicting fatigue damage regardless of the type of loading (mechanical, thermal, uniaxial, multiaxial proportional or non-proportional) [13,18,22,23]. Therefore energy-based models for fatigue damage predictions have been a good counterweight to equivalent prediction models. However, according to the available literature, there have been attempts to include the mean

http://dx.doi.org/10.1016/j.rinp.2016.06.007 2211-3797/® 2016 The Author(s). Published by Elsevier B.V.

This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/).

Nomenclature

a mean stress parameter Dijkl nonlinear elastic stiffness

S Prandtl density Dppqq linear elastic stiffness

D damage operator i, j ; k ,l second and fourth order tensor indices

DSa, p elastic principal strain range i stress history index (only in Appendix B)

Aeae ,p equivalent elastic principal strain range j index for number of reversal points in residuum (only in

AUae ,p total equivalent elastic principal energy range Appendix B)

cycle sa uniaxial experimental strain amplitude j stress index of the spring-slider model

ea ,p elastic principal strain amplitude k temperature index of the spring-slider model

sij, skl elastic strain ns number of time increments

sp; sq elastic principal strain nu number of fictive yield energies

ores Sj,P ea ,p residuum principal strain nT number of temperature divisions

elastic principal strain amplitude p index of principal components

sae p equivalent linear elastic principal strain amplitude Re strain ratio

sae p equivalent nonlinear elastic principal strain amplitude Rr stress ratio

emax ,p maximum absolute principal strain s time index

eo ,p origin of elastic principal strain S logical operator (only in Appendix B)

Ocaycle uniaxial experimental stress amplitude t time

Oa ,p principal stress amplitude T temperature

Oae,p equivalent linear principal stress amplitude u yield surface

Oae p equivalent nonlinear principal stress amplitude Ua total elastic energy amplitude

stress Ua, p total elastic principal energy amplitude

, Ocycle Om uniaxial experimental mean stress Uae total equivalent elastic energy amplitude

Om p principal mean stress Ucycle uae total experimental elastic energy amplitude

Omax p maximum principal stress Uae ,p total equivalent elastic principal energy amplitude

Oo p origin of principal stress Udj back stress of the spring-slider model

Op principal stress Ue p total equivalent elastic principal energy

C elastic complementary energy per unit volume Ue total equivalent elastic energy

Ca elastic complementary energy amplitude Um total elastic mean energy

Ca,p elastic principal complementary energy amplitude Um ,p total elastic principal mean energy

Cae p equivalent elastic principal complementary energy Uo ,p origin of total equivalent elastic principal energy

amplitude W elastic strain energy per unit volume

Cm elastic mean complementary energy Wa elastic strain energy amplitude

Cm p elastic principal mean complementary energy Wa ,p elastic principal strain energy amplitude

d ,df equivalent cycle damage Wae equivalent elastic strain energy amplitude

D, Df accumulated fatigue damage Wae p equivalent elastic principal strain energy amplitude

Da damage due to total elastic energy amplitude Wm elastic mean strain energy

Dae damage due to total equivalent elastic energy amplitude Wm p elastic principal mean strain energy

Dm damage due to total elastic mean energy

stress correction in the energy-based methods e.g. [18,22], but to date no established criterion has been accepted as e.g. are the Sm ith-Watson-Topper (SWT) or Bergmann mean stress criteria for the uniaxial stress-strain states [5,16,18,24,25].

Here we present how energy-based fatigue damage predictions can be applied to a given variable multiaxial thermomechanical loading and nonlinear elastic solid, and hence show how they could be applied to materials such as metals, rubbers, polymer networks, liquid crystal elastomers and new biological materials under large strains. Furthermore, the approach is extended to consider the mean energy influence of the load cycles which can exactly reproduce the well-known uniaxial SWT correction [2426] or can be adapted for another experimentally observed influence on the mean stress level, i.e. a Bergmann type correction [25] by introducing an additional parameter a. The approach presented here is incorporated into a robust method based on Prandtl operators [6,10,12,23] that estimates the accumulated thermome-chanical fatigue damage at any time instant during the load history by calculating the cyclic fatigue damage evolution.

Energy calculation

The material response under cyclic loading is assumed to be temperature dependent and nonlinear elastic. This means that it

is independent of the load history (path independent) and that the stress tensor Oj and strain tensor skl form a nonlinear constitutive law

Oij(ts) = 0-ij(Diju(Ts), £kl(ts))

which depends on a temperature dependent stiffness tensor Dijkt and temperature Ts = T(ts) for every time instant ts;s = 1,... ,ns. It will be assumed here that both stress and strain tensors Oij and ekl for every time instant have been determined in advance according to a nonlinear elastic model as the material response modelling is not the main focus of the paper. As the approach is general, the stress and strain tensor can be considered as multiaxial and nonproportional. In the equations below, the time and temperature dependence will be omitted for simplicity though they are considered throughout the calculation. Additionally, all the quantities in this paper apply only to elastic materials, e.g., ski = seJ ,Ue = Uf, hence ''el" superscripts will be omitted for clarity.

First, strain energy and its complementary energy must be defined. These phenomena are crucial for the calculations to follow. Eqs. 1-13 refer to [27] where the reader can find further details on the strain and complementary energies.

For a given stress state rij, an infinitesimal amount of strain energy per unit volume dW (referred to as strain energy hereafter) during a load cycle can be calculated as

dW = (2)

Integrating the contributions of strain energy dW over a range from the beginning of loading to a given strain e«, the total strain energy can be calculated as

dC = Ej-dCTj. (9)

By integrating the differentials of the complementary energy dC over a range from the beginning of loading to a given stress Cy, the total complementary energy can be calculated as

W = / Oydfiy.

C = SijdOij.

As the material is elastic, the path independence can be considered and hence the differential of the strain energy can be written as

dW = d£j

and inserting Eq. (2) into Eq. (4) we obtain

r« — d£j = 0

' de«.

Eq. (6) suggests that the stress tensor is a constitutive law dependent function of the strain energy.

For every strain energy, complementary energy per unit volume (referred to as complementary energy hereafter) can be introduced (see Fig. 1) whereby

C = OijEij — W.

The complementary energy also depends only on the observed state therefore its differential can be calculated as

dC = do-j£y + Cyd£y - dW.

Inserting Eqs. (4) and (6) into Eq. (8), the differential of the complementary energy is observed as

Fig. 1. Principal stress and strain tensor components r1 and e1 between time instants t1 and t3 under changeable mechanical and thermal conditions. Strain energy W and complementary energy C form a rectangle at any time instant. Every cycle can be decomposed into amplitude and mean segments (shown for an arbitrary point between time instants t2 and t3).

As the constitutive law is path independent, the differential of the complementary energy can then be written as

dC = -0 dOij

and by placing Eq. (9) into Eq. (11) one can obtain

en - -

do« = 0

Eq. (13) proves that the strain tensor is a constitutive law dependent function of the complementary energy. Adding both strain and complementary energies of the ith and jth tensor components together, they always form a rectangular shape in any time instant as shown in Fig. 1.

Stress and strain tensors at any time instant can be expressed in the coordinate system of principal directions as a principal stress tensor ap and principal strain tensor ep, where p = 1,2,3 represents the indices of the principal values. Despite the transformation both strain and complementary energies remain the same as they are invariant, i.e. the amount of energy cannot be influenced by the coordinate system [27].

Next, the mean stress effect and the multiaxiality of the stressstrain states are addressed. The components of the principal stress tensor rp which occur at any point in a material during loading can always be decomposed into mean rm,p (static) and amplitude (dynamic) segments ra,p as shown in Fig. 1. The same is valid for the strain tensor components. For every principal component at any instant in time, the origins of both stress ro,p and strain eo,p can be determined. Thus by knowing the current state of a material, the correct calculation of the amplitudes ra,p and ea,p and the mean rm,p are possible.

According to the SWT criterion [18,24,25], the mean stress effect of an arbitrary load cycle in an observed time instant ts can be described as

rae,prae,p = rmax,p ra,p (14)

where rmaxp = rmp + rap is the maximum stress tensor, rap is the amplitude tensor and o*aep is the equivalent linear amplitude tensor originating from a fully reversed cycle without mean stress and resulting in the same damage as the combination of rmp and rap. Obviously Eq. (14) can only be valid for positive values of the maximum stress tensor component. The relations between equivalent linear stress and strain tensors are

rae,p = Dppqqeae,q; rm,p = Dppqqem,q and ra,p = Dppqqea,q (15)

hence Eq. (14) can be rewritten as the following expressions

0ae,pDPPqqeae,q ~ rm,pDppqqea,q + 0a,pDppqqea,q D—pqq0ae ,pDppqqeae ,q = D—pqq°m ,pDppqqea ,q + D—pqq0a ,pDppqqea, q rae ,pdppqqeae,q = rm ,pdppqqea ,q + ra ,pdppqqea ,q rae ,peae ,p = rm ,pea ,p + ra ,pea ,p

( 1 6)

which confirm the observation, that the background of the mean stress correction according to the Smith-Watson-Topper criterion

Fig. 2. Accumulated energies for linear and nonlinear material responses.

is the energy calculation and consists of the mean energy contribution and the energy amplitude contribution [18].

Equating areas considering either linear or nonlinear material response, it can be seen in Fig. 2 that

rae,p£ae ,p = rae ,peae,p (17)

which enables to apply the energy calculation from Eq. (16) also to a nonlinear material response. A similar approach for equating linear and nonlinear energies can be found in e.g. [28]. Substituting the left side of Eq. (17) with Eq. (16) gives

peae p — rm ,pea p + ra,pea.

Furthermore, if the ranges Aeae;P = 2sae>p and Aep = 2sa>p are considered rather than the amplitudes sae>p and eap, Eq. (18) can be multiplied by a factor of 2. Hence the total principal energy is expressed as

°p + v a,pDGp , p + Cm , p + W a ,

+ Ca ,

Uae ,p — Um ,p + Ua ,p

as shown in Fig. 3. The quantities oaepAsaep, ompAsp and oapAsp represent the total equivalent principal energy amplitude Uaep, the total mean principal energy Ump and the total principal energy amplitude Uap, respectively, each one containing information on strain and complementary energies W and C. As is the case for Eq.

and corresponds to the equivalent energy state of the observed material point with stress tensor rij and strain tensor skl.

Finally, according to Eqs. (19) and (21) we can conclude that damage due to the total equivalent energy amplitude Dae — D(Uae) can be given in terms of damage due to the total mean energy Dm — D(Um) and damage due to the total energy amplitude Da — D(Ua),

D(Uae)—D(Um)+D(Ua).

Considering the composition of the total energy from the strain and the complementary part as shown in Eq. (19), Eq. (22) becomes

D(Uae)—D(Wm)+D(Cm)+D(W a)+D(Ca).

Now considering that the displacement and damage of material are caused by strain energy only [29-31], the damage contribution of the complementary energies equals 0, so

D(Uae)—D(Wm)+D(Wa).

Equating Eq. (22) with (24) shows that damage contributions of the total and strain energies are equal,

(14), Eq. (19) is valid only if ffm.p + ffa.p > 0, otherwise D(Uae)—D(Wae) , D(Um)—D(Wm) and D(Ua )—D(Wa).

Ua,p — Um,p — Uae,p —

Adding the range of the total equivalent principal energy amplitude AUaep — 2Uaep to the current origin of total equivalent principal energy Uo p, the total equivalent principal energy Ue p can be obtained as

Ue ,p — Uo ,p + DUae,p.

By summing all the total equivalent principal energies, the total equivalent energy Ue can then be obtained as

Hence the total equivalent energy according to Eq. (21) can be directly used to predict fatigue life. It is shown in Appendix A that Eq. (25) is in accordance with the conventional strain energy calculation considering the contributions of mean and amplitude load levels.

Additionally, if a mean stress parameter a is introduced to the total energy calculation in Eq. (19), it is possible to adjust the importance of the mean stress effect by the equation

Ue — Ue

Uae ,p — aUm ,p + Ua p.

If a — 1, the Smith-Watson-Topper's mean stress correction is considered, if a — 0, no mean stress correction will be calculated.

2rae p£ae p — 2ffm n£a p + 2r

m p a p

a p a p

ae p ae p

Wae p + Cae p — W

For 0 < a < 1, an arbitrary mean stress correction is considered, e.g. the Bergmann type [25] which can be based on experimental observations, e.g. by fatigue tests at different mean stress values.

Energy life curves

S - N and e - N curves are gained from high-cycle fatigue (HCF) and low-cycle fatigue (LCF) tests, respectively. The controlled variables, i.e. stress in HCF tests and strain in LCF tests, alternate between limit states continuously until failure of the specimen according to the stress ratio Rr or the strain ratio Re, respectively. As the material response is nonlinear elastic as well as different in tension and compression, the stress and strain ratios are usually not equal. Energy life curves are created from r - e curves and either e - N or S - N curves at distinct constant temperatures. If the e - N curve for a strain ratio Re and temperature T is available, the number of cycles to failure N can be obtained for a given strain e as shown in Fig. 4(a). From the r - e curve in Fig. 4(b), the stress r for a given strain e is retrieved and both the mean energy Um and the energy amplitude Ua are calculated according to Eq. (19). Considering Eqs. (19) and (20) for the experimental values of stress and strain as well as the chosen mean stress parameter a, the equivalent energy Uycle can be calculated from

- ^ Re = -1 (27)

raea + armea otherwise.

uaycle = 2

It turns out that U - N curves which follow from this calculation should be transformed into energy life U - d curves as

d = 1 =N (28)

to avoid numerical difficulties due to unlimited number of cycles to failure for equivalent energies beneath the endurance limit. The equivalent cycle damage d is always bounded between the values of 0 and 1 and is therefore more convenient than the number of cycles to failure for further analyses. The energy life U - d curve is presented in Fig. 4(c). However it should be noted that caution must be taken when determining ea, ra and rm. The mean stress parameter a must be the same as that used during the energy calculation otherwise erroneous fatigue damage predictions may occur. Eq. (27) is also valid only if rm + ra > 0, otherwise Uaycle = 0.

If the S - N curve for the stress ratio Rr at temperature T is available, the number of cycles to failure N is obtained for a given stress r. The strain e and strain ratio Re are then obtained from the experimental r - e curve. The rest of the transformation into a U - d curve is the same as for the e - N curve. The transformation shown in Fig. 4 is mandatory for each test temperature.

Continuous energy and damage calculations

For simple load cycles the calculations of the total mean and amplitude energies are fairly straightforward. But if real-life situations of arbitrarily loaded components are being analysed, their values are not so obvious. Therefore the approach has been incorporated into a method which enables the continuous calculation of rm,p and ra,p based on memory rules that provide information of the equivalent elastic stress and strain amplitude, and the total energy at every time instant. A flow chart illustrating the stages of the method is given in Fig. 5.

The method is similar to the one presented in [26] but with an additional outer loop to consider all principal directions and is then followed by a damage calculation based on Prandtl operators. A pseudo code illustrating the continuous total energy calculation is given in Appendix B.

The procedure depicted in Fig. 5 is carried out as follows:

1. Stress and strain tensors rj(ts) and ekl(ts), the mean stress parameter a as well as uniaxial r - e and S - N or e - N curves are required as input data.

2. U - d curves for test temperatures are created (as outlined in Section ''Energy life curves").

3. U - d curves are discretised into fictive yield energies uj = 1,..., nu which allows the utilisation of the Prandtl spring-slider model in Fig. 6 for the fatigue damage calculation. A similar procedure can be found in [6,12].

4. Prandtl densities dj(Tk); j = 1,..., nu; k = 1,...,nT are determined as

W-j^fP-X <^+1 - ui))

5. Principal stresses rp(ts) and principal strains ep(ts) are then searched for.

6. For every principal component p and time instant ts origins of both stress rop(ts) and strain eo ,p(ts), and the origin of total equivalent principal energy Uo p(ts) are determined, indicating the current cyclic state of the material based on the identification of rainflow and Clormann-Seeger load cycles [26]. This allows calculation of the amplitudes rap(ts) and eap(ts) and the mean rm p(ts).

7. The total equivalent principal energy amplitude Uaep(ts) is calculated from Eq. (26).

8. The total equivalent principal energy Ue p(ts) is calculated using Eq. (20).

Fig. 4. (a) e - N curve for Re = 0 [32], (b) r - e response [8] with calculated energies and (c) transformed U - d curve for a = 1.

Fig. 5. Incorporation of the approach into a method for continuous energy and damage calculations.

9. The total equivalent energy Ue(ts) is gained from Eq. (21).

10. The Prandtl spring-slider model is loaded by the total equivalent energy Ue(ts) and the back stresses U5j(ts) on individual springs are calculated as

USj(ts) — max jUe(ts) - Uj, min jUe(ts) + ujVjT)-U^)} J.

11. The contribution of each spring-slider segment to the cyclic damage evolution is given by

Dj(ts) — Sj(Ts)USj(ts) (31)

and the damage operator D(ts) representing the cyclic damage evolution as

D(ts) — X>(ts). (32)

Fig. 7. Uniaxial load history for node 321. A fragment of a mechanical component (shown as a meshed solid block in reversal points) is arbitrarily loaded in the vertical direction (uz). The total displacement |u| is given in colours. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

12. Accumulated fatigue damage Df (ts) is expressed as the variation of the damage operator over time as

Df (ts)^|D(ti)-D(ti_1)|.

Compared to the conventional energy calculation, the proposed method does not involve integration so it is more time efficient, which is particularly beneficial if FE models with a vast number of nodes and elements and longer load histories are being analysed. Considering the path independence of the energy, only the reversal points can be used for the energy calculation. In this case however, the nonlinear information between reversal points is irretrievable. As the strain and complementary energies always form a rectangular shape in space, the procedure is valid for any arbitrary nonlinear elastic material response.

Numerical examples and discussion

A solid block (Fig. 7), representing a fragment of an arbitrary mechanical component, consisting of 4500 brick elements and 5376 nodes, has been chosen to illustrate the method presented in Section ''Continuous energy and damage calculations". The material used for numerical examples is a blend of filled natural rubber, poly-butadiene rubber and styrene-butadiene rubber, syn-thesised via free-radical polymerisation as an emulsion in water. It is used in industry as a part of the composite that is used for producing air spring bellows. For further details on the material used in this study the reader is referred to [8,32]. A two-parameter Mooney-Rivlin material model [8] has been used to describe the nonlinear elastic stress-strain response of the rubber in a finite element calculation. The r - e curve and the e - N curve with its energy equivalent U - d curve are given in Fig. 4. Two load histories have been considered: a uniaxial load history (Fig. 7), to give a simple demonstration of how the method might be used and to illustrate the relation with the uniaxial SWT parameter; and a multiaxial load history (Fig. 8) to give a fuller numerical example demonstrating the ability of the method. Node 321 and element 286, both indicated in Fig. 7, have been chosen for fatigue analysis for both load histories. For simplicity, a constant temperature of 25 °C is assumed throughout.

The base of the block is fixed in its normal direction whereas the block can freely move in either transverse direction. The top side of the block is loaded with displacements in normal and/or transverse directions as shown in Figs. 7 and 8. The displacements of node 321 which represent the variable load history are also given in Figs. 7 and 8.

The calculated stress-strain histories, and the corresponding principal values and the energy calculation are given in Figs. 912. Whilst the convention on sorted principal stresses ffi P r2 P r3 is used here, it is not mandatory. Principal stresses and strains in element 286 are given in Fig. 9 for the uniaxial case and in Fig. 10 for the multiaxial case. The energy variation and the

principal direction 1

principal direction 2

principal direction 3

ra Û.

A lu 5P? \

1 0

0 1 2 3 4 5 7

0.8 0.6 0.4 0.2 0 -0.2 -0.4 -0.6 -0.8

¿0,1 0

1 T

u 6

0 12 3 4 5

8 6 4 2 0 -2 -4

8 4 0 -4 -8 -12 -16

i 1 — cra. »

— ~ — o;„ i \ \ \ V

/ /

/: "-i a M 1 X \ i

V V V/ \

0 1 2 3 4 5 6 7

i f/,i pxx

— -o— Ti

ra» ra» ra» ra» ¡roo« ra» ra» rax»

2 3 4 5 6 7

0.8 0.6 0.4 0.2 0 -0.2 -0.4 -0.6 -0.8

4 2 0 -2 -4

2 3 4 5 6 7

ri p „P

ÇDCC y

s es!* tro s ¡roo \ /

íux í \ L

1 0 1 2 3 4 5 6 7

1 1 - Ca.2

- — — - Olli.

1 t/e, 2

uro UDO« ra» ra» ra» rax* ra» rax»

2 0 -2 -4

- O— 33

0 i 2 3 4 5 6 7

0.8 0.6 0.4 0.2 0 -0.2 -0.4 -0.6 -0.8

rV > 0

/s A /

bn 1 A, L

¡roo;

0 1 2 3 4 5 6 7

6 4 2 0 -2 -4

2 1.6 1.2 0.8 0.4 0

-0.4 -0.8 -1.2 -1.6 -2

2 3 4 5 t[s]

2 3 4 5 t|s|

2 3 4 5 t[s]

■ i -3

— — - Cm.3

\ N 1

0 i 2 3 4 5 6 7

1 t/e, 3

0 ^0,3

< ra» ra» rax*

raxx icq

Fig. 9. Energy calculation for the uniaxial load history. Stresses, strains and their origins together with the mean stresses, stress amplitudes and the total principal energies with their origins are calculated for t = i; i = 1, 2,..., 8 s. The number of circles represents the number of steps performed in the analysis.

predicted damage of the chosen component are given in Fig. 11 for the uniaxial case and in Fig. 12 for the multiaxial case.

An initial comparison of Figs. 9 and 10 shows that more steps are needed for a multiaxial load history FE analysis to converge

than for the uniaxial case (the number of circles in Figs. 9 and 10 represents the number of steps performed in the analysis). Displacements which are transverse to the loading direction occur in the uniaxial load case as a consequence of Poisson's ratio, but

Fig. 10. Energy calculation for the multiaxial load history. Stresses, strains and their origins together with the mean stresses, stress amplitudes and the total principal energies with their origins are calculated for t — i; i — 1,2 ,... ,8 s. The number of circles represents the number of steps performed in the analysis.

produce no additional stresses in the material. In the multiaxial load case however, displacements in the lateral directions will produce additional normal and shear stresses. In the uniaxial case, only one principal stress is strongly expressed; that is the

maximum principal stress in tension between time intervals of 0 to 3 and 7 to 8 s and the minimum principal stress in compression between 3 to 7 s. In the multiaxial case all three principal stresses will simultaneously occur depending on the loading.

Fig. 11. Energy variation and fatigue damage calculation for the uniaxial load history.

The origins for stresses, strains and equivalent energies shown in Figs. 11 and 12, allow a visual check of the mean and amplitude calculation. The first principal stress r in Fig. 9 increases from 0 to

1 s whereas its origin 0"o>1 decreases for the same amount in the opposite direction as an alternating cycle is assumed at the beginning. To confirm this it can be seen that the amplitude ra1 is rising from 0 to 1 s whereas the mean stress rm1 remains at zero. There is a reversal point at 1 s so r1 starts decreasing from here on and the origin ro1 is set to the value of r1 at 1 s. At 2 s the origin ro1 is set to the value of r1 at 2 s as r1 starts increasing again. When the first principal stress r1 exceeds its value at 1 s, a closed loop starting at

2 s is obtained so the material follows the original path which started at the beginning (t = 0). Consequently the origin jumps to the opposite side, the stress amplitude increases and the mean stress value is set to 0. Similar behaviour is noticed throughout the history for the first, second and third principal directions. Furthermore, the stress, strain and equivalent energy calculations for the multiaxial load history given in Fig. 10 follow this same behaviour, i.e. bigger closed loops of the first principal stress are obtained just before 3, 6 and 8 s and bigger closed loops of the third principal stress are obtained just before 3 and 7 s. Jumps in origins of stresses and equivalent energies may occur during the amplitude and mean calculation as is observed for the first principal stress r1 between 4 and 5 s (Fig. 9). The reason for the jumps of ori-

gins lies in the strain-dependence of reversal points as the strain is their driving quantity [26], and it has been shown that by knowing the reversal points of the strain alone allows for the identification of not only the strain amplitude, but also the mean stress and the stress amplitude. Here this also applies to the calculation of the principal energy amplitude and principal mean energy. When e1 starts decreasing at 4 s in Fig. 9, its origin lies at the level of e1 at

4 s so the origins of r1 and Ue>1 also lie at their respective values of r1 and Ue1 at 4 s. But when e1 reaches a reversal point before

5 s, a small nested cycle is closed and the origin is therefore shifted to the previous reversal point, in this case to the one at 3 s. According to the strain-driven reversal points the origins of r1 and Ue1 at this particular moment also change to the values of the reversal points at 3 s.

From Figs. 11 and 12 it can be seen that the equivalent energy of the block increases or decreases not only depending on the product of the principal stresses and strains, but also the mean stress influence. On the contrary, the accumulated fatigue damage Df always increases and is therefore a good indicator of the damaged state of the block. In the uniaxial case where the second and for some time the third principal stresses remain at zero, though the principal strains vary throughout the analysed load history, the equivalent energy as a product of these quantities remains zero and does not contribute to the total damage. The fatigue damage in the uniaxial case is therefore only due to either the first principal stress and strain when the solid is loaded in tension; or the third principal stress and strain when loaded in compression. This is in accordance with uniaxial stress-strain theory.

In the multiaxial case however, the second and the third principal strains are not merely a consequence of the block contracting during loading but also actively producing stresses as the solid is loaded in all directions. Higher loads also cause higher damage which can easily be seen by comparing the stress, strain, energy and damage jumps between 6 and 7 s with the remainder of the loading for the uniaxial case in Figs. 9 and 11 and between 2 to 3 s and 6 to 7 s for the multiaxial case in Figs. 10 and 12.

The influence of the mean stress level can be noted in Figs. 11 and 12. Two limit contributions are depicted, with a = 0 (no mean stress correction, thin lines) and a = 1 (mean stress correction, thick lines). If no mean stress is considered, smaller fluctuations of the energy variations Ue are observed as compared to the considered mean stress correction. But if the mean stress is not considered for the determination of the U - d curve (Fig. 4), the amplitudes become more influential. Therefore more than 5 times higher fatigue damage is observed if the mean stress is not considered throughout the calculation in both, the uniaxial and the multiaxial case. Furthermore, the simulated multiaxial history would cause approximately 35% more damage compared to the simulated uniaxial history. Setting a between 0 and 1 changes the amount of accumulated fatigue damage.

Conclusions

The approach presented here enables energy-based fatigue damage calculations including the mean stress effect that resembles the equivalent uniaxial approach with the Smith-Watson-Top per or the Bergmann damage parameter. If the mean stress effect is not considered, the approach produces the same result as the conventional energy-based approach but in a significantly shorter time as the calculation rests on the subtraction of two energy states rather than the integration process between them. The mean stress effect is considered with the parameter a. Numerical examples show high dependence of the fatigue damage on the consideration of the mean stress effect in both uniaxial and multiaxial cases. Further validations of the method are on-going.

while p < 3 do

i = Uj = Oiemj,= 0 while i < n do

res _ res

j,P I" p^'i' ~J,p\

°i-\,P "J-2.P - \°j.p °j-l,P

if then

*max,p e0Jt,) = -s (0; = Uop(tl) = -2sign(ep{t,))ap(t,)sp(ti)

s (t ) = e"\ : a (t ) = <rrcs, ; U (t) = U"\

0,pV I s y-2,/7' 0 I / 7-2,/7' is 7-2,/7

end if

j = j-2;S = j>2 else

if 7 > 2 and e™ pe™ <0 and ft ) <0 and >|e"s, p| and |s"slj;| >| then

= |*,ft)|; ^ft) = -^ft); = -o-pft); ^ft) = -2sign(/,»<7, ft )£,(/,)

j = j-2; S = false else

if »hen

= MM; fio,„ft) = -Sft); CTo,,ft) = -Opft); ^,,ft) = -2signft,ft))<T,ft)^ft); -S = false

«„.,(/,) = C; CTo,„ft) = t^ft) = s=f^se

end if

^m»,i,ft)=^m>Pft)+KPft)|

if <wft)>0then

Um.p ft ) = 2fm>p ft ) ft )| f-j-ft) = sign(fia,„ft))(^ +ccUm„) else

end if

^ft) = ^ft) + 2t/„,„ft) if S = false then if f^ft) is reversal point then

7=7 + 1; «£=*,(/,); = MO; UZ =t/e,,ft) end if i = / +1 end if end while p = p +1 end while i = l;

while i<ndo

^eft) = ^,lft) + ^,2ft) + ^,3ft) /' = ;' +1 end while

Fig. 13. Pseudo code illustrating the continuous total energy calculation.

Acknowledgement

We thank the Faculty of Mechanical Engineering, University of Ljubljana, for financial support.

Appendix A. Equivalence with the SWT parameter

If the total equivalent principal strain energy amplitude of a load cycle Waep consists of the principal strain energy amplitude Wap and the mean principal strain energy Wmp, the following relation holds

W ae p — Wap + Wm>p.

By the definition of strain energy in Eq. (3) and assuming a constant temperature it can be shown that the principal strain energy amplitude can be calculated as

f&a q T-const f&a>q

Wap — Wa>p(£a,q) — Oa,pd£a,q ——ns/ Dpp„£a,qd£a,q

-1D -1

— 2 Dppqq£a,q£a,q — 2 Oa,P£a,q

and the mean principal strain energy as

1 /*£m q+£a q

Wm p — Wm ,p(£q) — -

J T— C£I

,p — Wm ,p(£q) — 4 I Opd£q — 4

^ £m, q—£a, q

„ 1 f£m,q+£a q

t—const 1 I D £ d£

A-/ppqqcqucq

' £m, q-£a,q

Wm ,p — g Dppqq(^ (£m, q + £a, q) (£m, q £a, q) ^ Wm, p — Dwqq( £m q + 2£m,q£a,q + £2! q — ( Ä q — 2£m,q£(

m ,p — g^ppqq^ °m ,q + ^ °m ,q + ^ ,q ^ °m ,q * °m ,q + ,q

W m ,p — Ö Dppqq£m, q£a, q

Wm ,p — 2 rm ,pea, q•

With dea p and dsa>q being interchangeable integration variables and p,q = 1,2 ,3, Eqs. (A.2) and (A.3) can then be rewritten as 1

Wa p = 2 ra ,p£a p

Wmp — 2 rmp£ap

The total equivalent principal strain energy amplitude according to the general expression in Eq. (A.1) is thus

Waep — 2 ra,P£aP ^ 2 rm,P£a,P

which is equal to the Smith-Watson-Topper mean stress criterion given by Eq. (16), rearranged and divided by 2.

Appendix B. Method for continuous total equivalent energy calculation

Steps 6-9 of the method for calculating the continuous total equivalent energy as presented in Section ''Continuous energy and damage calculations" are given in Fig. 13. The strain reversal points are calculated prior to initiation of the counters p, i, and j in lines 1 and 3. The counters represent the current principal direction, time instant and the number of reversal points in the residuum, respectively. The maximum absolute principal strain £maxp is also initiated.

Next, for every principal direction in line 2 all the data points in the load history in line 4 are processed. When a data point is processed it is moved to the residuum.

If the number of reversal points in residuum j > 2, line 5 is checked to see whether residuum strains sj"i:s2 p, jl p, ej^ and strain Sp(tj) form a rainflow cycle [33]. If a rainflow cycle is not found then the Clormann-Seeger cycle [34] is searched for in line 13. If j 6 2 nor a rainflow cycle, nor a Clormann-Seeger cycle is found, line 17 is checked to see if ep(tt) has exceeded the current maximum absolute principal strain - this condition is always fulfilled for the first point in the load history (emaxp = 0). If none of the above conditions is met, strain origin eop(tt), stress origin ffo>p(tj) and the origin of total equivalent principal energy Uo>p(tj) are set to the jth residuum in line 20. The amplitude and mean values of stress and strain are determined in line 22. If the maximum stress in line 23,

Omax p (ti) — Omp(ti) + |<7aj)(ti)|,

is greater than zero, the total equivalent principal energy amplitude in line 27 equals

Uae p(ti) — sign(£a,p(ti))(Ua,p + aUm p)

otherwise Uaep(ti) = 0 in line 29. The total equivalent principal energy in line 31 referring to time ti is then

Uep (ti) — Uop(ti) + 2Uaep(ti ).

Logical operator S is a switch that defines which point represents the origin in a given time instant. If it is set to false in line 32, line 33 is checked as to whether ep(tt) is a reversal point. The first and last points in the load history are also treated as reversal points. If ep(ti) is a reversal point, then j is incremented and residuum strain ejps, residuum stress cj^ and residuum total equivalent principal energy are set to the i - th values in line 34. After all the data points have been processed, the total equivalent energy is calculated in line 43.

References

[1] Shiri S, Yazdani M, Pourgol-Mohammad M. A fatigue damage accumulation model based on stiffness degradation of composite materials. Mater Des 2015;88:1290-5.

[2] Verron E, Andriyana A. Definition of a new predictor for multiaxial fatigue crack nucleation in rubber. J Mech Phys Solids 2008;56:417-43.

[3] Straffelini G, Benedetti M, Fontanari V. Damage evolution in sinter-hardening powder-metallurgy steels during tensile and fatigue loading. Mater Des 2014;61:101-8.

[4] Cho JR, Yoon YH, Seo CW, Kim YG. Fatigue life assessment of fabric braided composite rubber hose in complicated large deformation cyclic motion. Finite Elem Anal Des 2015;100:65-76.

[5] Liakat M, Khonsari MM. An experimental approach to estimate damage and remaining life of metals under uniaxial fatigue loading. Mater Des 2014;57:289-97.

[6] Seruga D, Hansenne E, Haesen V, Nagode M. Durability prediction of EN1.4512 exhaust mufflers under thermomechanical loading. Int J Mech Sci 2014;84:199-207.

[7] Bazant ZP, Hubler MH. Theory of cyclic creep of concrete based on Paris law for fatigue growth of subcritical microcracks. J Mech Phys Solids 2014;63:187-200.

[8] Oman S, Nagode M. On the influence of the cord angle on air-spring fatigue life. Eng Fail Anal 2013;27:61-73.

[9] Shiri S, Pourgol-Mohammad M, Yazdani M. Effect of strength dispersion on fatigue life prediction of composites under two-stage loading. Mater Des 2015;65:1189-95.

[10] Nagode M, Längler F, Hack M. A time-dependent damage operator approach to thermo-mechanical fatigue of Ni-resist D-5S. IntJ Fatigue 2011;33:692-9.

[11] Yazdanipour M, Pourgol-Mohammad M. Stochastic fatigue crack growth analysis of metallic structures under multiple thermal-mechanical stress levels. Mater Des 2016;95(5):599-611.

[12] Nagode M, Seruga D, Hack M, Hansenne E. Damage operator-based lifetime calculation under thermomechanical fatigue and creep for application on Uginox F12T EN 1.4512 exhaust downpipes. Strain 2012;48(3):198-207.

[13] Skelton RP. The energy density exhaustion method for assessing the creep-fatigue lives of specimens and components. Mater High Temp 2013;30 (3):183-201.

[14] Zhao T, Jiang Y. Fatigue of 7075-T651 aluminum alloy. Int J Fatigue 2008;30:834-49.

[15] Meggiolaro MA, de Castro JTP. An improved multiaxial rainflow algorithm for non-proportional stress or strain histories - Part I: enclosing surface methods. Int J Fatigue 2012;42:217-26.

[16] Meggiolaro MA, de Castro JTP. An improved multiaxial rainflow algorithm for non-proportional stress or strain histories - Part II: the modified Wang-Brown method. IntJ Fatigue 2012;42:194-206.

[17] Gates NR, Fatemi A. Fatigue life of 2024-T3 aluminum under variable amplitude multiaxial loadings: experimental results and predictions. Proc Eng 2015;101:159-68.

[18] Ince A, Glinka G. A modification of Morrow and Smith-Watson-Topper mean stress correction models. Fatigue Fract Eng Mater Struct 2011;34(11):854-67.

[19] Park J, Nelson D. Evaluation of an energy-based approach and a critical plane approach for predicting constant amplitude multiaxial fatigue life. Int J Fatigue 2000;22(1):23-39.

[20] Dowling NE. Mean stress effects in stress-life and strain-life fatigue, SAE Technical Paper 2004-01-2227, 2004.

[21] Dowling NE, Calhoun CA, Arcari A. Mean stress effects in stress-life fatigue and the Walker equation. Fatigue Fract Eng Mater Struct 2009;32: 163-79.

[22] Jahed H, Varvani-Farahani A, Noban M, Khalaji I. An energy-based fatigue life assessment model for various metallic materials under proportional and nonproportional loading conditions. Int J Fatigue 2007;29(4):647-55.

[23] Gosar A, Nagode M. Energy dissipation under thermomechanical fatigue loading. IntJ Fatigue 2012;43:160-7.

[24] Smith KN, Watson P, Topper TH. A stress-strain function for the fatigue of metals. J Mater 1970;5:767-78.

[25] Wehner T, Fatemi A. Effects of mean stress on fatigue behaviour of a hardened carbon steel. IntJ Fatigue 1991;13(3):241-8.

[26] Nagode M. Continuous damage parameter calculation under thermo-mechanical random loading. MethodsX 2014;1:81-9.

[27] Ottosen NS, Ristinmaa M, The mechanics of constitutive modelling, chapter 4, 2005.

[28] Neuber H. Theory of stress concentration for shear-strained prismatical bodies with arbitrary non-linear stress-strain law. J Appl Mech 1961;28:544-50.

[29] Minano M, Montäns FJ. A new approach to modeling isotropic damage for Mullins effect in hyperelastic materials. Int J Solids Struct 2015;7-68:272-82.

[30] Macha E, Sonsino CM. Energy criteria of multiaxial fatigue failure. Fatigue Fract Eng Mater Struct 1999;22:1053-70.

[31] Liu F, Li H, Li W, Wang B. Experimental study of improved modal strain energy method for damage localisation in jacket-type offshore wind turbines. Renewable Energy 2014;72:174-81.

[32] Oman S, Nagode M, Fajdiga M. The material characterization of the air spring bellow sealing layer. Mater Des 2009;30:1141-50.

[33] Amzallag C, Gerey JP, Robert JL, Bahuaud J. Standardization of the rainflow counting method for fatigue analysis. IntJ Fatigue 1994;16:287-93.

[34] Clormann UH, Seeger T. Rainflow - HCM - Ein Zählverfahren für Betriebsfestigkeitsnachweise auf werkstoffmechanischer Grundlage. Stahlbau 1986;55:65-71.