Step-by-Step Building of a Four Dimensional Fatigue Compatible Regression Model including Frequencies

Enrique Castillo^1,2, Alfonso Fernández-Canteli³, Sergio Blasón^4,5, Golta Khatibi⁶, Bernhard Czerny⁶, Mohammad Zareghomsheh⁶

¹Spanish Royal Academy of Engineering, Madrid, Spain.
²Spanish Royal Academy of Sciences, Madrid, Spain.
³Department of Construction and Engineering Manufacturing, Polytechnic School of Gijón, University of Oviedo, Campus de Viesques, Asturias, Spain.
⁴Bundesanstalt für Materialforschung und Prüfung (BAM), Berlin, Germany.
⁵Department 9-Component Safety, Division 9.4, Weld Mechanics, Berlin, Germany.
⁶Christian Doppler Laboratory for Lifetime and Reliability of Interfaces in Complex Multi-Material Electronics, Chemical Technologies and Analytics, Vienna, Austria.

DOI: 10.4236/ojs.2021.116064                

Abstract

The purpose of this research is to develop a model, with emphasis on compatibility conditions and model building, valid for high cycle fatigue design components such as wind turbines, automobiles, high speed railways and aeronautical material. In this work, we have added the frequency as one more variable to an existing fatigue model that already includes maximum stress, stress ratio and lifetime. As a result, a model and estimation method has been proposed and a random variable V has been identified, which, allows the accumulated damage and the probability of failure to be assessed for any load history in terms of stress levels, stress ranges and frequencies. Finally, the model is validated using a large set of real experimental data.

Keywords

Compatibility Conditions, Cumulative Damage, Frequency Effect, Functional Equations, Regression Model

1. Introduction

The current structural integrity concept in the fatigue design in aeronautics, electronics, railway material, rotating machines and many other fields demands the probabilistic fatigue characterization of the materials and the mechanical and structural components they are made off, so that a reliable safety level is ensured during their service life. With this aim, high performance test machines, from middle-high frequencies (200 to 1000Hz) [1], up to ultrasonic ones [2] [3], have been developed, which accomplish two objectives: First, to reach fatigue lifetimes in the very high cycle fatigue domain, i.e. till 10⁹, or even 10¹⁰ cycles, where new fatigue mechanisms may arise other than those observed in the high cycle fatigue domain. Second, a remarkable reduction of test duration and cost implied when conventional servo-hydraulic or electric fatigue test machines are used. In addition, the correspondence among the fatigue S-N fields for different test frequencies, ranging from low till ultrasonic frequencies must be established and validated if a high quality congruent fatigue design is desired.

Traditionally, most fatigue models have been considered for low frequencies and they are two-dimensional, that is, involving a stress, for example, the stress range, $\Delta \sigma$, and lifetime as the main variables. However, this was proved as not sufficient.

Fatigue has been considered from a long time ago, as one of the most important causes of failure, and being responsible for important accidents in machines, structures and bridges, for example. However, the appearance of models able to analyze the problem in detail has been slow and almost always below the expectations and, mainly, the needs. Today, it is known that a precise analysis of fatigue requires a random model including at least four variables: the stress variables, stress level and stress range, the frequency and the lifetime.

A simple revision of the existing literature shows that there is no unanimity in the selection of these two variables. Common models choose them from a set of six different candidate stress related variables: $\left\{{\sigma }_M\mathrm{,}{\sigma }_m\mathrm{,}\Delta \sigma \mathrm{,}{\sigma }_a\mathrm{,}R\mathrm{,}{\sigma }_{mean}\right\}$, which are the maximum stress, the minimum stress, the stress range, the stress amplitude, the stress ratio, and the mean stress, respectively. However, though in theory different choices are equivalent, as indicated in [4], maximum stress, ${\sigma }_M$ and stress ratio, R, are the most convenient stress variables.

Recently a three dimensional model has been given, see [4], which involves two stresses, the maximum stress ${\sigma }_M$ and the stress ratio R, what implied an important advance in order to predict fatigue for given load histories. In this paper a four dimensional fatigue model including the effects of the maximum stress ${\sigma }_M$, the stress ratio R, and the frequency f on the fatigue lifetime, N, is presented, as a step further. It is built based on the three-dimensional model in [4] to which the frequency, as a fourth and relevant variable, is incorporated. To understand the proposed model in depth and the important role of functional equations in building compatible models, it is necessary to explain in detail how the model arises from scratch. So, as a state of the art, a detailed description of all the steps followed to reach the final model is included. In addition, a functional equation is used later to derive the proposed cumulative damage formula to avoid incompatible formulas causing inconsistencies. Thus, in this paper functional equations are indirectly and strongly recommended for modeling, not only in fatigue but also in problems of general nature.

2. First Step: The S-N Bidimensional Models

Fatigue models started with the up-and-down method, that was designed for determining the endurance limit as its only variable. However, most fatigue models are bidimensional and include only one of the first two stress variables already mentioned and the lifetime, which are the models dealt with in this section. However, one must keep in mind that these models are valid only when the remaining two variables, the second stress and the frequency variables, are considered constant.

A significative advance in the development of fatigue models took place with the appearance of bidimensional S-N regression models, initially linear, later bilinear and finally of hyperbolic type. The selection of the mathematical structure of these models was initially based on convenience and simplicity of the models. This is the case of the linear models, which, when they were perceived as not adequate, a log transformation of one or the two axes was used to improve them. When this was not sufficient, the bilinear models were used as a simple and practical solution.

These regression curves, in the sense of being conditional means, were used extensively, but unfortunately its random character and the corresponding percentile curves were not fully introduced as a necessity to deal with the actual random nature of the problem.

Fortunately, curved trend models of random type were soon accepted and other models appeared. This was the case of the model in [5] where the mathematical structure of the regression curve and percentiles came from the solution of a functional equation expressing a compatibility condition that must be satisfied.

In this context, a special mention deserves the use of functional equations in the process of building models. These equations, as differential or integral equations, express properties that must be satisfied by the models being looked for. However, the unknown functions neither appear as derivatives, nor inside integrals, what makes the interpretation and handling of the problem more convenient, but without losses in its power.

As a reminder, in Figure 1 the compatibility condition used to derive some bidimensional model is illustrated graphically. In it, a $\Delta \sigma -N$ plot is presented, in which the virtual $\Delta \sigma -N$ curves associated with the different specimens used in fatigue tests are shown as continuous lines, in blue color. Since these fatigue tests are normally destructive, only one use of each specimen is possible, and then, only one point on that curve can be observed, but the existence of this curve needs to be emphasized and its conceptual importance to be pointed out.

One first comment on it consists in observing that these blue curves do not

Figure 1. Compatibility conditions: The S-N field showing the compatibility condition implied by the $\Delta \sigma -N$ curves associated with each specimen.

intersect, as it occurs in the plot. This means that one specimen is more resistant than other no matter the value of $\Delta \sigma$, what is a simplifying, but a reasonable assumption.

The first compatibility functional equation, used in [5], is based on the following observation: Any curve in the left plot of Figure 1 entering the square on the left-down rectangular area through out the horizontal line $\Delta \sigma =\Delta {\sigma }_i$ exits this area through out the vertical line $N=N_j$, so that the density functions of the variables lifetime, N, given $\Delta \sigma$, and stress range, $\Delta \sigma$, given N, defined all over the S-N field, permit plotting the density functions, as those shadowed in the figure, representing the fractions of curves, passing through the horizontal and vertical segments, which must be coincident. Note that all this implies that

$F_{\Delta \sigma }\left(\Delta \sigma \mathrm{<mo>\; |\; </mo>}N\right)=F_N\left(N\mathrm{<mo>\; |\; </mo>}\Delta \sigma \right)\mathrm{,}$ (1)

where $F_{\Delta \sigma }\left(\Delta \sigma \mathrm{<mo>\; |\; </mo>}N\right)$ and $F_N\left(N\mathrm{<mo>\; |\; </mo>}\Delta \sigma \right)$ are the cdf of $\Delta \sigma$ given N and the cdf of N given $\Delta \sigma$, respectively. The variability field of both variables in a correct model must fulfil this strong condition. Otherwise, the fatigue model would become not valid.

Next, we assume that the random variables implied in the fatigue behavior definition belong to the same family of distributions, such as one of the extreme value distributions (Weibull, Gumbel, Fréchet) or any other location-scale family. This is a weak and natural assumption, because it would have no sense if for given units the distribution belongs to the family but it fails to belong when changing units of measure. Similarly, when asymptotes are involved it is natural to admit location changes so that changes in the asymptotes are included in the selected family of distributions.

If the cdfs involved in (1) are assumed strictly increasing, the functional equation

${\left(\frac{\Delta \sigma -{\lambda }_1\left(N\right)}{{\delta }_1\left(N\right)}\right)}^{{\beta }_1\left(N\right)}={\left(\frac{N-{\lambda }_2\left(\Delta \sigma \right)}{{\delta }_2\left(\Delta \sigma \right)}\right)}^{{\beta }_2\left(\Delta \sigma \right)},$ (2)

with six unknown functions: ${\lambda }_1\left(N\right)$, ${\delta }_1\left(N\right)$, ${\beta }_1\left(N\right)$, ${\lambda }_2\left(\Delta \sigma \right)$, ${\delta }_2\left(\Delta \sigma \right)$ and ${\beta }_2\left(\Delta \sigma \right)$ results.

A very general solution of this equation was found in [5], but a rigorous solution of it was given in [6] and [7], in which it was shown that there is only one family of physically percentile compatible solutions, with $B\mathrm{,}C\mathrm{,}D$ and $E$ constant parameters, of the form:

$H\left(\Delta \sigma \mathrm{,}N\right)=\frac{\left(N-B\right)\left(\Delta \sigma -C\right)+E}{D}=\text{constant}\mathrm{.}$ (3)

Next, assuming validity of the weakest link principle for the lifetimes, the fatigue lifetime distribution represents a minimum problem asymptotically defined by the generalized extreme value family, see [8] [9] [10] [11] or [12]. Because the existence of a lower bound in the lifetime (it cannot be less than zero) the Weibull distribution becomes the only candidate, see [13] [14] and [15]. The Weibull distribution together with the percentile curves of the form in (3) complete the bidimensional model.

Consequently, these models were selected as the only appropriate ones for fatigue analysis when ${\sigma }_M$ is kept constant and when R is kept constant during the test, see [12] and [16].

Since most bidimensional fatigue models used in practice are written in terms of the lifetimes supplied for the stress range, for a fix value of the stress level, selected from the set of secondary variables $\left\{{\sigma }_M\mathrm{,}{\sigma }_m\mathrm{,}{\sigma }_{mean}\mathrm{,}R\right\}$, to reduce the effect of this limitation, a hypothetical equivalence among the different S-N fields for different stress levels is assumed, even though this is only a strong simplification. Some examples of this approach are the models of [17] - [25]. Some comparative studies among these models were reported in [26] - [32]. However, only very few works address the probabilistic analysis of models focused on the stress level effect, [33] and [34], or multiple explanatory variables, [35].

Unfortunately, these empirical type proposals are supported only by their simplicity and the associated fitting quality of the associated experimental results, but with no theoretical justification. A revision of some of these methods used in fatigue analysis can be seen in [12], together with some new and original proposals for fatigue models dealing with crack growth and damage accumulation.

The final conclusion is that bidimensional models are not satisfactory for fatigue analysis, but the compatibility condition must hold for all bidimensional models coming from higher dimension models when the remaining variables are kept constant.

Since the same considerations used for the $\Delta \sigma -N$ models for constant ${\sigma }_M$ are also valid for constant R, this implies that we have two possible compatible bidimensional models to be used:

$H\left(\Delta \sigma \mathrm{,}N\mathrm{<mo>\; ;\; </mo>}{\sigma }_M\right)=\frac{\left(N-B_M\left({\sigma }_M\right)\right)\left(\Delta \sigma -C_M\left({\sigma }_M\right)\right)+E_M\left({\sigma }_M\right)}{D_M\left({\sigma }_M\right)}$ (4)

and

$H\left(\Delta \sigma \mathrm{,}N\mathrm{<mo>\; ;\; </mo>}R\right)=\frac{\left(N-B_R\left(R\right)\right)\left(\Delta \sigma -C_R\left(R\right)\right)+E_R\left(R\right)}{D_R\left(R\right)}\mathrm{.}$ (5)

Even though these two similar models can be used, if not carefully chosen they will be incompatible. In the following section the compatibility problem of both bidimensional models is stated and solved. This permits identifying the only compatible three-dimensional models of the form $H\left({\sigma }_M\mathrm{,}R\mathrm{,}N\right)$.

3. Second Step: A Three-Dimensional Compatibility Model

In [36] and [37] the two bi-dimensional models in (4) and (5) were made consistent using a functional equation leading to the three-dimensional model ${\sigma }_M-R-N$.

Figure 2 shows some percentile curves of models (4) and (5) for two ${\sigma }_M$ constant values, ${\sigma }_{M_1}$ and ${\sigma }_{M_2}$, as continuous lines, and the cases of tests for two R constant values, $R_1$ and $R_2$, as dashed lines. The intersections of the percentile lines represent four tests, indicated by thick segments, where each one belongs to two sets of the previous four test types and then, they must provide exactly the same results. The horizontal intersections of the percentiles are not a coincidence, but a strong cross compatibility condition that need to be satisfied. This condition leads to the following functional equation, which is obtained from Equations (4) and (5) by eliminating $\Delta \sigma$ and replacing it by ${\sigma }_M\left(1-R\right)$, that is:

Figure 2. The cross compatibility of the probabilistic distributions arising from two bidimensional S-N fields referred to constant ${\sigma }_M$ and constant R.

$\begin{array}{c}H\left(N,{\sigma }_M,R\right)=\frac{\left(N-B_M\left({\sigma }_M\right)\right)\left({\sigma }_M\left(1-R\right)-C_M\left({\sigma }_M\right)\right)+E_M\left({\sigma }_M\right)}{D_M\left({\sigma }_M\right)}\\ =\frac{\left(N-B_R\left(R\right)\right)\left({\sigma }_M\left(1-R\right)-C_R\left(R\right)\right)+E_R\left(R\right)}{D_R\left(R\right)},\end{array}$ (6)

which depends on five unknown functions $B_M\left({\sigma }_M\right)$, $C_M\left({\sigma }_M\right)$, $D_M\left({\sigma }_M\right)$, $E_M\left({\sigma }_M\right)$, $B_R\left(R\right)$, $C_R\left(R\right)$, $D_R\left(R\right)$ and $E_R\left(R\right)$. Even though its non-consideration implies incompatibility, most new fatigue models proposed in the last thirty years do not pay attention to this property, as indicated in [38].

4. Third Step: Asymptotes Compatibility

In addition, some of the functions appearing in (6), which define model asymptotes, as indicated in [4], are known. In fact, we must have:

$C_M\left({\sigma }_M\right)=0;C_R\left(R\right)=\left(1-R\right){\sigma }_{M0}.$ (7)

The functional Equation (6) and the asymptotes in (7) lead to the following Weibull model:

$\boxed{p=1-\exp \left\{-{\left[\frac{\left({\log }_e\frac{N}{N_0}-{\delta }_1\frac{R}{1-R}\right)\left(1-R\right)\left(\frac{{\sigma }_M}{{\sigma }_{M0}}-1\right)-{\beta }_1\left(1-R\right)-{\gamma }_1}{a}\right]}^K\right\}}$ (8)

which depends on seven parameters: $a,{\sigma }_{M0},N_0,{\beta }_1,{\delta }_1,{\gamma }_1$ and K.

From (8), the following V Weibull random variable, which will play a key role in the cumulative damage evaluation, arises:

$\boxed{V\left(N;{\sigma }_M,R\right)=\left({\log }_e\frac{N}{N_0}-{\delta }_1\frac{R}{1-R}\right)\left(1-R\right)\left(\frac{{\sigma }_M}{{\sigma }_{M0}}-1\right)-{\beta }_1\left(1-R\right)}$ (9)

The physical meaning of these parameters is as follows: a is the scale parameter of the V Weibull distribution; ${\sigma }_{M0}$ is the endurance limit for ${\sigma }_M$ ; $N_0$ is the N asymptote for $R=0$ ; ${\beta }_1$ is the increment of the Weibull location parameter of the V variable for $R=0$ with respect to $R=1$ ; ${\gamma }_1$ is the Weibull location parameter of the V variable; ${\delta }_1$ is the displacement of the ${\log }_e{N}^{\mathrm{<mo>\; *\; </mo>}}$ asymptote for $R=1/2$ with respect to $R=0$ and K is the Weibull shape parameter of the V variable.

In addition, the associated regression models are:

${\sigma }_M={\sigma }_{M0}\left\{1+\frac{{\beta }_1\left(1-R\right)+h\left({\gamma }_1,a,p,K\right)}{\left({\log }_e\frac{N}{N_0}-{\delta }_1\frac{R}{1-R}\right)\left(1-R\right)}\right\};N\ge N_0\exp \left(\frac{{\delta }_{1}R}{1-R}\right),$ (10)

where

$h\left({\gamma }_1\mathrm{,}a\mathrm{,}p\mathrm{,}K\right)={\gamma }_1+a{\left(-{\log }_e\left(1-p\right)\right)}^{1/K}$ (11)

is the random part of the model, and letting $\Delta \sigma ={\sigma }_M\left(1-R\right)$ one gets:

$\Delta \sigma ={\sigma }_{M0}\left\{1-R+\frac{{\beta }_1\left(1-R\right)+h\left({\gamma }_1,a,p,K\right)}{\left({\log }_e\frac{N}{N_0}-{\delta }_1\frac{R}{1-R}\right)}\right\};N\ge N_0\exp \left(\frac{{\delta }_{1}R}{1-R}\right),$ (12)

which allows us to plot ${\sigma }_M$ and $\Delta \sigma$ versus ${\log }_e\frac{N}{N_0}$ for constant R values, and

$\Delta \sigma ={\sigma }_M\frac{{\delta }_1\left(\frac{{\sigma }_M}{{\sigma }_{M_0}}-1\right)+h\left({\gamma }_1,a,p,K\right)}{\left({\log }_e\frac{N}{N_0}+{\delta }_1\right)\left(\frac{{\sigma }_M}{{\sigma }_{M_0}}-1\right)-{\beta }_1};N\ge N_0\exp \left(\frac{{\beta }_1}{\frac{{\sigma }_M}{{\sigma }_{M_0}}-1}-{\delta }_1\right),$ (13)

which allows us to plot $\Delta \sigma$ versus ${\log }_e\frac{N}{N_0}$ for constant ${\sigma }_M$ values.

The asymptotes of the curves in (10), (12) and (13) are given in Table 1. The compatibility of the inequality relations in the plots in Figure 2, which must be taken into account for the model to be valid, implies the two constraints ${\delta }_1\ge 0$ and ${\beta }_1\ge 0$.

Next, a physical interpretation of the V variable defined in Expression (9) is given.

The V variable, which has arisen in a natural form, when evaluating the probability of failure of this model and imposing the compatibility constraints, can be interpreted as damage level associated with an effective load, the pair $\left({\sigma }_M\mathrm{,}R\right)$ and the lifetime N. In other words, it can be associated with the damage caused to the material when subjected to N cycles under a load $\left({\sigma }_M\mathrm{,}R\right)$ for constant frequency. Consequently, using the V variable in (9), a damage comparison can be made of two specimens subject to different load histories, when the frequency is constant. This is relevant, because later, based on this variable, the damage accumulation associated with any varying load history, will be calculated.

Table 1. Horizontal and vertical asymptotes of the curves in (10), (12) and (13).

5. Fourth Step: Reproducing the Frequency Effect

In this section the basic ideas and formulas used to reproduce the effect of frequency on fatigue tests are discussed. It is assumed that the influence of frequency on fatigue is manifested in the variation of the elastic modulus tending to an effective elastic modulus, $E_{ef}$. The observation is consistent with experimental experience that the speed in a static test influences the apparent modulus of elasticity.

The higher the test speed, the more rigid the material behaves. Therefore, the higher the frequency, the higher the effective elastic modulus, $E_{ef}$. Asymptotically, when the frequency is practically zero, the behavior of the material, i.e. the effective modulus of elasticity, is the nominal one, so $E_{ef}=E_0$. When the frequency is very high, the effective modulus of elasticity tends to an asymptotic value, $E_{\infty }$.

The effect of frequency on the S-N fields has been studied in [39], where a model based on a three-parameter Weibull cumulative distribution was proposed to analyze the influence of frequency and stress ratio on the fatigue life in concrete, both plain and reinforced with fibers, showing the reduction of effective stresses due to increased frequencies.

All this is in accordance with a viscoelastic simile of the material, in the sense that for high values of frequency, the behavior of the material tends to reach an asymptotic value according to a sigmoidal function, such as a biparametric Weibull cdf. If we observe a cycle in fatigue under the condition of test in load control, this assumes that the induced deformation in a load cycle decreases as the frequency increases, as it can be seen in Figure 3.

Since deformation is the determining mechanism in the fatigue process, the range of deformations deduced from the range of stresses, assuming that the static modulus of elasticity is the key factor. This implies that in reality, a lower deformation is being applied according to the relationship $E_0/E_{ef}$ and, accepting the linear relationship between strain and stress given by $E_0$, for the representation in the $S-N$ field, a shift to the right of the $S-N$ curves will be observed as the frequency increases.

Figure 3. Illustration of the role played by the Young modulus on the fatigue frequency effect and the stress ${\sigma }_M$ and ${\sigma }_m$ reduction due to frequency.

According to this, we can normalize the Young modulus and, since the normalized range becomes the closed interval $\left[\mathrm{0,1}\right]$ we can assume a Weibull cdf for the logarithm of the frequency, as follows:

$\bar{E}=\frac{E_{ef}-E_0}{E_{\infty }-E_0}=1-\exp \left[-{\left(\frac{{\log }_e\left(f/f_0\right)}{\tau }\right)}^{\beta }\right],$ (14)

where $f_0$ is a reference frequency and $\tau$ is a dimensionless scale parameter. Note that the Weibull distribution has been justified for the logarithm of $f/f_0$, a dimensionless variable, as it must always be for logarithm and trigonometric functions, a condition frequently ignored. Then, we have

$E_{ef}=E_0+\left(1-\exp \left[-{\left(\frac{{\log }_e\left(f/f_0\right)}{\tau }\right)}^{\beta }\right]\right)\left(E_{\infty }-E_0\right).$ (15)

From now on, we call

$a_{ef}=\frac{E_0}{E_{ef}}\le 1.$ (16)

which, from (15) gives

$a_{ef}=\frac{1}{1+\left(1-\exp \left[-{\left(\frac{{\log }_e\left(f/f_0\right)}{\tau }\right)}^{\beta }\right]\right)\left(E_{ratio}-1\right)}$ (17)

where

$E_{ratio}=\frac{E_{\infty }}{E_0}\ge 1.$ (18)

Accordingly, the value of ${\sigma }_M$ and the stress amplitude are affected by the factor $a_{ef}$, while the mean stress is not affected. This is due to the way of carrying out the test, in which the load is placed in the medium stress and then it is expanded until it reaches the desired amplitude. With which, it turns out that:

${\sigma }_{M,ef}={\sigma }_{mean,ef}+{\sigma }_{a,ef}={\sigma }_{mean,nom}+a_{ef}{\sigma }_{a,nom}$ (19)

$=\frac{{\sigma }_{M,nom}+{\sigma }_{m,nom}}{2}+a_{ef}\frac{{\sigma }_{M,nom}-{\sigma }_{m,nom}}{2}$ (20)

$=\frac{\left(1+a_{ef}\right){\sigma }_{M,nom}+\left(1-a_{ef}\right){\sigma }_{m,nom}}{2}$ (21)

$=\frac{\left[\left(1+a_{ef}\right)+\left(1-a_{ef}\right)R_{nom}\right]{\sigma }_{M,nom}}{2}$ (22)

and then,

$R_{ef}=\frac{{\sigma }_{m,ef}}{{\sigma }_{M,ef}}=\frac{{\sigma }_{mean,nom}-{\sigma }_{a.ef}}{{\sigma }_{mean,nom}+{\sigma }_{a.ef}}$ (23)

$=\frac{{\sigma }_{mean,nom}-a_{ef}{\sigma }_{a.nom}}{{\sigma }_{mean,nom}+a_{ef}{\sigma }_{a.nom}}=\frac{\frac{{\sigma }_M+{\sigma }_m}{2}-a_{ef}\frac{{\sigma }_M-{\sigma }_m}{2}}{\frac{{\sigma }_M+{\sigma }_m}{2}+a_{ef}\frac{{\sigma }_M-{\sigma }_m}{2}}$ (24)

$=\frac{\left(1-a_{ef}\right){\sigma }_M+\left(1+a_{ef}\right){\sigma }_m}{\left(1+a_{ef}\right){\sigma }_M+\left(1-a_{ef}\right){\sigma }_m}=\frac{\left(1-a_{ef}\right)+\left(1+a_{ef}\right)R_{nom}}{\left(1+a_{ef}\right)+\left(1-a_{ef}\right)R_{nom}}.$ (25)

Referring now to the nominal values, those intended in the tests, the following efective ${\sigma }_{M\mathrm{,}ef}$ and $R_{ef}$ values must replace the nominal ones, ${\sigma }_{M\mathrm{,}nom}$ and $R_{nom}$ in (10), (12) and (13):

${\sigma }_{M,ef}=\frac{\left[\left(1+a_{ef}\right)+\left(1-a_{ef}\right)R_{nom}\right]{\sigma }_{M,nom}}{2};$ (26)

$R_{ef}=\frac{\left(1-a_{ef}\right)+\left(1+a_{ef}\right)R_{nom}}{\left(1+a_{ef}\right)+\left(1-a_{ef}\right)R_{nom}}$ (27)

or their inverses:

${\sigma }_{M,nom}=\frac{\left(\left(1+a_{ef}\right)-\left(1-a_{ef}\right)R_{ef}\right){\sigma }_{M,ef}}{2a_{ef}};$ (28)

$R_{nom}=\frac{\left(1-a_{ef}\right)-\left(1+a_{ef}\right)R_{ef}}{\left(1-a_{ef}\right)R_{ef}-\left(1+a_{ef}\right)}.$ (29)

This means that the influence of the frequency on fatigue is equivalent to a change in the ${\sigma }_M$ and R values by their equivalent nominal ones ${\sigma }_{M\mathrm{,}ef}$ and $R_{ef}$ and that there exists a unique fatigue limit, independently of the frequency and R.

If the model turns out to be correct, the frequency characterization of a material, that is, the determination of the curve $\varphi =\varphi \left(f\right)$ would be easily determinable. It would only be a matter of applying a relatively high stress at different frequencies and measuring the deformation due to that stress, which would allow us to observe the evolution of the modulus of elasticity as a function of frequency.

6. Proposed Model

The model proposed in this paper starts from the model developed in [4], in which a three-dimensional compatible ${\sigma }_M-R-N$ model for constant frequency was derived using functional equations. In this section this model is extended to four dimensions by including the frequency variable, which plays a very important role in several applications, e. g. high speed train bridge safety.

It is clear that this model and consequently, the following equations must be applied to the effective values of ${\sigma }_M$ and R, because, even though the material is apparently subject to given nominal values of ${\sigma }_M$ and R, the frequencies modify these actual values to the effective ones. Thus, the following equations refer to effective values.

In the previous model it was demonstrated that the variable in (9) plays a crucial role, because all the remaining variables can be considered as deterministically dependent on this one. In other words, the model can be written in a form such that its deterministic and its random parts are identified, that is, they become separable. In addition, since the variable V is unidimensional, only the random parameters associated with it rule the random behavior. Obviously, the remaining model parameters are deterministic.

Replacing ${\sigma }_M,R$ and $\Delta \sigma$ by ${\sigma }_{M_{ef}},R_{ef}$ and $\Delta {\sigma }_{ef}$, respectively, in (10), (12) and (13) one gets:

1) ${\sigma }_{M_{ef}}$ for given $R_{ef}$ :

${\sigma }_{M_{ef}}={\sigma }_{M0}\left\{1+\frac{{\beta }_1\left(1-R_{ef}\right)+h\left({\gamma }_1\mathrm{,}a\mathrm{,}p\mathrm{,}K\right)}{\left({\log }_e\frac{N}{N_0}-{\delta }_1\frac{R_{ef}}{1\text{}-\text{}{R}_{ef}}\right)\left(1-v{R}_{ef}\right)}\right\}\mathrm{<mo>\; ;\; </mo>}N\ge N_0\exp \left(\frac{{\delta }_1{R}_{ef}}{1-R_{ef}}\right)\mathrm{.}$ (30)

2) $\Delta {\sigma }_{ef}={\sigma }_{M_{ef}}\left(1-R_{ef}\right)$ for given $R_{ef}$ :

$\Delta {\sigma }_{ef}={\sigma }_{M0}\left\{1-R_{ef}+\frac{{\beta }_1\left(1-R_{ef}\right)+h\left({\gamma }_1,a,p,K\right)}{\left({\log }_e\frac{N}{N_0}-{\delta }_1\frac{R_{ef}}{1-R_{ef}}\right)}\right\};N\ge N_0\exp \left(\frac{{\delta }_1{R}_{ef}}{1-R_{ef}}\right),$ (31)

3) $\Delta {\sigma }_{ef}$ for given ${\sigma }_{M_{ef}}$ :

$\Delta {\sigma }_{ef}={\sigma }_{M,ef}\frac{{\delta }_1\left(\frac{{\sigma }_{M_{ef}}}{{\sigma }_{M_0}}-1\right)+h\left({\gamma }_1,a,p,K\right)}{\left({\log }_e\frac{N}{N_0}+{\delta }_1\right)\left(\frac{{\sigma }_{M_{ef}}}{{\sigma }_{M_0}}-1\right)-{\beta }_1};N\ge N_0\exp \left(\frac{{\beta }_1}{\frac{{\sigma }_{M_{ef}}}{{\sigma }_{M_0}}-1}-{\delta }_1\right),$ (32)

where ${\delta }_1\mathrm{,}{\beta }_1\ge 0$ and ${\sigma }_{M_{ef}}>{\sigma }_{M0}$ for the model to be valid.

These expressions, allow us to plot ${\sigma }_{M_{ef}}$ and $\Delta {\sigma }_{ef}$ versus ${\log }_e\frac{N}{N_0}$ for constant $R_{ef}$ values and $\Delta {\sigma }_{ef}$ versus ${\log }_e\frac{N}{N_0}$ for constant ${\sigma }_{M_{ef}}$ values (some examples can be seen in the same paper [4] ).

It is important to realize that the common term $h\left(\gamma \mathrm{,}a\mathrm{,}p\mathrm{,}K\right)$ in (11) appearing in expressions (30) (31) and (32) refers to the Weibull distribution of the V in Expression (9), where $\gamma$, a and K are its location, scale and shape parameters, respectively. Alternatively, other family of distributions can be selected to reproduce the random behavior of V.

The asymptotes of the curves in (30), (31) and (32) are given in Table 2.

Note finally, that the role of the V damage variable remains valid, when frequencies are considered if the effective stresses are used instead of the nominal ones.

7. Parameter Estimation

To estimate the ten model parameters, $a,N_0,{\sigma }_{M_0},{\beta }_1,{\delta }_1,\tau ,{\gamma }_1,f_0,\beta$ and K, a two-step method is used. In the first step the seven deterministic parameters, $\left\{N_0\mathrm{,}{\sigma }_{M0}\mathrm{,}{\delta }_1\mathrm{,}\tau \mathrm{,}{\beta }_1,f_0\mathrm{,}\beta \right\}$, are estimated using a regression estimate based on

Table 2. Horizontal and vertical asymptotes of the curves in (30), (31) and (32).

least squares, which minimizes the sum of error terms, that is, the differences between the observed values, ${\sigma }_{M_{f_i}}$, and those given by the regression curves of Expression (10), which correspond to the V mean value $\bar{V}={\gamma }_1+a\Gamma \left(1+\frac{1}{K}\right)$. In other words, the following optimization problem is solved:

$\underset{N_0\mathrm{,}{\sigma }_{M0}\mathrm{,}{\delta }_1\mathrm{,}\tau \mathrm{,}{\beta }_1\mathrm{,}\bar{V}\mathrm{,}{f}_0\mathrm{,}\beta }{\text{Minimize}}Z={\displaystyle \underset{i=1}{\overset{n}{\sum }}}{\left[{\sigma }_{M_{e{f}_i}}-{\sigma }_{M0}\left(1+\frac{{\beta }_1\left(1-R_{e{f}_i}\right)+\hat{V}}{\text{}\text{}\left({\log }_e\frac{N_i}{N_0}-{\delta }_1\frac{R_{e{f}_i}}{1-R_{e{f}_i}}\right)\left(1-R_{e{f}_i}\right)}\right)\right]}^2$ (33)

where $M_{e{f}_i}$ and $R_{e{f}_i}$ depend on $f_0\mathrm{,}\tau$ and $\beta$, subject to

$N_i\ge N_0\exp \left(\frac{{\delta }_1{R}_{e{f}_i}}{1-R_{e{f}_i}}\right)\mathrm{,}$ (34)

${\delta }_1\mathrm{,}{\beta }_1\mathrm{,}{N}_0\ge \mathrm{0,}$ (35)

${\sigma }_{M_{el{l}_i}}\ge {\sigma }_{M0}\mathrm{,}$ (36)

where n is the sample size.

In the second step, and keeping $\bar{V}=\gamma +a\Gamma \left(1+\frac{1}{K}\right)$, that is, the mean of V,

constant, the random parameters, $\left\{\gamma \mathrm{,}a\mathrm{,}K\right\}$ are obtained by restricted maximum likelihood. This restriction guaranties that the previously estimated regression curves are not modified and also avoids overfitting. This means that the following optimization problem is solved:

$\underset{\gamma \mathrm{,}a\mathrm{,}K}{\text{Minimize}}Z={\displaystyle \underset{i=1}{\overset{n}{\sum }}}\left(-{\log }_{e}K+{\log }_{e}a+W_i^K-\left(K-1\right){\log }_e{W}_i\right)$ (37)

subject to $K\ge 0$, where the expression between parentheses in (37) is the log-likelihood, changed sign, and

$W_i=\frac{v_i-\gamma }{a};i=1,2,\cdots ,n.$ (38)

This optimization problem has been implemented in matlab, using the function “fmincon”. One example of application is given in Section 9.

8. Damage Accumulation

The damage interpretation of the V variable in (9) offers us the possibility of deriving from it a damage accumulation conversion law, which will allow us to calculate the probability of failure derived from any varying load history.

It has already been indicated that the V expression in (9) gives the cumulative damage associated with N cycles when the load ${\sigma }_{M_{ef}}$ and $R_{ef}$ is applied assuming that the initial damage is zero. In addition, its cdf gives the probability of failure, associated with N cycles when that load is applied assuming a threshold value, $V_0$, associated with a threshold number of load cycles, for which the probability of failure is zero. Thus, the damage level associated with N load cycles and zero initial damage, is:

$V\left(N\mathrm{,}{\sigma }_{M_{ef}}\mathrm{,}{R}_{ef}\right)=Q\left(\mathrm{0,}N\mathrm{,}{\sigma }_{M_{ef}}\mathrm{,}{R}_{ef}\right)\mathrm{.}$ (39)

However, we look for a more general formula

$V=Q\left(V_0,N,{\sigma }_{M_{ef}},R_{ef}\right),$ (40)

allowing us to consider how the damage level evolves after a given damage, $V_0$, is reached. This is one of the aims of this section.

We start by saying that the Expression in (40) cannot be given arbitrarily. Otherwise, using different threshold damages will lead to different final damage levels even though the loads and the number of cycles coincide (see, for example, [40] ).

To avoid this problem, the formula to obtain V in (40) must be subjected to a compatibility condition, which states that the same damage level must be obtained if (a) we start from a $V_0$ threshold damage and apply $N_0+N$ cycles of the load or (b) we start from the damage level associated with $N_0$ cycles and afterwards we apply N cycles of the same load. This, compatibility condition can be written as:

$Q\left(V_0\mathrm{,}{N}_0+N\mathrm{,}{\sigma }_{M_{ef}}\mathrm{,}{R}_{ef}\right)=Q\left(Q\left(V_0\mathrm{,}{N}_0\mathrm{,}{\sigma }_{M_{ef}}\mathrm{,}{R}_{ef}\right)\mathrm{,}N\mathrm{,}{\sigma }_{M_{ef}}\mathrm{,}{R}_{ef}\right)\mathrm{,}$ (41)

which is a functional equation, with its only unknown function being $Q\left(V_0\mathrm{,}N\mathrm{,}{\sigma }_{M_{ef}}\mathrm{,}{R}_{ef}\right)$ (see [41] [42] [43] [44] or [45] ). From (41) with $V_0=0$, we get

$Q\left(\mathrm{0,}{N}_0+N\mathrm{,}{\sigma }_{M_{ef}}\mathrm{,}{R}_{ef}\right)=Q\left(Q\left(\mathrm{0,}{N}_0\mathrm{,}{\sigma }_{M_{ef}}\mathrm{,}{R}_{ef}\right)\mathrm{,}N\mathrm{,}{\sigma }_{M_{ef}}\mathrm{,}{R}_{ef}\right)\mathrm{.}$ (42)

and using (39) twice one obtains:

$V\left(N_0+N\mathrm{,}{\sigma }_{M_{ef}}\mathrm{,}{R}_{ef}\right)=Q\left(V\left(N_0\mathrm{,}{\sigma }_{M_{ef}}\mathrm{,}{R}_{ef}\right)\mathrm{,}N\mathrm{,}{\sigma }_{M_{ef}}\mathrm{,}{R}_{ef}\right)\mathrm{.}$ (43)

and letting $V\left(N_0,{\sigma }_{M_{ef}},R_{ef}\right)=V_0$ in (43) the solution of Equation (41) is obtained:

$Q\left(V_0,N,{\sigma }_{M_{ef}},R_{ef}\right)=V\left(V^{-1}\left(V_0,{\sigma }_{M_{ef}},R_{ef}\right)+N,{\sigma }_{M_{ef}},R_{ef}\right),$ (44)

where the inverse of the V function refers to its first argument and it has been assumed that it exists.

If the load history is defined as a set of loads $\left\{\left({\sigma }_{M_{ef}}\left(j\right),R_{ef}\left(j\right)\right);j=1,2,\cdots ,m\right\}$, for given frequencies, repeated for $N\left(j\right)$ cycles, the partial and the final damage states can be calculated as follows:

$\begin{array}{l}{V}_0\left(j\right)=Q\left(V_0\left(j-1\right),N\left(j\right),{\sigma }_{M_{ef}}\left(j\right),R_{ef}\left(j\right)\right)\\ =V\left(V^{-1}\left(V_0\left(j-1\right),{\sigma }_{M_{ef}}\left(j\right),R_{ef}\left(j\right)\right)+N\left(j\right),{\sigma }_{M_{ef}}\left(j\right),R_{ef}\left(j\right)\right);j=1,2,\cdots ,m\end{array}$ (45)

where $V_0\left(j-1\right)$ is the threshold damage for load $\left\{{\sigma }_{M_{ef}}\left(j\right)\mathrm{,}{R}_{ef}\left(j\right)\right\}$, and the initial and final damage levels are $V_0=V_0\left(0\right)$ and $V_0\left(m\right)$, respectively.

It is obvious that these two formulas become crucial when one needs to evaluate damage states caused by different load histories, as it happens in real practice.

9. Example of Application

In this example we have a set of 243 fatigue test data. The samples consisted of commercial high purity Al wires (type Al-H11 from Heraeus in as-received condition with a diameter of 400 μm and a gauge length of 10 mm. A TA ElectroForce DMA 3200 testing machine with a load capacity of 500 N, a frequency range of 0.00001 - 300 Hz and a displacement range of ±6 mm was used for the experiments. The load controlled fatigue tests were carried out in tension-tension mode at a stress ratio of $R\approx 0.1$ by using 6 - 10 samples at each stress level. The S-N curves were obtained by collecting the data in a broad range of 10² up to about 10⁸ loading cycles at three testing frequencies of 2, 20 and 200 Hz, and the model parameters have been estimated using the methods described in Section 7 and commented below. E_ratio is assumed equal 3.

The resulting lifetimes are plotted in Figure 4, using three different colors, blue, for $f=2\text{Hz}$, green, for $f=20\text{Hz}$ and red, for $f=200\text{Hz}$, respectively.

A simple look to the data in Figure 4 reveals the presence of a few outliers, which are indicated using squares instead of circles for the data points. We

Figure 4. Plot of ${\sigma }_M$ versus lifetime N of the laboratory experimental data, with all data and outliers in the left plot and with removed outliers in the right plot.

advance that these five outliers, which are only 2% of the 243 sample data points, were selected using both, probability papers and analyzing to percentile curves of the fitted model including all data points (see Figures 4-8).

Table 3 shows the five outliers indicating their sample indices, their lifetimes and the corresponding values of ${\sigma }_M$ in MPa.

The two-step method described in detail in Section 6 has been used to estimate the parameters with all data points and the resulting fitted model plotted to analyze the possibility of some outliers, as revealed by our first look to the data points. The outliers were confirmed based on Figures 4-8, but this can also be done using well known methods, such as, for example, those given in [46] and [47].

Once this has been done and the outliers identified, they have been removed from the sample and the estimation and possible detection of more outliers repeated to select the final model for its use in practical applications, such as an analysis of safety against fatigue in high speed railway bridges.

For the sake of facilitating the comparison of both results, with and without outliers, and to analyze how them affect the final results, the plots corresponding to both estimation processes have been grouped, such that the left plots refer to the whole data and the right plots to the removed outliers case, respectively.

The first estimation step, in which the regression curves are obtained using least squares, provides the seven deterministic parameter estimates and the second step gives the three random ones, obtained by restricted maximum likelihood, assuming a Weibull distribution for the residuals. They are shown on the left and right parts of Table 4, respectively.

Figure 5 shows the data sample points and the ${\sigma }_M-N$ and ${\sigma }_{M_{ef}}-N$ resulting regression curves for the three different frequencies: 2, 20 and 200 Hz, which show a good fit.

Figure 6 shows the V values on a reverse Gumbel probability paper, revealing

Table 3. Outliers revealed after analyzing the information in Figures 4-8.

Table 4. Deterministic and random parameter estimates obtained with all data points and after removing the five outliers.

Figure 5. Data, ${\sigma }_M-N$ and ${\sigma }_{M_{eff}}-N$ regression curves for three different frequencies: 2, 20 and 200 Hz, fitted by least squares, using all data, on the left plots, and after removing the five outliers, on the right plots.

Figure 6. V resulting sample data plotted on a reverse Gumbel probability paper, corresponding to the whole sample and after removing the outliers.

a close to linear general trend but showing some outliers at both ends. As already indicated, the five outliers have been selected using this plot and the following percentile plots. They clearly show that we are in front of outliers, some of them confirming our first impression after looking to the data in Figure 4.

The right plot in Figure 6 corresponds to the data points after removing the outliers, which also shows a linear trend and suggests no more outliers present. The lower end, suggesting a vertical asymptote, indicates a Weibull distribution close to a Gumbel one ( $K=6.66$ ) (see [48] ).

The plots in Figure 7 show the empirical cdf of all the data and all data with outliers removed, together with the 0.01, 0.05, 0.50, 0.95 and 0.99 percentiles of the order statistics of a Weibull model, confirming once more the Weibull family as a reasonable choice.

It is clear that the Weibull distribution of V is justified, theoretically and physically, by extreme value theory, because fatigue failure is a process that takes place using the weakest link or weakest region principle, stating that failures occur at the weakest locations.

Alternatively, in the second step, the errors or residuals, $V_i$, are also calculated and plotted, in Figure 8, together with the “ksdensity” matlab kernel estimate, which is based on a normal kernel function and uses a window parameter (bandwidth) as explained in [49].

The plots in Figure 8 show the Weibull and kernel pdf estimates together with the data for the whole data and the removed outliers on the top and bottom plots, respectively, and the Weibull and kernel cases, on the left and right plots,

Figure 7. Weibull probability bands for the order statistics of the sample considering all data and with no outliers.

Figure 8. Weibull and non-parametric errors density estimates, left and right, respectively, with all sample data and after removing the outliers, upper and lower, respectively, together with the experimental data.

respectively. It is interesting to see how the top plots allow identifying the outliers on their left and right zones.

If the Weibull model is used, the resulting regression curves and the associated percentiles are shown in Figure 9 and Figure 10, for ${\sigma }_M$ and ${\sigma }_{M_{ef}}$, respectively.

Figure 9 shows the three data sets, in different colors, and the resulting 0.01, 0.05, 0.50, 0.95 and 0.99 percentiles of the ${\sigma }_M$ versus lifetime N regression curves using a Weibull parametric estimate, for outliers included (left) and removed (right) cases, respectively.

Finally, Figure 10 shows the three data sets, in different colors, and the resulting 0.01, 0.05, 0.50, 0.95 and 0.99 percentiles of the ${\sigma }_{M_{ef}}$ versus lifetime N regression curves using the Weibull parametric estimate. Note that the data are in accordance with the model showing a good fit and then, confirming the

Figure 9. The 0.01, 0.05, 0.50, 0.95 and 0.99 percentiles of the ${\sigma }_M$ versus lifetime N regression curves, for three frequencies of 2, 20 and 200 Hz and the resulting Weibull model, with (left) and without (right) outliers.

Figure 10. The 0.01, 0.05, 0.50, 0.95 and 0.99 percentiles of the ${\sigma }_{M_{ef}}$ versus lifetime N regression curves, for three frequencies of 2, 20 and 200 Hz and the resulting Weibull model, with (left) and without (right) outliers.

selected model. As it can be seen, this also reveals the presence of some outliers, corresponding to the two sides of the regression curve.

Note also that the outliers also appear in Figure 9 and Figure 10 on both sides of the mean, as already indicated by the density estimate in Figure 8, where the five outliers correspond to the five points on the left and right extremes.

The regression curves ${\sigma }_M$ versus N and ${\sigma }_{M_{ef}}$ versus N corresponding to the kernel estimates together with the associated percentiles, as those in Figure 9 and Figure 10, are omitted because they are practically identical to the Weibull model ones. This is due to the fact that the densities in Figure 8 are very similar for the Weibull and kernel cases.

Finally, the resulting diagram of $E_{ef}{}^{\ast }$ versus frequency is shown in Figure 11, showing its sigmoidal shape.

In summary, it can be concluded that, at least in this example including a

Figure 11. Sigmoidal diagram of $E_{e{f}^{\mathrm{<mo>\; *\; </mo>}}}$ -frequency, resulting after removing the outliers.

large size data set, the model behaves well and as expected. Future analysis will be required to confirm the use of this model. However, the possibility of a four dimensional analysis, that is, including frequencies as the fourth variable is now possible.

10. Conclusions

The most important conclusions from this paper are:

1) The necessity of replacing or improving existing models, broadly acknowledged and recommend by standards and guidelines, needs to be emphasized to avoid the consequences of application of invalid models as reported in the literature.

2) A new method for considering the joint effect of maximum stress, stress ratio and frequency on the lifetime is presented in Sections 0 to 1. The frequency is included as a new and relevant variable to be considered in the design of important high-risk components such as wind turbines and automobiles, high speed railways and aeronautical material.

3) It has been shown that compatibility conditions are important tools to build models without arbitrary assumptions. These conditions were easily stated as functional equations, in Sections 2 to 4, which when solved provide the unique models satisfying the corresponding conditions. To our knowledge, no other method exists satisfying these conditions.

4) The model shows, in Section 5, that the frequency effect is equivalent to a concurrent modification of the maximum stress-stress ratio pairs from each test, what has immediate consequences in how to design the tests conducting to a precise parameter estimation.

5) The model parameters are classified, in Section 7, into two distinct classes: those defining the deterministic model, which are estimated by least-squares regression, and those reproducing the random ones, estimated by a restricted maximum-likelihood method. An estimation method is proposed and applied independently to both classes to avoid overfitting, in Section 7, and tested later in Section 9.

6) A model damage measure variable, V, is identified in Section 4, which together with its cdf, permits the damage accumulation to be evaluated and the associated probability of failure to be determined subject to any varying load history. A general formula is given in Section 8, which allows the cumulative damage to be analytically calculated from the joint effects of the three changing variables, maximum stress, stress ratio and frequency.

7) Due to its integrated character, because it includes all variables, the use of the model, will certainly allow us to reduce the sample size, and also to optimize the test sequences when used in designing test studies.

8) A more general formula that provides the damage of any load when the initial damage is known is obtained in Section 8 as the result of a functional equation that guaranties the compatibility, that is, damage to become independent on the way the loads are considered.

9) An unusually, rarely available large sample data, from an experimental campaign in which 243 fatigue specimens were tested, was used to validate the model, in Section 9, resulting very good results, as shown in the figure examples.

Acknowledgements

G. Khatibi, B. Czerny, and M. Zareghomsheh gratefully acknowledge the financial support by the Austrian Federal Ministry for Digital and Economic Affairs and the National Foundation for Research, Technology and Development for the support provided to the experimental program in this work.

Conflicts of Interest

The authors declare no conflicts of interest.

References

[1]	Berchtold, M. and Klopfer, I. (2019) Fatigue Testing at 1000 Hz Testing Frequency. Procedia Structural Integrity, 18, 532-537. https://doi.org/10.1016/j.prostr.2019.08.197
[2]	Stanzl-Tschegg, S. (2014) Very High Cycle Fatigue Measuring Techniques. International Journal of Fatigue, 60, 2-17. https://doi.org/10.1016/j.ijfatigue.2012.11.016
[3]	Stanzl-Tschegg, S. (1999) Fracture Mechanisms and Fracture Mechanics at Ultrasonic Frequencies. Fatigue and Fracture of Engineering Materials and Structures, 22, 267-279. https://doi.org/10.1046/j.1460-2695.1999.00180.x
[4]	Castillo, E. and Fernndez-Canteli, A. (2022) A Compatible Regression Weibull Model for the Description of the Three-Dimensional Fatigue S-N-R Field as a Basis for Cumulative Damage Approach. International Journal of Fatigue, 155, 1-14. https://doi.org/10.1016/j.ijfatigue.2021.106596
[5]	Castillo, E., Fernández-Canteli, A., Esslinger, V. and Thürlimann, B. (1985) Statistical Models for Fatigue Analysis of Wires, Strands and Cables. IABSE Proceedings, Vol. 82, 1-140.
[6]	Castillo, E. and Galambos, J. (1987) Lifetime Regression Models Based on a Functional Equation of Physical Nature. Journal of Applied Probability, 24, 160-169. https://doi.org/10.2307/3214067
[7]	Castillo, E. and Hadi, A.S. (1995) Modeling Lifetime Data with Application to Fatigue Models. Journal of the American Statistical Association, 90, 1041-1054. https://doi.org/10.1080/01621459.1995.10476606
[8]	Freudenthal, G.E. and Gumbel, A.M. (1956) Physical and Statistical Aspects of Fatigue. In: Advances in Applied Mechanics, Volume 4, Academic Press, New York, 117-158. https://doi.org/10.1016/S0065-2156(08)70372-7
[9]	Bolotin, V.V. (1969) Statistical Methods in Structural Mechanics. 1em plus 0.5em minus 0.4em. Holden-Day, San Francisco.
[10]	Bolotin, V.V. (1981) Wahrscheinlichkeitsmethoden zur Berechnung von Konstruktionen. Verlag fur Bauwesen, Berlin.
[11]	Bolotin, V.V. (1998) Mechanics of Fatigue. CRC Press, Boca Ratón.
[12]	Castillo, E. and Fernández-Canteli, A. (2009) A Unified Statistical Methodology for Modeling Fatigue Damage. Springer, Berlin.
[13]	Castillo, E. (1988) Extreme Value Theory in Engineering. Academic Press, San Diego.
[14]	Castillo, E., Fernández-Canteli, A. and Hadi, A. (1999) On Fitting a Fatigue Model to Data. International Journal of Fatigue, 21, 97-106. https://doi.org/10.1016/S0142-1123(98)00048-6
[15]	Castillo, E., Hadi, A.S., Balakrishnan, N. and Sarabia, J.M. (2005) Extreme Value and Related Models with Applications in Engineering and Science. Wiley, New York.
[16]	Castillo, E., Fernández-Canteli, A., Hadi, A.S. and López-Aenlle, M. (2007) A Fatigue Model with Local Sensitivity Analysis. Fatigue and Fracture of Engineering Materials and Structures, 30, 149-168. https://doi.org/10.1111/j.1460-2695.2006.01099.x
[17]	Gerber, W. (1874) Bestimmung der zulässigen Spannungen in Eisenkonstruktionen. Zeitschrift der Bayerischen Architecten und Ingenieure Vereins, 6, 101-110.
[18]	Goodman, J. (1930) Mechanics Applied to Engineering. 9th Edition, Longmans, Green and Co., London, 1.
[19]	Manson, S. (1953) Behaviour of Materials under Conditions of Thermal Stress. NACA Technical Note, Tech. Rep. 2933.
[20]	Coffin, L. (1954) A Study of the Effects of Cyclic Thermal Stresses on a Ductile Metal. Transactions of the American Society of Mechanical Engineers, 76, 931-950.
[21]	Morrow, J. (1968) Fatigue Properties of Metals. In: Fatigue Design Handbook, Soc. of Automotive Engineers, Warrendale, Pub. No. AE-4, ch. Section 3.2., 1-132.
[22]	Bagci, C. (1981) Fatigue Design of Machine Elements Using the “Bagci Line” Defining the Fatigue Failure Surface Line (Mean Stress Diagram). Mechanics and Machine Theory, 16, 339-359. https://doi.org/10.1016/0094-114X(81)90009-4
[23]	Wang, S., Dixon, M., Huey Jr., C. and Chen, S. (2000) The Clemson Limit Stress Diagram for Ductile Parts Subjected to Positive Mean Fatigue Loading. Journal of Mechanical Design, 122, 143-146. https://doi.org/10.1115/1.533557
[24]	Sekercioglu, T. (2009) A New Approach to the Positive Mean Stress Diagram in Mechanical Design. Materialwiss Werkstofftech, 40, 713-717. https://doi.org/10.1002/mawe.200900509
[25]	Nieslony, A. and Bhm, M. (2013) Mean Stress Effect Correction Using Constant Stress Ratio S-N Curves. International Journal of Fatigue, 52, 49-56. https://www.sciencedirect.com/science/article/pii/S0142112313000686 https://doi.org/10.1016/j.ijfatigue.2013.02.019
[26]	Kwofie, S. (2001) An Exponential Stress Function for Predicting Fatigue Strength and Life Due to Mean Stresses. International Journal of Fatigue, 23, 829-836. https://www.sciencedirect.com/science/article/pii/S0142112301000445 https://doi.org/10.1016/S0142-1123(01)00044-5
[27]	Dowling, N., Calhoun, C. and Arcari, N. (2009) Mean Stress Effects in Stress-Life Fatigue and the Walker Equation. Fatigue and Fracture of Engineering Materials and Structures, 32, 163-179. https://doi.org/10.1111/j.1460-2695.2008.01322.x
[28]	Papuga, J., Vizkova, I., Lutovinov, M. and Nesladek, M. (2018) Mean Stress Effect in Stress-Life Fatigue Prediction Re-Evaluated. MATEC Web of Conferences, 165, Article No. 10018. https://doi.org/10.1051/matecconf/201816510018
[29]	Pallarés-Santasmartas, L., Albizuri, J., Leguinagoicoa, N., Saintier, N. and Merzeau, J. (2019) The Effect of Mean Axial and Torsional Stresses in the Fatigue Strength of 34CrNiMo6 High Strength Steel. MATEC Web of Conferences, 300, Article No. 16004. https://doi.org/10.1051/matecconf/201930016004
[30]	Dowling, N. (1988) Estimation and Correlation of Fatigue Lives for Random Loading. International Journal of Fatigue, 10, 179-185. https://www.sciencedirect.com/science/article/pii/0142112388900606 https://doi.org/10.1016/0142-1123(88)90060-6
[31]	Dowling, N. (1972) Fatigue Failure Predictions for Complicated Stress-Strain Histories. Journal of Materials, 7, 271-287.
[32]	Arcari, A., Apetre, N., Sarkar, S., Iyyer, N., Dowling, N., Stanzl-Tschegg, S., Vasudevan, A., Phan, N. and Kang, P. (2012) Influence of Superimposed VCHF Loadings in Cyclic Fatigue of 7075-T6 Aluminum Alloy. 53rd AIAA/ASME/AHS/ASC Structures, Structural Dynamics and Materials Conference, Honolulu, 23-26 April 2012, 1-11. https://doi.org/10.2514/6.2012-1732
[33]	Apetre, N., Arcari, A., Dowling, N., Iyyer, N. and Phan, N. (2015) Probabilistic Model of Mean Stress Effects in Strain-Life Fatigue. Procedia Engineering, 114, 538-545. https://www.sciencedirect.com/science/article/pii/S1877705815017427 https://doi.org/10.1016/j.proeng.2015.08.103
[34]	Lu, S., Su, Y., Yang, M. and Li, Y. (2018) A Modified Walker Model Dealing with Mean Stress Effect in Fatigue Life Prediction for Aeroengine Disks. Mathematical Problems of Engineering, 2018, Article ID: 51482278. https://doi.org/10.1155/2018/5148278
[35]	Fine, J.P. and Gray, R.J. (1999) A Proportional Hazards Model for the Subdistribution of a Competing Risk. Journal of the American Statistical Association, 94, 496-509. https://doi.org/10.1080/01621459.1999.10474144
[36]	Castillo, E. and Fernández-Canteli, A. (2006) A Parametric Lifetime Model for the Prediction of High Cycle Fatigue Based on Stress Level and Amplitude. Fatigue and Fracture of Engineering Materials and Structures, 29, 1031-1038. https://doi.org/10.1111/j.1460-2695.2006.01068.x
[37]	Castillo, E., Fernández-Canteli, A. and Ruiz-Ripoll, M.L. (2008) A General Model for Fatigue Damage Due to Any Stress History. International Journal of Fatigue, 30, 150-164. https://doi.org/10.1016/j.ijfatigue.2007.02.011
[38]	Fernández-Canteli, A., Castillo, E., Blasón, S., Correia, J. and de Jesus, A. (2021) Generalization of the Weibull Probabilistic Model to Assess Fatigue Data into the Three Domains LCF, HCF and VHCF. International Journal of Fatigue, Under Review.
[39]	Saucedo, L., Yu, R.C., Medeiros, A., Zhang, X. and Ruiz, G. (2013) A Probabilistic Fatigue Model Based on the Initial Distribution to Consider Frequency Effect in Plain and Fiber Reinforced Concrete. International Journal of Fatigue, 48, 308-318. https://doi.org/10.1016/j.ijfatigue.2012.11.013
[40]	Melby, A. and Kobayashi, N. (1999) Damage Progression and Variability on Breakwater Trunks. In: Proceedings of Coastal Structures 1999, Balkema, Rotterdam, 309-316. https://doi.org/10.1061/9780784404119.141
[41]	Aczél, J. (1966) Lectures on Functional Equations and Their Applications. Mathematics in Science and Engineering. Academic Press, Cambridge, 19.
[42]	Aczél, J. (1987) A Short Course on Functional Equations Based upon Recent Applications to the Social and Behavioral Sciences. Springer, Berlin. https://doi.org/10.1007/978-94-009-3749-9
[43]	Eichhorn, W. (1978) Functional Equations in Economics, Ser. Applied Mathematics and Computation. Addison-Wesley Publishing Co., London, 11.
[44]	Castillo, E. and Ruiz Cobo, M. (1992) Functional Equations and Modelling in Science and Engineering, ser. Monographs and Textbooks in Pure and Applied Mathematics. Marcel Dekker Inc., New York, 161.
[45]	Castillo, E., Iglesias, A. and Ruz-Cobo, R. (2005) Functional Equations in Applied Sciences, ser. Mathematics in Science and Engineering. Elsevier B. V., Amsterdam, 199.
[46]	Hadi, A.S. and Simonoff, J.S. (1993) Procedures for the Identification of Multiple Outliers in Linear Models. Journal of the American Statistical Association, 88, 1264-1272. https://doi.org/10.1080/01621459.1993.10476407
[47]	Hadi, A.S. and Simonoff, J.S. (1997) A More Robust Outlier Identifier for Regression Data. Bulletin of the International Statistical Institute, 88, 281-282.
[48]	Castillo, E. and Sarabia, J.M. (1992) Engineering Analysis of Extreme Value Data. Journal of Waterway, Port, Coastal, and Ocean Engineering, ASCE, 118, 129-146. https://doi.org/10.1061/(ASCE)0733-950X(1992)118:2(129)
[49]	Bowman, A.W. and Azzalini, A. (1997) Applied Smoothing Techniques for Data Analysis. Oxford University Press, Oxford.