Fatigue Analysis of a Nonproportionally Loaded Shaft with a Fillet

This benchmark model is based on the example found in section 5.4.3 of Ref. 1. It shows how to perform a high-cycle fatigue analysis for nonproportional loading using critical plane methods.

Model Definition

The geometry is a circular shaft with two different diameters, 10 mm and 16 mm. At the transition between the two diameters there is a fillet with a radius of 2 mm.

Figure 1: The notched shaft.

Two time-dependent loads are applied at the small end of the shaft: a transverse force and a twisting moment. The force varies between 0 and 1.94 kN and the torque varies between −28.7 and +28.7 Nm. Figure 2 shows the history of one loading cycle.

Figure 2: Load history.

The big end of the shaft is fixed. The material is Elastic with E = 100 GPa and ν = 0.

In Ref. 1 it is stated that the fatigue limit for completely reversed axial tension is 700 MPa, while the fatigue limit for pure tension is 560 MPa. These values are the stress amplitudes. In uniaxial loading, the Findley criterion can be written as

where σa is the stress amplitude and σmax is the maximum stress experienced in a fatigue cycle. This means that you have to solve the simultaneous equations

to get the Findley parameters f and k. The result is f = 440 MPa and k = 0.23.

The Matake criterion is similar to the Findley criterion, with the difference that the critical plane is defined solely by the maximum shear stress. For a uniaxial case, the Matake expression is

which gives the corresponding system of equations as

The solution is f = 467 MPa and k = 0.33 as parameters for the Matake case.

The Dang Van criterion utilizes the stress history on a plane that has undergone an elastic shakedown. It evaluates fatigue using the expression

where τmax is the maximum shear stress and σH is the hydrostatic stress, while a and b are material parameters. For the fatigue tests above, the Dang Van relation transforms to

where σa is the stress amplitude and σmax is the maximum stress. The material constants are obtained from the following set of equations

The material parameters of the Dang Van model are b = 467 MPa and a = 0.5.

Results and Discussion

Figure 3 and Figure 4 show stress distributions from the two basic load cases. The location for the maximum equivalent stress is at the surface of the fillet, at a radius slightly larger than the minimum radius of the shaft.

In Figure 5 the equivalent stress from the combined load case with transverse force and positive torque is shown. It is symmetric with respect to the XY-plane and is identical also for the case when the torque is reversed.

Figure 3: Axial stress from transverse force.

Figure 4: Equivalent stress from torque.

Figure 5: Equivalent stress distribution for one of the combined load cases.

The results from the fatigue evaluation is shown in Figure 6, Figure 7, and Figure 8. With the Findley criterion, the fatigue usage factor is computed to 0.98, in perfect agreement with Ref. 1.

Using the Matake criterion, the fatigue usage factor decreases to 0.90, which shows that there can be significant differences between results from seemingly similar models. The critical plane computed in the Matake model differs from the one used in the Findley model. As a consequence, the maximum normal stress on the critical plane can be significantly lower in the Matake case than in the Findley case.

The Dang Van model predicts the fatigue usage factor to 0.94. This is slightly larger than 0.92 as predicted in Ref. 1. The difference arises from the fact that in the reference the shear stress arising from the constant bending force is approximated to zero. In reality, this shear force is 16% as compared with the shear force arising from the twisting moment and has an impact on the results. Moreover, when a comparison between different models is made, attention must be made to the discretization of the geometry. The critical point of each model is found in different locations. Thus, a prefect agreement in one model does not mean that there will be a perfect agreement in another model. Based on the discretization, the evaluation for certain models is made closer to its critical point and thus results are in better agreement with the theoretical value. A stress state that is only a few percent from its theoretical value can have a larger impact on the fatigue prediction, since the fatigue model requires further data manipulation.

Figure 6: Fatigue usage factor using the Findley criterion.

Figure 7: Fatigue usage factor using the Matake criterion.

Figure 8: Fatigue usage factor using the Dang Van criterion.

There is a large difference in the fatigue usage factor between the top (tension) and bottom (compression) side of the bar (see Figure 9), even though the equivalent stress is the same at both locations. This shows how the criterion captures the difference between the predominantly tensile stress states at the critical spot and the compressive stress states on the other side. Moreover, there is a clear difference between the models. In tension, the Dang Van model predicts higher fatigue usage factor than the Matake model, while in compression, this is reversed. Therefore, the computed values should not be taken as an exact prediction of the fatigue but rather as a probable outcome. Each model has a preferred area of use and application in a slightly different domain can result in less accurate fatigue prediction.

Figure 9: Fatigue usage factor prediction in tension and compression.

Notes About the COMSOL Implementation

In this model, you use the load case functionality in COMSOL Multiphysics to produce the load cycle. In the first study, the two basic load cases are analyzed. This study is not essential for the analysis, but it allows you to inspect the results of the individual basic load cases.

Reference

1. D.F. Socie and G.B. Marquis, Multiaxial Fatigue, SAE, 1999.

Fatigue_Module/Stress_Based/shaft_with_fillet

Modeling Instructions

From the menu, choose .

New

In the window, click

Model Wizard

In the window, click

In the tree, select .

Click .

Click

In the tree, select .

Click

Geometry 1

In the window, under click .

In the window for , locate the section.

From the list, choose .

Work Plane 1 (wp1)

In the toolbar, click

In the window for , click

Work Plane 1 (wp1)>Plane Geometry

In the window, click .

Work Plane 1 (wp1)>Polygon 1 (pol1)

In the toolbar, click

In the window for , locate the section.

From the list, choose .

In the text field, type 0 0 32 32 50 50 0.

In the text field, type 0 5 5 8 8 0 0.

Work Plane 1 (wp1)>Fillet 1 (fil1)

In the toolbar, click

On the object , select Point 3 only.

In the window for , locate the section.

In the text field, type 2.

Click

Revolve 1 (rev1)

In the window, right-click and choose .

In the window for , locate the section.

Clear the check box.

Locate the section. Find the subsection. In the text field, type 1.

In the text field, type 0.

Click

Solid Mechanics (solid)

Fixed Constraint 1

In the window, under right-click and choose .

Select Boundaries 21–24 only.

Rigid Connector 1

In the toolbar, click

and choose .

Select Boundaries 1, 3, 5, and 7 only.

Applied Force 1

In the toolbar, click

and choose .

In the window for , locate the section.

Specify the F vector as


0	x
0	y
-1.94[kN]	z

In the toolbar, click

and choose .

Rigid Connector 1

In the window, click .

Applied Moment 1

In the toolbar, click

and choose .

In the window for , locate the section.

Specify the M vector as


28.7[N*m]	x
0	y
0	z

In the toolbar, click

and choose .

Global Definitions

Transverse force

In the window, under click .

In the window for , type Transverse force in the text field.

In the text field, type lgF.

Twisting moment

In the window, under click .

In the window for , type Twisting moment in the text field.

In the text field, type lgM.

Materials

Material 1 (mat1)

In the window, under right-click and choose .

In the window for , locate the section.

In the table, enter the following settings:


Property	Variable	Value	Unit	Property group
Young’s modulus	E	100[GPa]	Pa	Basic
Poisson’s ratio	nu	0	1	Basic
Density	rho	0	kg/m³	Basic

Mesh 1

In the window, under click .

In the window for , locate the section.

From the list, choose .

Click

A finer mesh is needed in the fillet in order to resolve the stress concentration.

Locate the section. From the list, choose .

Size 1

In the window, right-click and choose .

In the window for , locate the section.

From the list, choose .

Size 2

Right-click and choose .

In the window for , locate the section.

From the list, choose .

Select Edges 13, 14, 16, and 18 only.

Locate the section. Click the button.

Locate the section. Select the check box.

In the associated text field, type 0.5.

Select the check box.

In the associated text field, type 1.2.

Click

Study 1

Step 1: Stationary

In the window, under click .

In the window for , click to expand the section.

Select the check box.

Click two times.

In the table, enter the following settings:


Load case	lgF	Weight	lgM	Weight
Transverse force	√	1.0		1.0
Twisting moment		1.0	√	1.0

In the window, click .

In the window for , type Study 1 (Basic load cases) in the text field.

In the toolbar, click

Results

Stress (solid)

Visualize the difference between the tension and compression on the opposite sides of the shaft.

In the toolbar, click

In the window for , locate the section.

From the list, choose .

Surface 1

In the window, expand the node, then click .

In the window for , click in the upper-right corner of the section. From the menu, choose .

In the toolbar, click

Add Study

In the toolbar, click

to open the window.

Go to the window.

Find the subsection. In the tree, select .

Click in the window toolbar.

In the toolbar, click

to close the window.

Study 2

Step 1: Stationary

In the window for , locate the section.

Select the check box.

Click three times.

In the table, enter the following settings:


Load case	lgF	Weight	lgM	Weight
No load		1.0		1.0
+F -M	√	1.0	√	-1.0
+F +M	√	1.0	√	1.0

In the window, click .

In the window for , type Study 2 (Combined load cases) in the text field.

In the toolbar, click

Results

Stress (solid) 1

In the window, expand the node, then click .

In the toolbar, click

As a last step, perform a fatigue analysis on the load cycle.

Add Physics

In the toolbar, click

to open the window.

Go to the window.

In the tree, select .

Find the subsection. In the table, clear the check boxes for and .

Click in the window toolbar.

In the toolbar, click

to close the window.

Fatigue (ftg)

Findley

Right-click and choose the boundary evaluation .

In the window for , type Findley in the text field.

Locate the section. From the list, choose .

Locate the section. In the Q text field, type 16.

Matake

In the toolbar, click

and choose .

In the window for , type Matake in the text field.

Locate the section. From the list, choose .

Locate the section. In the Q text field, type 16.

Dang Van

In the toolbar, click

and choose .

In the window for , type Dang Van in the text field.

Locate the section. From the list, choose .

Materials

Because the fatigue model is active only on the boundaries, you need to define a material on the boundaries.

Material 2 (mat2)

In the window, under right-click and choose .

In the window for , locate the section.

From the list, choose .

Locate the section. In the table, enter the following settings:


Property	Variable	Value	Unit	Property group
Normal stress sensitivity coefficient	k_Findley	0.23	1	Findley
Limit factor	f_Findley	440[MPa]	Pa	Findley
Normal stress sensitivity coefficient	k_Matake	0.33	1	Matake
Limit factor	f_Matake	467[MPa]	Pa	Matake
Hydrostatic stress sensitivity coefficient	a_DangVan	0.5	1	Dang Van
Limit factor	b_DangVan	467[MPa]	Pa	Dang Van