Uncertainty analysis of torque-induced bending deformations in ball screw systems

Abstract

In order to study the torque-induced bending deformation in ball screw systems with uncertain-but-bounded parameters, uncertainty analysis of torque-induced bending deformations in ball screw systems based on an numerical method combined second-order Taylor series expansion with difference of convex functions algorithm is conducted in this research. The maximum torque-induced bending deformation in ball screw systems with certain parameters is solved by modeling the ball screw as a simple Euler–Bernoulli beam. On the basis of the model, the general expression of maximum deformation of ball screw systems in terms of uncertain parameters is derived based on the second-order Taylor series expansion. Then the problem of determining bounds of the maximum deformation of uncertain systems is transformed into a series of box constrained quadratic programming problems. Box constrained quadratic programming problems are solved by difference of convex functions algorithm. The effectiveness of the analysis is validated by Monte Carlo simulation method. Through demonstration, the conclusion can be drawn that the torque-induced bending deformation in ball screw systems with uncertain-but-bounded parameters can be analyzed explicitly by the numerical method combined second-order Taylor series expansion with difference of convex functions algorithm, which offers a new approach for estimating the torque-induced bending deformation in ball screw systems.

Introduction

As an important force and motion transfer device, ball screw mechanism (BSM) has the capability of converting rotation into linear motion, or vice versa. It is widely applied to machine tools, robots, measuring instruments, and automatic equipment for its advantages of precise positioning and high efficiency. ¹ However, the deformations and vibrations of the ball screw drive components adversely affect the positioning accuracy and fatigue life of the drive. ²

Many works have been conducted on the deformations of the ball screw drive components. Zhang and Sun ³ simplified the screw as a Timoshenko beam, and the screw–nut interface is simulated by the stiffness and damp in six directions considering the effects of pretension, shear, and moment of inertia on the axial, lateral, and torsional vibration. Zaeh et al. ⁴ established a finite element model of ball screw drive systems considering the lateral deformation of the ball screw, but both the models used for the screw–nut interface stiffness matrix fail to capture the coupling between the axial, torsional, and lateral dynamics of ball screw drives. Okwudire and colleagues^5,6 demonstrated theoretically and experimentally that the nut introduces a coupling between the dynamics in the axial/torsional directions and the lateral (bending) direction in BSM. Torque-induced bending deformations were optimized to minimize the static coupling between the applied torque and lateral deformations of the screw in Okwudire. ⁷ However, from the survey in ball screw manufacturers, the authors found that the preload of the ball screws always depends on the experiences of workers. What is more, wear will occur in the ball screw in the course of service. The preload will degrade gradually. The uncertainty of preload breeds the uncertainty of axial stiffness. The wear and preload degradation result in change in the contact position, which gives rise to the uncertainty of contact angle. But there are few researches on uncertainty analysis of torque-induced bending deformations in ball screw systems.

Uncertainty analysis of torque-induced bending deformations in ball screw systems based on a numerical method combined second-order Taylor series expansion with difference of convex functions algorithm (DCA) is conducted in this research. The maximum torque-induced bending deformation in ball screw systems with certain parameters is solved by modeling the ball screw as a simple Euler–Bernoulli beam. On the basis of the model, the general expression of maximum deformation of ball screw systems in terms of uncertain parameters is derived based on the second-order Taylor series expansion. The problem of determining bounds of the maximum deformation of uncertain systems is transformed into a series of box constrained quadratic programming (QB) problems. The problems are solved by DCA. The effectiveness of the analysis is validated by Monte Carlo simulation (MCS) method. The bounds of torque-induced bending deformations in a certain ball screw system are discussed in this article as well.

Maximum torque-induced bending deformation in ball screw systems

A detailed model for maximum torque-induced bending deformation in ball screw systems is presented in Okwudire. ⁷ The result is a dimensionless displacement in lateral direction. This section summarizes the derivation of the maximum torque-induced bending deformation in ball screw systems to provide the reader with some context needed for the rest of the article.

Screw–nut interface stiffness model

Screw–nut interface stiffness modeling is the key problem in deformation analysis of ball screw systems. Therefore, screw–nut interface stiffness model should be established to obtain the accurate model of ball screw systems.

Axial preload is usually applied to two nuts in the ball screw by inserting a shim between them. Thus, balls contacting with the two nuts have different contact status, one is L-configuration and the other is U-configuration, which is depicted in Figure 1.

Figure 1.

Schematic diagram of preload in ball screw.

In the process of screw–nut interface stiffness modeling, the mass of balls is neglected. The contact stiffness is considered as $k_{\mathrm{ball}}$ , and its direction is along the contact normal direction, as shown in Figure 2.

Figure 2.

Simplification of balls.

In order to convert the contact stiffness of balls to the screw–nut interface stiffness, two steps are conducted. First, the contact stiffness of balls is converted to the global coordinate. Second, the stiffness of a single ball obtained by the first step is accumulated into the screw–nut interface stiffness. Detail solution process is listed as follows.

Contact stiffness conversion of a single ball

According to the geometrical relationship of the ball screws, equation (1) can be obtained

$$p=\tan \alpha ·\pi d_{\mathrm{p}}$$

(1)

where $p$ is the pitch, $\alpha$ is the pitch angle of the ball screws, and $d_{\mathrm{p}}$ is the pitch diameter of the ball screws.

Two local coordinates CS_1L and CS_2L are established as shown in Figure 3. CS_1U and CS_2U are also established in the similar way. Where the z-axis of CS_1L (CS_1U ) coordinate points along the contact direction, and the z-axis of CS_2L (CS_2U ) coordinate (global coordinate) points along the axial direction of the ball screw shaft. The subscripts L and U represent different contact conditions.

Figure 3.

Schematic diagram of local coordinates in L-configuration.

Taking L-configuration, for example, $T_{2L-1L}$ represents the transformation matrix from CS_1L to CS_2L , and it can be expressed as

$$T_{2L-1L}={\mathrm{rot}}_y(-\alpha )·{\mathrm{rot}}_x(\beta )$$

(2)

where

$$\begin{array}{cc}\multicolumn{1}{c}{}& {\mathrm{rot}}_x(\theta )=\left[\begin{array}{ccc}\hfill 1\hfill & \hfill 0\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{0}& \hfill \cos \theta \hfill & \hfill -\sin \theta \hfill \\ \multicolumn{1}{c}{0}& \hfill \sin \theta \hfill & \hfill \cos \theta \hfill \end{array}\right],\hfill \\ \multicolumn{1}{c}{}& {\mathrm{rot}}_y(\theta )=\left[\begin{array}{ccc}\hfill \cos \theta \hfill & \hfill 0\hfill & \hfill \sin \theta \hfill \\ \multicolumn{1}{c}{0}& \hfill 1\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{-\sin \theta }& \hfill 0\hfill & \hfill \cos \theta \hfill \end{array}\right],\hfill \\ \multicolumn{1}{c}{}& {\mathrm{rot}}_z(\theta )=\left[\begin{array}{ccc}\hfill \cos \theta \hfill & \hfill -\sin \theta \hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{\sin \theta }& \hfill \cos \theta \hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{0}& \hfill 0\hfill & \hfill 1\hfill \end{array}\right]\hfill \end{array}$$

Similarly, for balls on the U-configuration, $T_{2U-1U}$ can be expressed as

$$T_{2U-1U}={\mathrm{rot}}_y(\alpha )·{\mathrm{rot}}_x(-\beta )$$

(3)

$T_{2L-z1L}(T_{2U-z1U})$ is defined as $T_{2L-1L}(T_{2U-1U})$ postmultiplied by a unit vector in the direction of z1-axis

$$\begin{array}{cc}\multicolumn{1}{c}{}& T_{2L-z1L}=T_{2L-1L}·{[\begin{array}{ccc}\hfill 0\hfill & \hfill 0\hfill & \hfill 1\hfill \end{array}]}^T;\hfill \\ \multicolumn{1}{c}{}& T_{2U-z1U}=T_{2U-1U}·{[\begin{array}{ccc}\hfill 0\hfill & \hfill 0\hfill & \hfill 1\hfill \end{array}]}^T\hfill \end{array}$$

(4)

CS₃ coordinate is established with its z-axis parallel to the axial direction of the ball screw shaft. Its origin is located at the center of each ball. Its orientation depends on the azimuth angle $\phi$ .The azimuth angle $\phi$ represents the angle measured in the counterclockwise direction from the x-axis in CS₃ coordinate to the x-axis in CS_2L (CS_2U ) coordinate.

The transformation matrix $T_{3-2L}(T_{3-2U})$ from CS_2L (CS_2U ) to CS₃ can be expressed as

$$T_{3-2L}={\mathrm{rot}}_z(\phi ),T_{3-2U}={\mathrm{rot}}_z(\phi )·\left[\begin{array}{ccc}\hfill 1\hfill & \hfill 0\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{0}& \hfill -1\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{0}& \hfill 0\hfill & \hfill -1\hfill \end{array}\right]$$

(5)

Therefore, the transformation matrix $T_{3-z1L}(T_{3-z1U})$ from each local coordinate to global coordinate can be listed as follows

$$T_{3-z1L}=T_{3-2L}·T_{2L-z1L}{T}_{3-z1U}=T_{3-2U}·T_{2U-z1U}$$

(6)

The matrix above is orthogonal, so the inverse matrix $T_{z1L-3}(T_{z1U-3})$ could be obtained by transposing the matrix $T_{3-z1L}(T_{3-z1U})$ .

Accumulation of the contact stiffness of balls

CS coordinate is established by considering the center of mass of nuts, which is marked as P, as its origin. Z-axis of the coordinate is along the axial direction of the ball screw shafts. The position of each ball could be represented as $\mathbf{r}$ , which in L-configuration is shown in Figure 4

$$\mathbf{r}=\left\{\begin{array}{c}\hfill r_x\hfill \\ \multicolumn{1}{c}{r_y}\\ \multicolumn{1}{c}{r_z}\end{array}\right\}=\left\{\begin{array}{c}\hfill R\sin \phi \hfill \\ \multicolumn{1}{c}{-R\cos \phi }\\ \multicolumn{1}{c}{r_g\phi }\end{array}\right\}$$

(7)

Figure 4.

Vector diagram of balls in L-configuration.

where $R$ is the constant radius measured from the axis of the screw to the center of each ball, $r_g=p/(2\pi )$ .

The contact point between the ball screw shaft and the ball is marked as ${\mathrm{P}}_{\mathrm{BS}}$ . The contact point between the nut and the ball is marked as ${\mathrm{P}}_{\mathrm{N}}$ ; ${\mathbf{r}}_{\mathrm{BS}}$ and ${\mathbf{r}}_{\mathrm{N}}$ are the position vectors from P to ${\mathrm{P}}_{\mathrm{BS}}$ and ${\mathrm{P}}_{\mathrm{N}}$ , respectively. ${\mathbf{r}}_{\mathrm{BS}}$ and ${\mathbf{r}}_{\mathrm{N}}$ can be expressed as

$${\mathbf{r}}_{\mathrm{BS}}=\mathbf{r}-{\mathbf{r}}_{\mathbf{ball}},{\mathbf{r}}_{\mathrm{N}}=\mathbf{r}+{\mathbf{r}}_{\mathbf{ball}}$$

(8)

where

$${\mathbf{r}}_{\mathbf{ball}}=T_{3-z1L}·R_{\mathrm{ball}}(L\mathrm{-}\mathrm{configuration})$$

(9a)

$${\mathbf{r}}_{\mathbf{ball}}=T_{3-z1U}·R_{\mathrm{ball}}(U\mathrm{-}\mathrm{configuration})$$

(9b)

According to the theory in Okwudire and Altintas, ⁵ the relationship below can be obtained

$$\begin{array}{cc}\multicolumn{1}{c}{}& {\mathbf{u}}_{3\mathrm{BS}}=\left[\begin{array}{cc}\hfill {\mathbf{I}}_{3\times 3}\hfill & \hfill -\mathbf{S}({\mathbf{r}}_{\mathrm{BS}})\hfill \end{array}\right]\left\{\begin{array}{c}\hfill {\mathit{u}}_{\mathrm{BS}}\hfill \\ \multicolumn{1}{c}{{\mathbf{\theta }}_{\mathrm{BS}}}\end{array}\right\},\hfill \\ \multicolumn{1}{c}{}& {\mathbf{u}}_{3\mathrm{N}}=\left[\begin{array}{cc}\hfill {\mathbf{I}}_{3\times 3}\hfill & \hfill -\mathbf{S}({\mathbf{r}}_{\mathrm{N}})\hfill \end{array}\right]\left\{\begin{array}{c}\hfill {\mathit{u}}_N\hfill \\ \multicolumn{1}{c}{{\mathbf{\theta }}_N}\end{array}\right\}\hfill \end{array}$$

(10)

where $[\begin{array}{cc}\hfill {\mathbf{I}}_{3\times 3}\hfill & \hfill -\mathbf{S}\hfill \end{array}({\mathbf{r}}_{\mathrm{BS}})]$ can be replaced by ${\mathbf{T}}_{3-\mathrm{BS}}$ , and $[\begin{array}{cc}\hfill {\mathbf{I}}_{3\times 3}\hfill & \hfill -\mathbf{S}\hfill \end{array}({\mathbf{r}}_{\mathrm{N}})]$ can be replaced by ${\mathbf{T}}_{3-\mathrm{N}}$ . ${\mathbf{u}}_{\mathrm{BS}}$ and ${\mathbf{\theta }}_{\mathrm{BS}}$ represent the displacement and angle displacement of the ball screw shaft, respectively. ${\mathbf{u}}_{\mathrm{N}}$ and ${\mathbf{\theta }}_{\mathrm{N}}$ represent the displacement and angle displacement of the nut, respectively.

The function $\mathbf{S}(\mathbf{r})$ could be expressed as follows

$$\mathbf{S}(\mathbf{r})=\left[\begin{array}{ccc}\hfill 0\hfill & \hfill -r_z\hfill & \hfill r_y\hfill \\ \multicolumn{1}{c}{r_z}& \hfill 0\hfill & \hfill -r_x\hfill \\ \multicolumn{1}{c}{-r_y}& \hfill r_x\hfill & \hfill 0\hfill \end{array}\right]$$

(11)

The axial stiffness of the screw–nut interface can be expressed as follows ⁸

$$k_{\mathrm{ax}}=0.8\times K\left[\frac{F}{\epsilon c}\right]$$

(12)

where $K$ is the stiffness that can be searched in manufacturers’ catalogs, and the unit is N/µN; F is the preload, the unit is N; $\epsilon$ is the factor for computing stiffness, and generally the value is 0.1; $c$ is the basic rated dynamical load, and the unit is N.

The contact stiffness of the balls is

$$k_{\mathrm{ball}}=\frac{k_{\mathrm{ax}}}{N_{\mathrm{ball}}{\cos }^2\alpha {\cos }^2\beta }$$

(13)

where $N_{\mathrm{ball}}$ is the number of balls in contact status.

The stiffness between ${\mathrm{P}}_{\mathrm{BS}}$ and ${\mathrm{P}}_{\mathrm{N}}$ is

$${\mathbf{K}}_{\mathrm{ball}}=k_{\mathrm{ball}}·\left[\begin{array}{cc}\hfill 1\hfill & \hfill -1\hfill \\ \multicolumn{1}{c}{-1}& \hfill 1\hfill \end{array}\right]$$

(14)

The displacement of ${\mathrm{P}}_{\mathrm{BS}}$ and ${\mathrm{P}}_{\mathrm{N}}$ is transformed to CS_1L (CS_1U ) coordinate. The transformation matrix is

$$T_{z1L-\mathrm{BS}}=T_{z1L-3}·T_{3-\mathrm{BS}},T_{z1L-N}=T_{z1L-3}·T_{3-N}$$

(15a)

$$T_{z1U-\mathrm{BS}}=T_{z1U-3}·T_{3-\mathrm{BS}},T_{z1U-N}=T_{z1U-3}·T_{3-N}$$

(15b)

A matrix composed by the transformation matrix of the nut and the ball screw is

$$\begin{array}{c}\multicolumn{1}{c}{T_{z1L-\mathrm{BSN}}=\left[\begin{array}{cc}\hfill T_{z1L-\mathrm{BS}}\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{0}& \hfill T_{z1L-N}\hfill \end{array}\right],}\\ \multicolumn{1}{c}{T_{z1U-\mathrm{BSN}}=\left[\begin{array}{cc}\hfill T_{z1U-\mathrm{BS}}\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{0}& \hfill T_{z1U-N}\hfill \end{array}\right]}\end{array}$$

(16)

Contribution of each ball to the stiffness of the screw–nut interface is

$${\mathit{K}}_L={\mathit{T}}_{z1L-\mathrm{BSN}}^T·{\mathbf{K}}_{\mathrm{ball}}·T_{z1L-\mathrm{BSN}}$$

(17a)

$${\mathit{K}}_U={\mathit{T}}_{z1U-\mathrm{BSN}}^T·{\mathbf{K}}_{\mathrm{ball}}·T_{z1U-\mathrm{BSN}}$$

(17b)

The stiffness of the screw–nut interface can be accumulated as follows

$${\mathit{K}}_{\mathrm{SN}}={\displaystyle \sum _{k=1}^{N_{\mathrm{ball}}}}(\frac{1}{{\phi }_{k,\mathrm{end}}-{\phi }_{k,\mathrm{st}}}{\displaystyle \underset{{\phi }_{k,\mathrm{st}}}{\overset{{\phi }_{k,\mathrm{end}}}{\int }}}{\mathit{K}}_{L/U}(\phi )d\phi )$$

(18)

where ${\mathit{K}}_{L/U}$ is determined by the contact status of the kth ball whether ${\mathit{K}}_{L/U}$ is chosen ${\mathit{K}}_L$ or ${\mathit{K}}_U$ . ${\phi }_{k,\mathrm{st}}$ and ${\phi }_{k,\mathrm{end}}$ are the entry and exit angles, respectively.

${\mathrm{K}}_{\mathrm{SN}}$ is a 12 × 12 stiffness matrix which connects the 6 degrees of freedom (DOFs) of the nut to the 6 DOFs of the screw. The stiffness is determined by its first six rows and columns

$${\mathit{K}}_{\mathrm{SN}1-6,1-6}=\left[\begin{array}{ccc}\hfill \begin{array}{cc}\hfill k_{x,x}\hfill & \hfill \hfill \\ \multicolumn{1}{c}{0}& \hfill k_{y,y}\hfill \end{array}\hfill & \hfill \begin{array}{cc}\hfill \hfill & \hfill \hfill \\ \multicolumn{1}{c}{}& \hfill \mathrm{sym}\hfill \end{array}\hfill & \hfill \begin{array}{cc}\hfill \hfill & \hfill \hfill \\ \multicolumn{1}{c}{}& \hfill \hfill \end{array}\hfill \\ \multicolumn{1}{c}{\begin{array}{cc}\hfill k_{x,z}\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{k_{x,\theta x}}& \hfill 0\hfill \end{array}}& \hfill \begin{array}{cc}\hfill k_{z,z}\hfill & \hfill \hfill \\ \multicolumn{1}{c}{k_{z,\theta x}}& \hfill k_{\theta x,\theta x}\hfill \end{array}\hfill & \hfill \begin{array}{cc}\hfill \hfill & \hfill \hfill \\ \multicolumn{1}{c}{}& \hfill \hfill \end{array}\hfill \\ \multicolumn{1}{c}{\begin{array}{cc}\hfill 0\hfill & \hfill k_{y,\theta y}\hfill \\ \multicolumn{1}{c}{k_{x,\theta z}}& \hfill 0\hfill \end{array}}& \hfill \begin{array}{cc}\hfill 0\hfill & \hfill 0\hfill \\ \multicolumn{1}{c}{k_{z,\theta z}}& \hfill k_{\theta x,\theta z}\hfill \end{array}\hfill & \hfill \begin{array}{cc}\hfill k_{\theta y,\theta y}\hfill & \hfill \hfill \\ \multicolumn{1}{c}{0}& \hfill k_{\theta z,\theta z}\hfill \end{array}\hfill \end{array}\right]$$

(19)

Maximum torque-induced bending deformation in ball screw systems

From the stiffness of the screw–nut interface obtained in section “Screw–nut interface stiffness model,” the conclusion can be drawn that deformation not only in torsional direction but also in axial and bending direction occurs when a torque is applied to one end of the screw shaft.

In this research, the situation is considered that the nut has one loaded circuit. This means that all the balls have the same entry and exit angles. Consequently, ${\phi }_{\mathrm{end}}=-{\phi }_{\mathrm{st}}=\mathrm{\Phi }$ .

The number of balls per unit angular rotation can be expressed as

$$c_{\mathrm{ball}}=\frac{N_{\mathrm{ball}}}{2\mathrm{\Phi }}$$

(20)

From the coupling stiffness matrix obtained in section “Screw–nut interface stiffness model,” the third and sixth rows of equation (19) are linearly dependent; the sixth row is equal to the third row scaled by a factor of − $r_{\mathrm{g}}$ .

$k_{x,z}$ can be expressed as

$$k_{x,z}=-2c_{\mathrm{ball}}{k}_{\mathrm{ball}}\sin 2\mathrm{\alpha }{\cos }^2\beta \sin \mathrm{\Phi }$$

(21)

And $k_{z,\theta x}$ can be expressed as

$$k_{z,\theta x}=-4c_{\mathrm{ball}}{k}_{\mathrm{ball}}R{\cos }^2\beta [{\sin }^2\alpha (\mathrm{\Phi }\cos \mathrm{\Phi }-2\sin \mathrm{\Phi })+\sin \mathrm{\Phi }]$$

(22)

In this research, ball screw system is modeled as a simple Euler–Bernoulli beam as shown in Figure 5.

Figure 5.

Model of ball screw systems.

In this model, the stiffness of bearings is considered much larger than that of others, and the ball screw is simplified as simply supported beam. The distance between the nut and the left end of the screw is $\xi L$ , where $L$ is the length of the ball screw. The value of $\xi$ is between 0 and 1. A unit torque is applied at the left end of the ball screw. Because of the coupling stiffness, the unit torque can cause the deformation in bending direction.

There is no coupling stiffness between y-direction and θ_y -direction. Therefore, the deformation in the two directions is ignored. The equilibrium of node at the position of the nut can be described as

$$\mathbf{Ku}=\mathbf{F}$$

(23)

where

$$K=\left[\begin{array}{cc}\begin{array}{cc}{k}_{x,x}+k_{bx}& \\ k_{x,z}& k_{z,z}+k_{bz}\end{array}& \begin{array}{cc}& \\ sym& \end{array}\\ \begin{array}{cc}{k}_{x,\theta x}& k_{z,\theta x}\\ -r_g{k}_{x,z}& -r_g{k}_{z,z}\end{array}& \begin{array}{cc}{k}_{\theta x,\theta x}+k_{b\theta x}& \\ -r_g{k}_{z,\theta x}& r_g{}^2{k}_{z,z}\end{array}\end{array}\right];\mathrm{\hspace{1em}}u={\left\{\begin{array}{cc}\begin{array}{cc}{u}_x& u_z\end{array}& \begin{array}{cc}{\theta }_x& {\theta }_z\end{array}\end{array}\right\}}^T;\mathrm{\hspace{1em}}F={\left\{\begin{array}{cc}\begin{array}{cc}0& 0\end{array}& \begin{array}{cc}0& 1\end{array}\end{array}\right\}}^T$$

(23)

In the stiffness matrix, $k_{\mathrm{bx}}$ , $k_{\mathrm{bz}}$ , and $k_{b\theta x}$ are the stiffness of the beam in the position of the nut. They can be expressed as

$$k_{\mathrm{bx}}=\frac{3\mathrm{EI}}{L^3{(1-\xi )}^2{\xi }^2}$$

(24a)

$$k_{\mathrm{bz}}=\frac{\mathrm{EA}}{L\xi }$$

(24b)

$$k_{b\theta x}=\frac{3\mathrm{EI}}{L(1-3\xi +3{\xi }^2)}$$

(24c)

where E, A, and I are the elastic modulus of the material, cross-sectional area of the ball screw shaft, and moment of inertia, respectively.

The relationship between maximum bending deformation of the system ( ${\mathit{u}}_{\mathit{x},\max }$ and ${\mathit{u}}_{\mathit{y},\mathit{max}}$ ) and the deformation of the beam on nut position is

$${\mathit{u}}_{\mathit{x},\mathit{max}}={\mathit{W}}_{\mathit{x}}{\mathit{u}}_{\mathit{x}};{\mathit{u}}_{\mathit{y},\mathit{max}}={\mathit{W}}_{\mathit{y}}L{\mathbf{\theta }}_{\mathit{x}}$$

(25)

where $W_x$ and $W_y$ are piecewise functions of $\xi$

$${\mathit{W}}_{\mathit{x}}=\{\begin{array}{c}\hfill \sqrt{\frac{{(1+\mathit{\xi })}^3}{27(1-\mathit{\xi }){\mathit{\xi }}^2}}0\le \mathit{\xi }\le 0.5\hfill \\ \multicolumn{1}{c}{\sqrt{\frac{{(2-\mathit{\xi })}^3}{27{(1-\mathit{\xi })}^2\mathit{\xi }}}0.5\le \mathit{\xi }\le 1}\end{array}$$

(26a)

$${\mathit{W}}_{\mathit{y}}=\{\begin{array}{c}\hfill \frac{1}{3\sqrt{3}}\frac{{(1-3{\mathit{\xi }}^2)}^{3/2}}{1-3\mathit{\xi }+3{\mathit{\xi }}^2}0\le \mathit{\xi }\le 0.5\hfill \\ \multicolumn{1}{c}{\frac{1}{3\sqrt{3}}\frac{{[3\mathit{\xi }(2-\mathit{\xi })-2]}^{3/2}}{1-3\mathit{\xi }+3{\mathit{\xi }}^2}0.5\le \mathit{\xi }\le 1}\end{array}$$

(26b)

Dimensionless maximum bending deformation of the system among different nut positions is defined as

$$J=\underset{0.1\le \mathrm{\xi }\le 0.9}{max}[{\left(\frac{u_{\mathit{x},\max }}{r_g{\theta }_z}\right)}^2+{\left(\frac{u_{\mathit{y},\max }}{r_g{\theta }_z}\right)}^2]$$

(27)

Uncertainty analysis of torque-induced bending deformations in ball screw systems

Based on the model in section “Maximum torque-induced bending deformation in ball screw systems,” the stiffness matrix $\mathit{K}$ can be expressed as a function of the parameter vector $\mathit{b}=[k_{\mathrm{ax}},\beta ]$

$$\mathit{K}=\mathit{K}(\mathit{b})$$

(28)

The parameters of $\mathit{b}$ are usually assumed to be deterministic for identical ball screw systems. However, these parameters may be uncertain due to wear, preload degradation, and installation errors. Therefore, there are uncertainties in the parameter vector $\mathit{b}$ . The thought of interval analysis^9,10 is followed, and it is assumed that the parameter vector $\mathit{b}=[k_{\mathrm{ax}},\beta ]$ belongs to a interval vector ${\mathit{b}}^I=[b_1^I,b_2^I]$ . The elements of the vector $\mathit{b}$ hold

$$b_j\in b_j^I=[\underset{¯}{b_j},\overline{b_j}],j=1,2$$

(29)

where $b_j^I$ is an interval number, and $\underset{¯}{b_j},\overline{b_j}$ are the lower bound and upper bound of $b_j^I$ , respectively. Then the linear algebraic equations of static response of structures with uncertain parameters may be rewritten as follows

$$\mathit{K}(\mathit{b})\mathit{x}(\mathit{b})=\mathit{F}$$

(30)

As the uncertainties in $\mathit{K}(\mathit{b})$ are caused by the interval parameter $\mathit{b}$ , then the elements of $\mathit{K}(\mathit{b})$ vary within a deterministic interval. And they can be written in the following inequality form

$$\min (\mathit{K}(\mathit{b}))=\underset{¯}{\mathit{K}}\le \mathit{K}(\mathit{b})\le \overline{\mathbf{K}}=\max (\mathit{K}(\mathit{b}))$$

(31)

Thus, the deformation $\mathit{x}(\mathit{b})$ of equation (30), which is a function of parameter vector $\mathit{b}$ , also varies within a deterministic interval $\mathit{x}{}^I(\mathit{b})=[\underset{¯}{\mathit{x}}(\mathit{b}),\overline{\mathit{x}}(\mathit{b})]$ , where

$$\underset{¯}{\mathit{x}}(\mathit{b})=\mathrm{Gmin}\{\mathit{x}:\mathit{K}(\mathit{b})\mathit{x}(\mathit{b})=\mathit{f}(\mathit{b})\},\mathit{b}\in {\mathit{b}}^I$$

(32a)

$$\overline{\mathit{x}}(\mathit{b})=\mathrm{Gmax}\{\mathit{x}:\mathit{K}(\mathit{b})\mathit{x}(\mathit{b})=\mathit{f}(\mathit{b})\},\mathit{b}\in {\mathit{b}}^I$$

(32b)

and $\mathrm{Gmax}\{\mathit{x}(\mathit{b})\}$ and $\mathrm{Gmin}\{\mathit{x}(\mathit{b})\}$ are the global maximum value and the global minimum value of $\mathit{x}(\mathit{b})$ , respectively.

Based on interval mathematics, the nominal value vector of the interval parameter vector ${\mathit{b}}^I=[b_1^I,b_2^I]$ is denoted by

$${\mathit{b}}^c=[b_1^c,b_2^c]$$

(33)

where $b_j^c=0.5({\underset{¯}{b}}_j+{\overline{b}}_j),j=1,2$ . And the uncertain radius vector of the interval parameter vector ${\mathit{b}}^I$ is denoted by

$$\Delta \mathit{b}=(\Delta b_1,\Delta b_2)$$

(34)

where $\Delta b_j=0.5({\overline{b}}_j-{\underset{¯}{b}}_j),j=1,2$ . Therefore, the parameter vector $\mathit{b}$ can be decomposed into the sum of the nominal value vector ${\mathit{b}}^c=[b_1^c,b_2^c]$ and the deviation vector $\mathit{\delta }=[{\delta }_1,{\delta }_2]$

$$\mathit{b}={\mathit{b}}^c+\mathit{\delta }$$

(35)

where ${\delta }_j\in {\delta }_j^I=[-\Delta b_j,\Delta b_j]$ is the uncertainty of $b_j$ . Based on the assumption that there are fairly large uncertainties in the parameters and using the second-order Taylor series expansion and expanding the response $\mathit{x}(\mathit{b})$ about ${\mathit{b}}^c=[b_1^c,b_2^c]$

$$\mathit{x}(\mathit{b})=\mathit{x}({\mathit{b}}^c)+{\displaystyle \sum _{j=1}^2}\frac{\partial \mathit{x}({\mathit{b}}^c)}{\partial b_j}{\delta }_j+\frac{1}{2}{\displaystyle \sum _{j=1}^2}{\displaystyle \sum _{k=1}^2}\frac{\partial \mathit{x}({\mathit{b}}^c)}{\partial b_j\partial b_k}{\delta }_j{\delta }_k$$

(36)

where $\mathit{x}({\mathit{b}}^c)$ can be obtained by substituting $\mathit{b}={\mathit{b}}^c$ into equation (30). The first-order partial derivative $\partial \mathit{x}({\mathit{b}}^c)/\partial b_j$ and the second-order partial derivative ${\partial }^2\mathit{x}({\mathit{b}}^c)/(\partial b_j\partial b_k)$ can be derived as follows. Differentiating both sides of equation (30) with respect to $b_j$ yielded

$$\mathit{K}(\mathit{b})\frac{\partial \mathit{x}(\mathit{b})}{\partial b_j}=-\frac{\partial \mathit{K}(\mathit{b})}{\partial b_j}\mathit{x}(\mathit{b})$$

(37)

Differentiating both sides of equation (36) with respect to $b_k$ yielded

$$\frac{\partial \mathit{K}(\mathit{b})}{\partial b_k}\frac{\partial \mathit{x}(\mathit{b})}{\partial b_j}+\mathit{K}(\mathit{b})\frac{{\partial }^2\mathit{x}(\mathit{b})}{\partial b_j\partial b_k}=-\frac{\partial \mathit{K}(\mathit{b})}{\partial b_j}\frac{\partial \mathit{x}(\mathit{b})}{\partial b_k}-\frac{{\partial }^2\mathit{K}(\mathit{b})}{\partial b_j\partial b_k}\mathit{x}$$

(38)

Equation (37) can be rewritten as

$$\mathit{K}(\mathit{b})\frac{{\partial }^2\mathit{x}(\mathit{b})}{\partial b_j\partial b_k}=-\frac{\partial \mathit{K}(\mathit{b})}{\partial b_j}\frac{\partial \mathit{x}(\mathit{b})}{\partial b_k}-\frac{{\partial }^2\mathit{K}(\mathit{b})}{\partial b_j\partial b_k}-\frac{\partial \mathit{K}(\mathit{b})}{\partial b_k}\frac{\partial \mathit{x}(\mathit{b})}{\partial b_j}$$

(39)

Substituting $\mathit{b}={\mathit{b}}^c$ into equations (37) and (39), $\partial \mathit{x}({\mathit{b}}^c)/\partial b_j$ and ${\partial }^2\mathit{x}({\mathit{b}}^c)/(\partial b_j\partial b_k)$ can be obtained by solving these algebraic linear equations.

When considering the static response analysis of ball screw systems with interval parameters, the primary interest lies in the bounds of the response. Making use of optimization method, from equation (36), the lower and upper bounds of the static response (30) can be determined as follows

$$\begin{array}{cc}\multicolumn{1}{c}{\underset{¯}{\mathit{x}}=}& \mathrm{Gmin}(\mathit{x}(\mathit{b}))=\mathit{x}({\mathit{b}}^c)+\mathrm{Gmin}\hfill \\ \multicolumn{1}{c}{}& ({\displaystyle \sum _{j=1}^2}\frac{\partial x_i({\mathit{b}}^c)}{\partial b_j}{\delta }_j+\frac{1}{2}{\displaystyle \sum _{j=1}^2}{\displaystyle \sum _{k=1}^2}\frac{{\partial }^2{x}_i({\mathit{b}}^c)}{\partial b_j\partial b_k}{\delta }_j{\delta }_k)\hfill \end{array}$$

(40a)

$$\overline{\mathit{x}}=\mathrm{Gmax}(\mathit{x}(\mathit{b}))=\mathit{x}({\mathit{b}}^c)+\mathrm{Gmax}({\displaystyle \sum _{j=1}^2}\frac{\partial x_i({\mathit{b}}^c)}{\partial b_j}{\delta }_j+\frac{1}{2}{\displaystyle \sum _{j=1}^2}{\displaystyle \sum _{k=1}^2}\frac{{\partial }^2{x}_i({\mathit{b}}^c)}{\partial b_j\partial b_k}{\delta }_j{\delta }_k)$$

(40b)

The optimization problems can be transformed into the following QB problems

$$\begin{array}{cc}\multicolumn{1}{c}{}& \mathrm{Gmin}(\mathrm{Gmax})f(\mathit{Y})=0.5〈\mathit{Y},\mathit{AY}〉+〈\mathit{d},\mathit{Y}〉\hfill \\ \multicolumn{1}{c}{}& s.t.-\Delta b_j\le y_j\le \Delta b_j,j=1,2\hfill \end{array}$$

(41)

where $\mathit{d}=[(\partial x_i({\mathit{b}}^c))/(\partial b_1),(\partial x_i({\mathit{b}}^c))/\partial b_2]$ is a two-dimensional vector, and $\mathit{A}=(a_{\mathrm{jk}})$ is a real symmetry matrix, and the value of element $a_{\mathrm{jk}}$ is $({\partial }^2{x}_i({\mathit{b}}^c))/(\partial b_j\partial b_k),j,k=1,2$ .

The elements of $\mathit{A}$ are problem dependent, and therefore, it is difficult to check the positive semi-definiteness of the matrix $\mathit{A}$ . The optimization problems (41) are two nonconvex QB problems. Due to its computational complexity, the method of DCA ¹¹ is conducted to solve these QB problems.

DCA proposed by Tao and Souad ¹² is a continuous optimization approach based on the duality for solving programming in a general form of difference of convex functions, and the DCA only requires matrix-vector products. For the QB problems (41), the DCA converges to a global optimal solution in polynomial time.

DCA for the QB problems or problems (41) is well established and described in detail in the literature. ¹¹ The QB problem can be written in the following form

$$\begin{array}{cc}\multicolumn{1}{c}{}& (\mathrm{QB})\mathrm{Gmin}f(\mathit{Y})=0.5〈\mathit{Y},\mathit{AY}〉+〈\mathit{d},\mathit{Y}〉\hfill \\ \multicolumn{1}{c}{}& s.t.l_i\le {\mathrm{y}}_i\le u_i,i=1,2\hfill \end{array}$$

(42)

where $\mathit{d}$ is a four-dimensional vector, and $\mathit{A}$ is just a real symmetry matrix; $l_i,u_i$ are the bounds of design variable $y_i\in \mathit{Y}$ , and $\mathit{Y}$ is a four-dimensional vector. Because the maximization problem can be transformed into minimization problem, the DCA also can obtain the global maximum solution.

The circumscribed ball of design region of the ball constrained quadratic programming (QBA) problem can be described by

$$\mathit{C}=\{\mathit{Y}\in R^2:\Vert \mathit{Y}-{\mathit{c}}_e\Vert \le r\}$$

(43)

where ${\mathit{c}}_e=[c_1,c_2]$ is the center of the circumscribed ball of the box region $\mathbf{D}=\{\mathit{Y}\in R^2:u_i{<\mathrm{y}}_i<l_i,i=1,2\}$ , and $c_i=0.5(l_i+u_i),i=1,2$ , and the radius of this ball is $r=0.5\sqrt{{\displaystyle \sum _{i=1}^2}{(u_i-l_i)}^2}$ . Substituting equation $\mathit{Y}={\mathit{c}}_e+\mathit{Y}\prime$ into equation (42)

$$\begin{array}{cc}\multicolumn{1}{c}{}& minf(\mathit{Y}\prime )=0.5〈\mathit{Y}\prime ,\mathit{QY}\prime 〉+〈\mathit{q},\mathit{Y}\prime 〉\hfill \\ \multicolumn{1}{c}{}& +〈{\mathit{c}}_e,0.5{\mathit{Ac}}_e+\mathit{d}〉\hfill \\ \multicolumn{1}{c}{}& s.t.\Vert \mathit{Y}\prime \Vert \le r\hfill \end{array}$$

(44)

where $\mathit{Q}=\mathit{A}$ and $\mathit{q}={\mathit{Ac}}_e+\mathit{d}$ .

By solving the QB problems (41) based on the DCA for QB problems, we obtained

$${\underset{¯}{x}}_i=x_i({\mathit{b}}^c)+\mathrm{Gmin}(f(\mathit{Y})),i=1,2,\dots ,4$$

(45a)

$${\overline{x}}_i=x_i({\mathit{b}}^c)+\mathrm{Gmax}(f(\mathit{Y})),i=1,2,\dots ,4$$

(45b)

where $\mathrm{Gmin}(f(\mathit{Y}))$ and $\mathrm{Gmax}(f(\mathit{Y}))$ are the global minimum value and the global maximum value of optimization formulation (41), respectively.

Validation

In order to validate the effectiveness of the analysis, MCS method is conducted and compared with the method in this article.

• Case 1. A lead of 20 mm is used for the simulations. This corresponds to $\alpha =17.7\circ$ , ${\phi }_{\mathrm{end}}=-{\phi }_{\mathrm{st}}=162.5\circ$

• Case 2. A lead of 5 mm is used for the simulations. This corresponds to $\alpha =4.5\circ$ , ${\phi }_{\mathrm{end}}=-{\phi }_{\mathrm{st}}=180\circ$

Both the two cases have the same uncertainty. $k_{\mathrm{ax}}$ is considered an interval parameter, and it is expressed by $k_{\mathrm{ax}}\in k_{\mathrm{ax}}^I=[1-u_k,1+u_k]k_{\mathrm{axc}}$ , where $k_{\mathrm{axc}}=3.43\times {10}^8\mathrm{N}/\mathrm{m}$ ; $u_k$ is the uncertain factor, and its value is 15%. Contact angle $\beta$ is another interval parameter, and it is expressed as $\beta \in {\beta }^I=[1-u_{\beta },1+u_{\beta }]{\beta }_c$ , where ${\beta }_c=75\circ$ ; and $u_{\beta }$ is the uncertain factor; its value is 2%. The result of second-order Taylor series expansion and first-order Taylor series expansion ¹³ is compared with that of MCS method, which is listed in Tables 1 and 2.

Table 1.

Comparison among different methods for Case 1.

Method	Upper bound	Error (%)	Lower bound	Error (%)
Second-order Taylor series expansion	1.58e−004	0.637	7.79e−005	0.129
First-order Taylor series expansion	1.53e−004	2.548	7.37e−005	5.270
Monte Carlo simulation method	1.57e−004		7.78e−005

Table 2.

Comparison among different methods for Case 2.

Method	Upper bound	Error (%)	Lower bound	Error (%)
Second-order Taylor series expansion	9.13e−003	0.220	6.36e−003	0.469
First-order Taylor series expansion	9.13e−003	0.220	6.34e−003	0.782
Monte Carlo simulation method	9.11e−003		6.39e−003

MCS method is the numerical method for estimating the upper and lower bounds in a box region $\mathit{D}=\{(k_{\mathrm{ax},}\beta )\in R^2:{\underset{¯}{k}}_{\mathrm{ax}<}{k}_{\mathrm{ax}}<{\overline{k}}_{\mathrm{ax}},\underset{¯}{\beta }<\beta <\overline{\beta }\}$ . The sample for computation is randomly chosen. Therefore, the result obtained by MCS is accurate enough if the sample is adequate, and the result obtained by other methods is compared with MCS.

From the comparison with MCS method in Tables 1 and 2, the errors of the two methods are relatively small, and the result of second-order Taylor series expansion is more accurate than that of first-order Taylor series expansion. The method in this research can be employed to conduct the uncertainty analysis of torque-induced bending deformation in ball screw systems.

Results and discussions

Detail uncertainty analysis for the two cases above is conducted. From the bounds of the maximum torque-induced bending deformations in ball screw systems depicted in Figure 6, the conclusion can be drawn that the upper and lower bounds of maximum deformation of ball screw systems studied in this article are approximately proportional to the uncertain factor $u_k$ and $u_{\beta }$ . Thus, the bounds of maximum deformation can be fitted by Figure 6 without calculating.

Figure 6.

Bounds of the maximum torque-induced bending deformations in ball screw systems: (a) Case 1 ( $\mathit{D}=\{(u_k,u_{\beta })\in R^2:0\le u_k\le 0.15,0\le u_{\beta }\le 0.02\}$ ), (b) Case 2 ( $\mathit{D}=\{(u_k,u_{\beta })\in R^2:0\le u_k\le 0.15,0\le u_{\beta }\le 0.02\}$ ), (c) Case 1 (u_β  = 0), (d) Case 1 (u_k  = 0), (e) Case 2 (u_β  = 0), and (f) Case 2 (u_k  = 0).

From the comparison in Figure 6, the larger the uncertain factors, the larger the errors between the results of second-order Taylor series expansion and that of first-order Taylor series expansion. Therefore, much larger error will be introduced if first-order Taylor series expansion is applied to analyze the bounds of the maximum deformation in ball screw systems.

According to the uncertain analysis for the two cases in Figure 6, the conclusion can be drawn that the bound of torque-induced bending deformations in ball screw systems can be obtained via second-order Taylor series expansion. If the uncertainty of $k_{\mathrm{ax}}$ and $\beta$ can be confirmed, the maximum torque-induced bending deformations in the system can be estimated.

Conclusion

The following conclusions can be drawn from the uncertainty analysis of torque-induced bending deformations in ball screw systems based on an numerical method combined second-order Taylor series expansion with DCA:

1. For ball screw system studied in this article, the upper and lower bounds of maximum deformation of ball screw systems are approximately proportional to the uncertain factor $u_k$ and $u_{\beta }$ . Thus, the bounds of maximum deformation can be fitted by the result obtained in this research without calculating.

2. By comparison among the results obtained by second-order Taylor series expansion, first-order Taylor series expansion, and MCS, it is concluded that second-order Taylor series expansion has high precision in predicting the bounds of the maximum torque-induced bending deformations in ball screw systems. Uncertain analysis of maximum deformation in ball screw systems can be conducted by the method combining second-order Taylor series expansion with DCA, which provides a new approach for estimating the maximum torque-induced bending deformations in ball screw systems.