## 3 Linear Stability Analysis  
The dynamic equations of the drilling system (Equations (14)-(17)) have a stable solution, i.e., \(z(\tau)\ =\phi(\tau)=0\) (drilling without vibration), which is the equilibrium point of the system. The stability of the equilibrium point is analyzed in this section. The linearized homogeneous state vector equation can be obtained from Equations (30) and (31) as:  
\[\dot{\mathbf{x}}(\tau)-\mathbf{A}\mathbf{x}(\tau)-\mathbf{B}\mathbf{x}(\tau-\tau_{n0})=\mathbf{0} \tag{32}\]  
This is a delayed differential equation, and its characteristic equation is:  
\[\Big{|}\lambda\mathbf{I}-\mathbf{A}-e^{-\lambda\tau_{n0}}\mathbf{B}\Big{|}=\mathbf{0}, \tag{33}\]  
where \(\lambda\) denotes an eigenvalue of the linearized system and the exponential term appears due to the time delay. There are an infinite number of eigenvalues in Equation (33). The drill string system is stable only when all the eigenvalues have negative real part \((Re(\lambda)<0)\), otherwise the system is unstable. Pure imaginary eigenvalues \((Re(\lambda)=0)\) corresponding to a specific condition determine the stability boundaries that divide the stable and unstable regions. In this case, the responses of the drilling system exhibit periodic oscillations with frequency \(\omega_{c}\). \[\mathbf{X}=[\lambda,1,\lambda k,k]^{\mathrm{T}} \tag{34}\]  
\[k_{SDD}=\frac{\lambda^{2}+v^{-1}Q_{A}^{-1}\lambda+v^{-2}}{\rho\psi k_{\tau}v^{ -2}\big{(}1-e^{-\lambda\tau_{n0}}\big{)}}+\rho^{-1} \tag{35}\]  
\[k_{CD}=\frac{-\psi\big{(}1-e^{-\lambda\tau_{n0}}\big{)}}{\lambda^{2}+Q_{T}^{-1 }\lambda+1} \tag{36}\]  
Assuming that the eigenvalue \(\lambda\) is equal to \(\omega_{c}\ (Re(\lambda)=0)\), the characteristic Equation (33) could be converted into a frequency equation, and the following four equations can be obtained by separating the real part and imaginary part. Using \(\omega_{c}\) as the parameter, the stability diagrams of the SDD model and CD model, respectively, could be drawn (see Figure 5). The abscissa is the non-dimensional rotational speed \(\Omega\) and the ordinate is the non-dimensional axial cutting depth \(\psi\). \[\Omega_{SDD}=\frac{\pi\omega_{c}}{N\Big{[}k\pi+\arctan\Big{(}\frac{C_{SDD}C_{ 1}+B_{SDD}C_{2}}{C_{SDD}C_{2}-B_{SDD}C_{1}}\Big{)}\Big{]}} \tag{37}\]\[\psi_{SDD}=\frac{C_{1}}{C_{SDD}[1-\cos(\tau_{n0}\omega_{c})]-B_{SDD}\sin(\tau_{n0} \omega_{c})} \tag{38}\]  
\[\Omega_{CD}=\frac{\pi\omega_{c}}{N\Big{[}k\pi+\arctan\Big{(}\frac{C_{CD}C_{1}+B _{CD}C_{2}}{C_{CD}C_{2}-B_{CD}C_{1}}\Big{)}\Big{]}} \tag{39}\]  
\[\psi_{CD}=\frac{C_{1}}{C_{CD}[1-\cos(\tau_{n0}\omega_{c})]-B_{CD}\sin(\tau_{n0 }\omega_{c})}, \tag{40}\]  
where \(k\) = 1, 2, 3,..., and represents the number of lobes. \[\left\{ \begin{array}{l}B_{SDD}=k_{r}v^{-2}\omega_{c}Q_{T}^{-1}+\rho v^{-1}Q _{A}^{-1}\omega_{c},B_{CD}=k_{r}v^{-2}\omega_{c}Q_{T}^{-1}\\ \\ C_{SDD}=k_{r}v^{-2}-\rho v^{-2}-k_{r}v^{-2}\omega_{c}^{2}+\rho\omega_{c}^{2},C _{CD}=k_{r}v^{-2}-k_{r}v^{-2}\omega_{c}^{2}\\ \\ C_{1}=-v^{-2}+v^{-1}\omega_{c}^{2}Q_{T}^{-1}Q_{A}^{-1}+v^{-2}\omega_{c}^{2}+ \omega_{c}^{2}-\omega_{c}^{4}\\ \\ C_{2}=-v^{-2}Q_{T}^{-1}\omega_{c}-v^{-1}Q_{A}^{-1}\omega_{c}+Q_{T}^{-1}\omega_ {c}^{3}+v^{-1}Q_{A}^{-1}\omega_{c}^{3}\end{array} \right. \tag{41}\]  
There are seven lobes in the Figure 5; the upper part of the lobe is the instability of the drilling system, the lower part of the lobe is the stable interval, and the lobe represents the critical cutting depth. In the case of drilling, the appropriate control parameters \(\Omega\) and \(\rho\) can be chosen through the stability lobes to realize a working condition without stick-slip vibration and improve drilling efficiency. From the figure, it can be seen that the system stability interval increases with non-dimensional rotational speed. On the other hand, the stability boundaries are higher for the state-dependent delay model than for the constant delay model. For the SDD model, the system stability interval increases with the control parameters \(\rho\).