Efficient frequency-domain based fatigue life estimation of spot welds in vehicle components

The estimation of the damage and fatigue life of spot welds in vehicle components with hundreds or thousands thereof within reasonable computation time is a challenge even for recent computers in automotive industry. Hence, an advanced procedure in the frequency-domain is proposed. Whereas the general frequency-domain approach is well-established, it nevertheless becomes efficient for a large number of spot welds only with a package of measures. Among others, a two-scale FE-method using a pre-processed data base for various spot weld/sheet geometries and specified basic load cases, application of the cutting-plane method, and a priori stress estimations to reduce the number of spot welds to be checked completely are suggested. By the example of a tank of a cooling medium in a vehicle it is shown that the proposed algorithm may reduce the computation time as compared to time-domain calculations substantially.


Introduction
During driving, the components of a vehicle permanently are subject to various excitations which originate not only from road irregularities but also from interior sources like engine vibrations, and they lie within a wide frequencyrange. Hence, the estimation of the fatigue life of vehicle Eray Arslan eray.arslan@tuwien.ac.at 1 MAGNA STEYR Engineering, Liebenauer Hauptstraße 317, 8041 Graz, Austria 2 Institute of Mechanics and Mechatronics, TU Wien, Getreidemarkt 9, 1060 Vienna, Austria components is an important task in automotive industry; for a quite recent example see [1]. Since for advanced models with large numbers of degree of freedom fatigue analyses by FE-simulations in the time-domain are very timeconsuming, frequency-domain methods are applied more and more during the early design process. Surveys of these methods and underlying basic ideas can be found, e.g., in [2][3][4], and a comprehensive literature survey particularly of frequency-domain multi-axial dynamic fatigue analysis is given in the state-of-the-art paper [5]. Essentially, the method is based on the categorization of random loading and the resulting stresses -in the nodes of an FE-meshby power spectral density (PSD) functions and modelling the structure under consideration by linear transfer functions [3]; in order to make the paper self-contained, a short summary of the general approach is given in Sect. 2.
Whereas the computation time can be reduced significantly by using this procedure as compared to time-domain calculations, in complex vehicle structures with a huge number of FE-nodes the computational effort still may constitute some problem. In this regard, a particular challenge is the treatment of spot welds since a typical modern vehicle body-in-white contains about 4000-6000 thereof. Indeed, spot welds are the dominant joining elements, there, and essentially two failure modes may occur: Fracture of the weld nugget through the plane of the weld or nugget pullout, i.e. fracture of the sheet around the weld while the nugget remains intact [6]. Of course, a failure of spot welds affects the performance and stability of a structure significantly, and therefore the strength and/or fatigue life of spot welds as well as the relations between their failure and the dynamic behavior of the structure have been the subject of several investigations. For example, the latter issue is in the focus of the studies [6][7][8][9][10], and static strength and failure analyses of spot welds can be found, e.g., in [11,12]. The influence of the sheet thickness(es) and the spot diameter on the fatigue life is a topic of [13,14], while in [15,16] both various FE-models and the effects of different load cases are studied.
Indeed, with respect to the modelling of a spot weld, in essence there are two approaches in the literature: On the one hand comparatively simple FE-models with only a few elements [15] and on the other hand very detailed models with hundreds or even thousands thereof [11]. Both of these model-types have their benefits and drawbacks. Whereas the simple models are computationally very efficient, they provide a rather gross approximation of the stress state within the nugget and its immediate vicinity only. On the contrary, complex models are able to capture the stress distribution necessary for accurate fatigue life predictions quite exactly, but the effort for extremely fine meshes at all the spot welds in a vehicle body could hardly be justified from an engineering point of view [6]. Although the latter issue is of paramount relevance particularly in time-domain calculations, it nevertheless cannot be neglected in frequency-domain analyses, too. Moreover, in the hitherto published literature on welded joints fatigue problems mostly time-domain approaches have been taken. Some of the few exceptions thereof are the study [10] on the effects of spot welds damage, and a recent investigation of weld seams [17].
Hence, in the present study on spot welds several measures are proposed to speed up the computations significantly. For example, an efficient superposition of the stresses from a modal analysis, using a pre-processed database, with the PSD input may increase the performance and allow fast consideration of variants of the fatigue life estimations both for different PSD input curves and different spot weld and sheet geometries. Moreover, the so-called cutting plane method used hitherto in the timedomain proves advantageous in the frequency-domain, too. Additionally, by appropriate a priori stress estimations the number of spot welds to be considered further can be reduced essentially. To show the benefits of these measures, a fatigue life estimation of an actual vehicle component both with the suggested method in the frequency-domain and with standard time-domain approach is presented and discussed.
Thus, the focus of the present study is rather on the computational efficiency of the software for spot weld fatigue life estimation than on a discussion of specific fatigue criteria from a materials science point of view. In this study, the well-established (maximum) normal stress criterion is used. Indeed, application of the normal stress criterion proves quite beneficial with respect to the computation time. Hence, this widely used engineering fatigue criterion is applied here, irrespective of its simplicity; for a discussion of several classical and recent approaches see [5,[17][18][19][20][21] and the references quoted therein.
The paper is organized as follows: In Sect. 2, the basic approach in frequency-domain calculations is summarized. Subsequently, in Sect. 3, the computational procedure of spot weld fatigue analysis in the frequency-domain is discussed in detail, and both the analogies and differences to the corresponding time-domain calculations are addressed. Then, in Sect. 4 the computational efficiency of both approaches is compared. Finally, some concluding remarks are made in Sect. 5.

Basic approach in frequency-domain fatigue calculations
As mentioned in the Introduction, here a brief recapitulation (essentially based on [2,3]) of the underlying ideas for the frequency-domain fatigue analysis shall be given. It is presumed that both the (random) loading and the response of a structure can be categorized by PSD functions. The PSD function represents the Fourier transform of the autocorrelation function of a signal, and it does not contain any phase data. Moreover, the signals are understood to be ergodic stationary Gaussian random processes. In a stationary Gaussian random process the statistics are not influenced by a shift of the origin, and the ensemble probability density function is a Gaussian one. Ergodicity of a stationary process then means the equality of the statistics of one sample with those of the entire signal.
Next, the concept of frequency-domain analysis requires the determination of the transfer functions for the model of the structure. By a transfer function -which for complex structures must be determined by FE-methods -the ratio of output (stress) amplitude to input amplitude can be found for each frequency of a harmonic excitation. Input quantities here may be forces or accelerations (see the subsequent Sect. 3.1, Step b). Then, one obtains the PSD response by multiplying the transfer function with the input loading PSD; of course, this approach is based on linearity of the equations governing the structure's response. Now, from the stress PSD an estimate of the probability density function (PDF) of the occurring stress ranges must be found. To this aim, several procedures for its derivation from the stress PSD were proposed (see, e.g., [2,4,22]). In the present paper, the widely-used and well-established Dirlik's method for broad band processes is applied (see Sect. 3.2 below), in which the computation of the spectral moments of the stress PSD up to the fourth order is necessary. Then, from the PDF the numbers of stress cycles with certain stress amplitudes to be expected during the exposure time of the structure are calculated. Nevertheless, note that a determination of non-zero mean values of stress cycles from the PSD remains out of the scope of Dirlik's method.
Thereafter, based on the notion that each stress cycle causes a certain fatigue damage, the total damage caused by the entire action of the excitation is calculated by linear superposition of these partial damage values according to the Palmgren-Miner rule. Finally, the damage is determined with reference to the (given) S-N-curve of the material under consideration, which shows the relation between stress cycles to failure N f and stress range S.
In contrast to this, in the well-known standard procedure in the time-domain the numbers of stress cycles with certain amplitudes and mean stresses usually are obtained from the time history of the stresses by a rainflow count (e.g., [23]).

Computational procedure of spot weld fatigue analysis: detailed outline
Prior to the discussion of the steps in the spot-weld fatigue life estimation it shall be mentioned that the following algorithm is implemented in the program FEMSITE which has been developed at MAGNA STEYR company and works quite efficiently as will be shown below (see also [24]). Figure 1 shows the schematic of the overall procedure. Whereas the general algorithm is applicable to the entire structure, of course, the specific features in several steps of the proposed approach for spot weld analysis will be pointed out in the following. Moreover, the flowchart gives a unified representation for both frequency-domain and time-domain calculations, so that both the analogies and differences can be seen clearly.

Calculation of the stress-PSDs (frequencydomain) or the stresses (time-domain)
In Step (a), a modal analysis of the entire structure is performed by an FE software (NASTRAN or ABAQUS) to obtain the modal stresses Modal i (e.g., [25]). In particular, at the surfaces of the (usually thin) sheets connected by a spot weld a plane state of stress in its immediate vicinity is considered [26, p. 102]. Hence, the components of the modal stresses in a node have the indices i = x; y; xy, there, and the dimension of each modal stress component is identical with the number j of considered modes. It can be specified how many modes shall be taken into account or until which frequency level the modes are requested.
A part of the essential measures to speed up the calculations for spot welds in the subsequent steps is summarized in the orange boxes of Fig. 1, see Steps (a spotwelds,1 ) and (a spotwelds,2 ). As a pre-processing an FE-calculation of 48 basic load cases with unit forces and unit moments at a spot weld model with fine mesh (element size approximately 0.25 mm) is performed, see Fig. 2. In the FEMSITE database, results for more than 900 such models for various sheet thicknesses are available. To accelerate the computations further, just the stress tensors at the four circular contours of a spot weld on both sides of each sheet are saved; since at each circle there are 48 nodes, totally 192 evaluation nodes for each spot weld are considered in the following. The stress tensors at these nodes for some specific spot weld then are loaded to the main program. Thereafter, to obtain the actual stresses at the spot weld, a superposition of the stresses caused by appropriate basic load cases is performed: To this aim, the respective unit forces and moments are replaced by the actual ones. These actual forces and moments are obtained from the aforementioned modal analysis of the entire structure with coarse mesh (Fig. 3a), where the spot welds are represented by simple hexahedron elements (CHEXA) with beam elements (CBAR) at the edges perpendicular to the shell elements, see Fig. 3b (the connection is realized as mesh independent by the multi-point constraint equations, MPC). By proper transformations the forces in the nodes of this coarse spot weld model then are converted to equivalent actual grid point forces and moments acting as inputs in appropriate basic load cases with fine mesh, and the stresses therefore are obtained as multiples of the (unit load case) basic ones. Hence, essentially the whole procedure corresponds merely to a multiplication of pre-calculated stresses with appropriate factors to obtain the actual ones and thus is computationally very efficient. In Step (b), in the frequency-domain a frequency response analysis (analogous in time-domain: transient response analysis) with unit acceleration (e.g., 1 g) for one or more principal directions is performed; in the present analysis, as is common practice in most vehicle component fatigue studies, the focus is on an excitation in one principal direction, i.e., one channel excitation. Thus, transfer function T j .f / (where f denotes the frequency) and phase angle ' j .f / for each considered mode are obtained. The sizes of both the transfer function array and the phase angle array for each mode depend on the size Á f r of the frequency array which is determined by the chosen frequency range and the fineness of its discretization; usually, the latter is Step (b) can be performed in NASTRAN with "Solution sequence SOL111-Modal frequency response" or in ABAQUS with "Steady State Dynamics (SSD)". In Step (c), as already mentioned in the Introduction, the cutting-plane method established in the time-domain now also is applied in the frequency-domain. There are two rea-sons for this choice. The most important is the following one: If the normal stress fatigue criterion is applied, the extremum normal stress could be determined by calculating the eigenvalues of the stress tensor at each node, in principle [3]. However, for thousands of spot welds with a large number of FE-nodes each this would be computationally hardly feasible within reasonable CPU time, even with the acceleration measures mentioned above; hence, the cuttingplane method proves very efficient from the computational point of view. And a further advantage is the possibility of comparing the results of frequency-domain and timedomain calculations directly without any influence of a selected equivalent stress; nevertheless, if desired, in FEM-SITE the equivalent v. Mises stress (or the shear stress) can also be taken without difficulty.
The method is based on the calculation of the modal normal stress Modal In Step (d), in the frequency-domain now the superposition of the results of the frequency response analysis with the given acceleration PSD [Input (e)], G a .f /, is performed according to where G .f / denotes the stress PSD and summation has to be done for all modes considered (this equation is valid for one channel excitation). A particular advantage of this method is the possibility of a fast investigation of the fatigue behavior for various different loadings -corresponding, e.g., to specific road profiles -which are characterized by different PSD-input-functions G a .f /. This is due to the fact that for once determined normal stresses at unit acceleration input the stress PSDs then simply follow by a multiplication. In analogous time-domain computations, acceleration for a spot weld can be achieved by considering the normal stresses at critical time points only. These critical time points are defined as the time points where for each mode the maximum of the respective modal normal stress component occurs, and then the values of the modal normal stress components of the other modes at the same time point are superposed. This is done consecutively for all critical time points as long as the subsequently given inequality is fulfilled: n;max is compared with the fatigue limit f of the material under consideration, and if n;max < .2=3/ f for all of the 18 cutting planes, the further steps are skipped for this node, and the next one is considered (of course, also some other factor than 2=3 could be chosen).

Calculation of the fatigue damage and the estimated fatigue life
Based on the thus determined stress PSD, the damage and hence the fatigue life estimates can be calculated. Although this part of the procedure is essentially the same as in other frequency-domain investigations, for the global spot welds analysis in the full vehicle (or large parts thereof) special attention must be paid to the choice of the integration rou-K tine for the computation of the moments of the stress PSD as will be discussed below.
In Step (g), for the application of Dirlik's method to obtain the PDF, the spectral moments (k th moment of area) m k of the stress PSD, must be calculated. Here, f 1 and f n f r mean the smallest and largest frequency taken into account in the frequency array, respectively. Since the data are available in tabulated form only, a numerical integration technique for discrete and arbitrarily spaced abscissa values must be applied. Due to the large number of FE-nodes to be considered, the routine should be well-balanced between the requirements of sufficient accuracy and short necessary CPUtime. Hence, four different integration methods, i.e. trapezoidal rule, Simpson's rule (with two sub-variants) and the commercially available routine Davint [27] were tested for this purpose. As a result, the elementary trapezoidal rule according to proved the most efficient one. In particular, this holds true also with respect to numerical stability if relatively large differences between the abscissas of the data points occur (for other integrators, this requires additional CPU-timeconsuming error estimation routines). Now, Dirlik's method (e.g., [3]) is used to generate an estimation of the stress PDF. Since the original Dirlik's equations [28] are however formulated with respect to the stress range S, for a simpler implementation in FEMSITE they are re-written using the stress amplitude a = S=2. Thus, the probability density function reads with the abbreviations In the above relations, the irregularity factor I = E OE0 =E OEP is the ratio of the expected number of zeros, E OE0 = p m 2 =m 0 , Then, the rounded number of stress cycles for a certain stress amplitude, N i . a /, to be expected in an exposure time T (which is a user-defined input quantity) is calculated from the relation where the probability of an occurrence of a stress amplitude a;i is approximated by p a;i d a (see Fig. 5). Since the value of the maximum stress amplitude with non-negligible probability is not known a priori, the total number of discrete amplitude levels to be used also must be specified by Fig. 6 Example of a PSWT curve with illustration of the meaning of N i and N f;i the user, and by an appropriate algorithm then a;max and thus also the value of d a is calculated iteratively.
In contrast to the above approach, time-domain calculations follow the well-known procedure where for the prescribed acceleration input the actual normal stresses in a node are obtained by superposition of the modal ones (at all time points). Thereafter, as already mentioned in Sect. 2, one obtains the number of stress cycles for a certain stress amplitude by a rainflow count, which -unlike Dirlik's method -also yields occurring non-zero mean values. In Step (h), finally the results for the damage and the estimated fatigue life for a node are obtained. To this aim, the S-N-curve in the particular form of a PSWT-curve ( [29], [30,Chapt. 3 There, 0 f means the fatigue strength coefficient, " 0 f the fatigue ductility coefficient, b the fatigue strength exponent, and c the fatigue ductility exponent. Based on this relation, damage D is calculated by Palmgren-Miner's damage accumulation hypothesis [3], for a graphic representation of the approach see Fig. 6. Finally, the estimated fatigue life T f l is found as the ratio of exposure time to damage, i.e., T f l = T =D.
A short remark on the choice of the PSWT-curve for damage calculation should however be made. It was proposed in [31] to include the effect of non-zero values of the mean stress, which in the frequency-domain however cannot be obtained by Dirlik's method. Nevertheless, in the time-domain an appropriate mean stress correction for any cycle with different ratio of maximum to minimum stress can be taken into account [24], and in frequency-domain calculations in FEMSITE it is possible to add some constant stress values caused, e.g., by gravity or prestressing forces.

Comparison of computational efficiency of frequency-domain and time-domain approach
An assessment of the proposed approach now shall be given by the example of a bracket for an expansion tank of a cooling medium. This expansion tank is mounted at the front of a passenger car next to the radiator, see Fig. 7a; the mass of the tank with the cooling medium is 1.9 kg. In the bracket, during durability test runs cracks occurred (see Fig. 7b), but in the spot welds failure did not occur: They a priori were recognized as possibly sensitive parts and dimensioned appropriately.

Fig. 9 Calculated damage in the bracket
For the numerical simulation a measured acceleration signal was applied to the carrier. This was done on the one hand by the acceleration time history in the time-domain and on the other hand by the corresponding PSD in the frequency-domain (Fig. 8). In the latter domain, first the FEanalysis with unit acceleration for the frequency response was performed with ABAQUS "Steady State Dynamics (SSD)", and subsequently the damage estimation was done with FEMSITE. In the time-domain, for the FE-analysis ABAQUS "Modal Dynamics" was used, and thereafter also FEMSITE was applied, where time-domain routines are implemented, too (for details of the time-domain approach see [24]). Figure 9 shows the comparison of frequency-domain and time-domain results for the bracket. As one observes, the general distribution of critical damage is almost indistinguishable. The (left) most vulnerable position found by the calculations is perfectly coinciding with the position of crack occurrence in tests shown in Fig. 7b, thus validat- ing the simulations. With regard to the spot welds, they are confirmed to be in a safe state (within the load history), and the difference between the D− values obtained for them by both methods is only about 11% (see also Table 1). As compared to the time-domain approach, the frequency-domain method predicts a slightly smaller damage. Whereas time-domain and frequency-domain calculations yield absolutely comparable results, particularly for the spot welds, Table 1 reveals that the differences in the necessary CPU-times are enormous. In ABAQUS, a reduction by a factor of about 20 and in FEMSITE of about 4 was possible for the present example. Moreover, it should be noted that here a comparatively short load history was used and that for longer load histories still larger factors will occur since the data volume of the PSD is (practically) independent of the length of the load history.
Furthermore, it shall be pointed out that for vehicle components with a much larger number of spot welds the reduction in CPU-time as compared to time-domain calculations will be still more pronounced; summing up, a substantial reduction will be possible by application of the above method.

Concluding remarks
In this paper, an advanced algorithm for the damage and fatigue life estimation for spot welds was proposed. Whereas the general frequency-domain approach is well-established, of course, the following special measures and features give rise to a superior computational performance for spot weld calculations: (i) Initially, a comparatively coarse FE-mesh for the entire vehicle component is used, and the grid point forces at the spot welds then are transformed to input forces for spot weld models with fine mesh, which are available in a pre-processed data base. These models (in FEMSITE more than 900) take various geometrical properties and basic unit load cases into account, and the actual modal stresses then simply are obtained as multiples of the ones caused by the unit loads. (ii) The cutting-plane method, hitherto used primarily in time-domain calculations, is applied and proves very efficient in the frequency-domain, too, since, e.g., time-consuming eigenvalue calculations are avoided; moreover, alternatively normal stress or shear stress (or also v. Mises equivalent stress) can be used in the further steps. (iii) By an appropriate a priori estimation of the ratio of maximum occurring (normal) stress to fatigue limit a substantial reduction of the spot welds to be considered in the further steps is achieved. (iv) Applying the classic trapezoidal rule in the integrations necessary for Dirlik's method yields the optimum balance between speed of computation, accuracy, and numerical stability. (v) Based on the above approach, the influences of qualitatively different load histories caused by, e.g., various road types can be studied within short time. This essentially is due to the fact that the calculation of the stress PSDs is performed by a computationally fast multiplication process using the input acceleration PSDs. Furthermore, since the latter may be considered independent of the length of the load history, the data volume to be handled remains more or less constant.
Thus, the proposed procedure gives the automotive engineer a computationally very efficient tool to estimate the damage and fatigue life of vehicle components with a large number of spot welds. Of course, it would be worthwhile to extend the computational scheme also to the significantly more intricate case of multi-channel excitation; this however remains a task of future studies.