CN103823787B - Multi-turning-tool parallel turning stability judgment method based on differential quadrature method - Google Patents

Abstract

The invention provides a multi-turning-tool parallel turning stability judgment method based on a differential quadrature method. The multi-turning-tool parallel turning stability judgment method comprises the steps that dynamics modeling is carried out on a multi-turning-tool parallel turning system, and a multi-delay second order differential equation is built; a normalized state space equation is built and obtained; the second kind of Chebyshev points is adopted as discrete points in adjacent unit intervals (0, 1) and (-1, 0); speed items are expressed by displacement items at the discrete points based on a Lagrange's interpolation function by using the differential quadrature method; the intervals where the discrete points of delay items are located are judged, and the delay items are expressed by the second kind of Chebyshev points of the intervals where the discrete points of the delay items are located; state transfer matrixes between the adjacent unit intervals are constructed, and the stability of an original system is judged based on the Floquet theory. Compared with traditional single-turning-tool turning, the dynamic characteristic of the multi-turning-tool parallel turning system is analyzed by adopting the differential quadrature method, optimized cutting parameters are obtained, and machining efficiency is greatly improved.

Description

Multi-lathe parallel turning stability judgment method based on differential quadrature method

Technical Field

The invention relates to a novel Method for judging the stability of multi-lathe parallel turning, namely a Differential Quadrature Method, in particular to a Method for reasonably selecting cutting parameters by utilizing the Differential Quadrature Method to avoid the influence of regenerative chatter on multi-lathe parallel turning and carry out high-efficiency high-quality machining.

Background

In the field of machining, turning is one of the most common and common machining methods. The traditional single lathe tool machining research is mature, and the machining mode of parallel turning of a plurality of lathe tools is a new concept which is only proposed by academia in recent years, particularly the research on the stability of parallel turning of a plurality of lathe tools is still in a starting stage. Theoretically, the processing efficiency of the multi-lathe tool parallel turning is far higher than that of the traditional turning, but the processing mechanism of the multi-lathe tool parallel turning is much more complicated than that of the traditional turning. The main factor influencing the processing quality of the multi-lathe tool parallel turning is the processing stability, so that the dynamic mechanism of the multi-lathe tool parallel turning is researched, the processing parameters are reasonably selected, and the occurrence of regenerative chatter vibration can be effectively avoided, so that the stable operation of the processing process is realized, and the high-efficiency turning of a workpiece is realized while the processing quality is ensured.

Document 1 "e.budak, e.ozturn, Dynamics and stability of parallel turning-Manufacturing Technology 60(2011) 383-. The multiple cutting tools in the literature are fixed on different tool holders, and no coupling effect exists between the cutting tools. The method comprises the following steps:

(1) establishing a dynamic model for parallel turning of multiple lathes;

(2) obtaining a kinetic equation at a stable boundary;

(3) obtaining a stable region of the multi-lathe parallel turning in a cutting parameter space by using a search method;

(4) and performing time domain simulation verification and cutting experiment verification.

Document 2 "E.Ozturk, E.Budak, Modeling dynamics of parallel turning operations. proceedings of 4^thCIRP International Conference on high performance Cutting,2010 "discusses two cases of parallel turning with multiple tools: a plurality of turning tools are arranged on different tool rests; a plurality of turning tools are arranged on the same tool rest. The analytical calculation methods in these two cases are similar and the same as the procedure of document 1.

The frequency domain method has simple steps and high calculation speed, but is not beneficial to considering various complex working conditions. The time domain simulation method has visual results, but the calculation speed is low, the processing stability under the condition of single cutting parameter combination can be obtained, and the stability map of the processing parameter space is not favorable to be drawn. Compared with the method, the differential integration method proposed by Bellman and the like in the seventies of the twentieth century has the advantages of high calculation speed, high result precision and the like, and is widely applied to various engineering technical fields. The dynamic equation of the multi-lathe-cutter parallel turning is a multi-time-lag differential equation, and the differential quadrature method is properly popularized and used for judging the processing stability of the multi-lathe-cutter parallel turning, so that the method has important production value and significance.

Disclosure of Invention

Aiming at the defects in the prior art, the invention aims to provide a novel method for judging the stability of the cutting process, namely a differential quadrature method, according to the dynamic characteristics of the multi-lathe tool parallel turning, so as to provide an effective basis for selecting the multi-lathe tool parallel turning processing parameters, obtain the highest processing efficiency on the premise of ensuring the high-quality processing without regenerative chatter and generate good economic benefits.

The invention provides a method for judging the stability of multi-lathe parallel turning based on a differential quadrature method, which comprises the following steps:

and performing dynamic modeling on the multi-lathe parallel turning system, and establishing a multi-time-lag second-order differential equation. Carrying out state space change, establishing a state space equation, and then carrying out normalization processing to obtain a normalized standard equation;

equivalently dispersing the standard equation into a group of algebraic equations by taking Chebyshev-Gauss-Lobatto Points as discrete Points on adjacent unit intervals [0,1] and [ -1,0 ];

representing a speed term by using a displacement term at a discrete point based on a Lagrange interpolation function by using a differential quadrature method; judging the interval where the discrete point of the time lag term is located, and representing the time lag term by using a second class Chebyshev point of the interval where the discrete point of the time lag term is located based on a Lagrange interpolation function;

constructing a state transition matrix between two adjacent unit intervals [0,1] and [ -1,0], and judging the stability of the original system according to the Floquet theory; if the moduli of all the eigenvalues of the state transition matrix are less than 1, the system is stable; if the modulus of any eigenvalue of the state transition matrix is larger than 1, the system is unstable; therefore, the stability spectrogram of the multi-lathe parallel turning system in the cutting parameter space can be obtained.

Specifically, the method for judging the stability of the multi-lathe parallel turning based on the differential quadrature method provided by the invention comprises the following steps:

step 1: establishing a multi-lathe tool parallel turning kinetic equation, and sorting the multi-lathe tool parallel turning kinetic equation to obtain a multi-time-lag second-order differential equation;

step 2: carrying out state space transformation on the multi-time-lag second-order differential equation to obtain a state space equation;

and step 3: carrying out normalization processing on the state space equation to obtain a state space equation in a standard form;

and 4, step 4: carrying out periodic dispersion on a state space equation in a standard form, and equivalently converting the state space equation into a group of algebraic equations serving as a state space equation expression;

and 5: based on a Lagrange interpolation function, taking a second class Chebyshev point on a unit interval [0,1] as a discrete point, and expressing a derivative term in a state space equation expression by using a displacement term by using a differential integration method;

step 6: judging the interval where the time lag term in the state space equation expression is located based on a Lagrange interpolation function, and if the discrete point of the time lag term belongs to the interval [0,1], representing the time lag term by taking the second class Chebyshev point on [0,1] as the discrete point; if the discrete point of the time lag term belongs to the interval [ -1,0], the time lag term is represented by taking the second type Chebyshev point on [ -1,0] as the discrete point;

and 7: and constructing a state transition matrix between adjacent unit intervals of [0,1] and [ -1,0], and judging the stability of the original system according to the Floquet theory.

Preferably, in step 1, the equation of dynamics of the multi-lathe tool parallel turning is shown in formula (1):

$\begin{array}{c}\left[\begin{array}{c}{\stackrel{\·\·}{z}}_1\left(t\right)\\ {\stackrel{\·\·}{z}}_2\left(t\right)\end{array}\right]+\left[\begin{array}{cc}2{\mathrm{\ζ}}_1{\mathrm{\ω}}_{n1}& 0\\ 0& 2{\mathrm{\ζ}}_2{\mathrm{\ω}}_{n2}\end{array}\right]\left[\begin{array}{c}{\stackrel{\·}{z}}_1\left(t\right)\\ {\stackrel{\·}{z}}_2\left(t\right)\end{array}\right]+\left[\begin{array}{cc}{\mathrm{\ω}}_{n1}^2& 0\\ 0& {\mathrm{\ω}}_{n2}^2\end{array}\right]\left[\begin{array}{c}{z}_1\left(t\right)\\ z_2\left(t\right)\end{array}\right]=\\ K_f\left[\begin{array}{c}\frac{{\mathrm{\ω}}_{n1}^2}{k_1}{a}_1(-z_1(t)+z_2(t-\frac{\mathrm{\τ}}{2}\left)\right)\\ \frac{{\mathrm{\ω}}_{n2}^2}{k_2}{a}_1(-z_2(t)+z_1(t-\frac{\mathrm{\τ}}{2}\left)\right)+\frac{{\mathrm{\ω}}_{n2}^2}{k_2}(a_2-a_1)(-z_2(t)+z_2(t-\mathrm{\τ}\left)\right)\end{array}\right]\end{array}---\left(1\right)$

wherein subscript 1 represents a first turning tool, subscript 2 represents a second turning tool, z (t) represents a kinematically responsive displacement, ω_nRepresenting natural circular frequency, K_fDenotes a cutting force coefficient in a feed direction, τ denotes a workpiece rotation period, ζ denotes a relative damping ratio, k denotes a stiffness coefficient, a₁Indicating the depth of cut of the first turning tool, a₂Indicating the depth of cut of the second turning tool,it is indicative of the dynamic acceleration response,representing the dynamic velocity response, t represents time;

the multi-lag second-order differential equation is shown in formula (2):

$\begin{array}{c}\left[\begin{array}{c}{\stackrel{\·\·}{z}}_1\left(t\right)\\ {\stackrel{\·\·}{z}}_2\left(t\right)\end{array}\right]+\left[\begin{array}{cc}2{\mathrm{\ζ}}_1{\mathrm{\ω}}_{1n}& 0\\ 0& 2{\mathrm{\ζ}}_2{\mathrm{\ω}}_{n2}\end{array}\right]\left[\begin{array}{c}{\stackrel{\·}{z}}_1\left(t\right)\\ {\stackrel{\·}{z}}_2\left(t\right)\end{array}\right]+\left[\begin{array}{cc}{\mathrm{\ω}}_{n1}^2+K_f\frac{{\mathrm{\ω}}_{n1}^2}{k_1}{a}_1& 0\\ 0& {\mathrm{\ω}}_{n2}^2+K_f\frac{{\mathrm{\ω}}_{n2}^2}{k_2}{a}_2\end{array}\right]\left[\begin{array}{c}{z}_1\left(t\right)\\ z_2\left(t\right)\end{array}\right]=\\ \left[\begin{array}{cc}0& K_f\frac{{\mathrm{\ω}}_{n1}^2}{k_1}{a}_1\\ K_f\frac{{\mathrm{\ω}}_{n2}^2}{k_2}{a}_1& 0\end{array}\right]\left[\begin{array}{c}{z}_1(t-\frac{\mathrm{\τ}}{2})\\ z_2(t-\frac{\mathrm{\τ}}{2})\end{array}\right]+\left[\begin{array}{cc}0& 0\\ 0& K_f\frac{{\mathrm{\ω}}_{n2}^2}{k_2}(a_2-a_1)\end{array}\right]\left[\begin{array}{c}{z}_1(t-\mathrm{\τ})\\ z_2(t-\mathrm{\τ})\end{array}\right]\end{array}---\left(2\right).$

preferably, in step 2, the state space equation is as shown in equation (3):

$\begin{array}{c}\left[\begin{array}{c}{\stackrel{\·}{z}}_1\left(t\right)\\ {\stackrel{\·}{z}}_2\left(t\right)\\ {\stackrel{\·\·}{z}}_1\left(t\right)\\ {\stackrel{\·\·}{z}}_2\left(t\right)\end{array}\right]=\left[\begin{array}{cccc}0& 0& 1& 0\\ 0& 0& 0& 1\\ -{\mathrm{\ω}}_{n1}^2-K_f\frac{{\mathrm{\ω}}_{n1}^2}{k_1}{a}_1& 0& -2{\mathrm{\ζ}}_1{\mathrm{\ω}}_{n1}& 0\\ 0& -{\mathrm{\ω}}_{n2}^2-K_f\frac{{\mathrm{\ω}}_{n2}^2}{k_2}{a}_2& 0& -2{\mathrm{\ζ}}_2{\mathrm{\ω}}_{n2}\end{array}\right]\left[\begin{array}{c}{z}_1\left(t\right)\\ z_2\left(t\right)\\ {\stackrel{\·}{z}}_1\left(t\right)\\ {\stackrel{\·}{z}}_2\left(t\right)\end{array}\right]+\\ \left[\begin{array}{cccc}0& 0& 0& 0\\ 0& 0& 0& 0\\ 0& K_f\frac{{\mathrm{\ω}}_{n1}^2}{k_1}{a}_1& 0& 0\\ K_f\frac{{\mathrm{\ω}}_{n2}^2}{k_2}{a}_1& 0& 0& 0\end{array}\right]\left[\begin{array}{c}{z}_1(t-\frac{\mathrm{\τ}}{2})\\ z_2(t-\frac{\mathrm{\τ}}{2})\\ {\stackrel{\·}{z}}_1(t-\frac{\mathrm{\τ}}{2})\\ {\stackrel{\·}{z}}_2(t-\frac{\mathrm{\τ}}{2})\end{array}\right]+\left[\begin{array}{cccc}0& 0& 0& 0\\ 0& 0& 0& 0\\ 0& 0& 0& 0\\ 0& K_f\frac{{\mathrm{\ω}}_{n2}^2}{k_2}(a_2-a_1)& 0& 0\end{array}\right]\left[\begin{array}{c}{z}_1(t-\mathrm{\τ})\\ z_2(t-\mathrm{\τ})\\ {\stackrel{\·}{z}}_1(t-\mathrm{\τ})\\ {\stackrel{\·}{z}}_2(t-\mathrm{\τ})\end{array}\right]\end{array}---\left(3\right)$

order toEquation (3) becomes:

$\stackrel{\·}{Z}\left(t\right)=AZ\left(t\right)+B_{1}Z(t-\frac{\mathrm{\τ}}{2})+B_{2}Z(t-\mathrm{\τ})---\left(4\right)$

wherein,representing the velocity of the state, A representing the coefficient matrix of the state displacement term, B₁Representing a matrix of time-lapse state shift terms, B₂Representing a time-lag state shift term coefficient matrix.

Preferably, the step 3 specifically comprises:

let t equal ξ. tau, then

The normalized state space equation becomes the standard form, as shown in equation (5):

$\stackrel{\·}{Z}\left(\mathrm{\ξ}\right)=\mathrm{\τ}AZ\left(\mathrm{\ξ}\right)+{\mathrm{\τB}}_{1}Z(\mathrm{\ξ}-0.5)+{\mathrm{\τB}}_{2}Z(\mathrm{\ξ}-1)---\left(5\right).$

preferably, the step 4 specifically includes:

taking n +1 second class Chebyshev discrete points ξ in the interval_iAs shown in equation (6):

${\mathrm{\ξ}}_i=\frac{1}{2}\[1-cos\left(\frac{i\mathrm{\π}}{n}\right)\],i=0,...,n---\left(6\right)$

then at each discrete point equation (5) needs to be satisfied, namely:

$\stackrel{\·}{Z}\left({\mathrm{\ξ}}_i\right)=\mathrm{\τ}AZ\left({\mathrm{\ξ}}_i\right)+{\mathrm{\τB}}_{1}Z({\mathrm{\ξ}}_i-0.5)+{\mathrm{\τB}}_{2}Z({\mathrm{\ξ}}_i-1)---\left(7\right)$

then, over the interval [0,1],

$\left[\begin{array}{c}\stackrel{\·}{Z}\left({\mathrm{\ξ}}_0\right)\\ \stackrel{\·}{Z}\left({\mathrm{\ξ}}_1\right)\\ .\\ .\\ .\\ \stackrel{\·}{Z}\left({\mathrm{\ξ}}_n\right)\end{array}\right]=\mathrm{\τ}\·I_{n+1}\⊗A\left[\begin{array}{c}Z\left({\mathrm{\ξ}}_0\right)\\ Z\left({\mathrm{\ξ}}_1\right)\\ .\\ .\\ .\\ Z\left({\mathrm{\ξ}}_n\right)\end{array}\right]+\mathrm{\τ}\·I_{n+1}\⊗B_1\left[\begin{array}{c}Z({\mathrm{\ξ}}_0-0.5)\\ Z({\mathrm{\ξ}}_1-0.5)\\ .\\ .\\ .\\ Z({\mathrm{\ξ}}_n-0.5)\end{array}\right]+\mathrm{\τ}\·I_{n+1}\⊗B_2\left[\begin{array}{c}Z({\mathrm{\ξ}}_0-1)\\ Z({\mathrm{\ξ}}_1-1)\\ .\\ .\\ .\\ Z({\mathrm{\ξ}}_n-1)\end{array}\right]---\left(8\right)$

wherein,representing the Kronecker product, i.e.I_n+1Represents an (n +1) -dimensional square matrix.

Preferably, the step 5 specifically includes:

expressing a velocity term at the left end of the formula (8) by using a displacement term by using a differential quadrature method; for convenience of presentation, noteIs composed off(ξ_i) A scalar expression symbol representing a state displacement term;

first using (n +1) points (ξ)₀,f(ξ₀)),(ξ₁,f(ξ₁)),…,(ξ_n,f(ξ_n) Lagrange interpolation is performed, the result is shown in equation (9):

$f\left(\mathrm{\ξ}\right)=\underset{j=0}{\overset{n}{\Σ}}{l}_j\left(\mathrm{\ξ}\right)f\left({\mathrm{\ξ}}_j\right)---\left(9\right)$

wherein the interpolation basis function is:

$l_i\left(\mathrm{\ξ}\right)=\underset{k\≠i}{\underset{k=0}{\overset{n}{\Π}}}\frac{(\mathrm{\ξ}-{\mathrm{\ξ}}_k)}{({\mathrm{\ξ}}_i-{\mathrm{\ξ}}_k)}---\left(10\right)$

derivation is performed on the lagrange interpolation function, and then each time point of the left-end speed term of the formula (8) is substituted, and the result is shown in the formula (11):

wherein, the expression of H is shown as formula (12):

$H_{ij}=\frac{{\mathrm{dl}}_j\left(\mathrm{\ξ}\right)}{d\mathrm{\ξ}}{|}_{\mathrm{\ξ}={\mathrm{\ξ}}_i}=\left\{\begin{array}{cc}\frac{\underset{k\≠i,j}{\underset{k=0}{\overset{n}{\Π}}}({\mathrm{\ξ}}_i-{\mathrm{\ξ}}_k)}{\underset{k\≠j}{\underset{k=0}{\overset{n}{\Π}}}({\mathrm{\ξ}}_j-{\mathrm{\ξ}}_k)},& j\≠i\\ \underset{k\≠i}{\underset{k=0}{\overset{n}{\Σ}}}\frac{1}{{\mathrm{\ξ}}_j-{\mathrm{\ξ}}_k},& j=i\end{array}\right.---\left(12\right).$

preferably, the step 6 specifically includes:

the time lag shift term Z (ξ -0.5) is [ -1,0 [ -1]And [0,1]]The second Chebyshev point representation on the interval, find i so that ξ_i-0.5 < 0 and ξ_i+1-0.5 > 0(i ═ 0,1, …, n); with [ -1,0 [ ]]And [0,1]]And taking the second class of Chebyshev points on the interval as interpolation points to obtain a Lagrange interpolation function, and substituting time points of time-lag displacement terms to obtain:

wherein,

$T_{ij}^d=\left\{\begin{array}{cc}\frac{\underset{k\&NotEqual;j}{\underset{k=0}{\overset{n}{\&Pi;}}}\&lsqb;({\mathrm{\&xi;}}_i-0.5)-({\mathrm{\&xi;}}_k-1)\&rsqb;}{\underset{k\&NotEqual;j}{\underset{k=0}{\overset{n}{\&Pi;}}}\&lsqb;({\mathrm{\&xi;}}_j-1)-({\mathrm{\&xi;}}_k-1)\&rsqb;},& ({\mathrm{\&xi;}}_i-0.5)\&NotElement;\{{\mathrm{\&xi;}}_k-1|k=0,1,\mathrm{...},n\}\\ 0,& ({\mathrm{\&xi;}}_i-0.5)\&Element;\{{\mathrm{\&xi;}}_k-1|k=0,1,\mathrm{...},j-1,j+1,\mathrm{...},n\}\\ 1,& {\mathrm{\&xi;}}_i-0.5={\mathrm{\&xi;}}_j-1\end{array}\right.---\left(14\right)$

wherein,

$T_{ij}^p=\left\{\begin{array}{cc}\frac{\underset{k\&NotEqual;j}{\underset{k=0}{\overset{n}{\&Pi;}}}\&lsqb;({\mathrm{\&xi;}}_i-0.5)-{\mathrm{\&xi;}}_k\&rsqb;}{\underset{k\&NotEqual;j}{\underset{k=0}{\overset{n}{\&Pi;}}}({\mathrm{\&xi;}}_i-{\mathrm{\&xi;}}_k)},& ({\mathrm{\&xi;}}_i-0.5)\&NotElement;\{{\mathrm{\&xi;}}_k-1|k=0,1,\mathrm{...},n\}\\ 0,& ({\mathrm{\&xi;}}_i-{0.5}_1)\&Element;\left\{{\mathrm{\&xi;}}_k\right|k=0,1,\mathrm{...},j-1,j+1,\mathrm{...},n\}\\ 1,& {\mathrm{\&xi;}}_i-0.5={\mathrm{\&xi;}}_j\end{array}\right.---\left(16\right)$

the dispersion of the time lag term Z (ξ -1) over the interval [ -1,0], as shown in equation (17):

can obtain T ═ I_n+1，

Wherein T represents

Preferably, the step 7 specifically includes:

constructing a Floquet transfer matrix:

matrix H, T^d,T^pAnd after T disappears the first line is respectively marked as

Wherein, T^dTo representT^pTo represent

By combining the above formulas, the state transition relationship between two adjacent intervals [ -1,0] and [0,1] can be obtained as follows:

$P\left[\begin{array}{c}Z\left({\mathrm{\&xi;}}_0\right)\\ Z\left({\mathrm{\&xi;}}_1\right)\\ .\\ .\\ .\\ Z\left({\mathrm{\&xi;}}_n\right)\end{array}\right]=Q\left[\begin{array}{c}Z({\mathrm{\&xi;}}_0-1)\\ Z({\mathrm{\&xi;}}_1-1)\\ .\\ .\\ .\\ Z({\mathrm{\&xi;}}_n-1)\end{array}\right]---\left(18\right)$

wherein,

the Floquet state transition matrix phi is shown as the formula:

wherein,representing the Penrose-Moor generalized inverse;

according to the Floquet theory, if the moduli of all the characteristic values of phi are less than 1, the system is stable; if the modulus of any characteristic value in phi is larger than 1, the system is unstable.

Preferably, the method further comprises the following steps:

and 8: and drawing a stability map of the system in a time-lag parameter space.

Compared with the prior art, the invention has the following beneficial effects:

compared with the traditional single turning tool turning processing, the method adopts a differential quadrature method to analyze the dynamic characteristics of a multi-tool parallel turning system, obtains optimized cutting parameters, and greatly improves the processing efficiency.

Drawings

Other features, objects and advantages of the invention will become more apparent upon reading of the detailed description of non-limiting embodiments with reference to the following drawings:

FIG. 1 is a schematic view of a multi-tool parallel turning system, with two tools mounted on different tool holders.

FIG. 2 is a view of multi-lathe parallel turning in lathe cutting depth parameter space (a)₁,a₂) Wherein the black parts represent unstable regions.

In FIG. 1, the formulae have the following meanings:

a denotes a first turning tool, B denotes a second turning tool, C_wExpressing damping coefficient, K, of the workpiece_wRepresenting the stiffness coefficient, C, of the workpiece₁Expressing the damping coefficient, K, of the first tool₁Representing the stiffness coefficient, C, of the first turning tool₂Expressing the damping coefficient, K, of the second tool₂Representing the stiffness coefficient of the second tool, a₁Indicating the depth of cut of the first turning tool, a₂Indicating the depth of cut, z, of the second turning tool₁Representing the kinematically responsive displacement of the first turning tool, z₂And the dynamic response displacement of the second turning tool is shown, and Z represents the flexibility direction of the turning tool.

In FIG. 2, a₁Indicating the depth of cut of the first turning tool, a₂Indicating the depth of cut of the second turning tool.

Detailed Description

The present invention will be described in detail with reference to specific examples. The following examples will assist those skilled in the art in further understanding the invention, but are not intended to limit the invention in any way. It should be noted that variations and modifications can be made by persons skilled in the art without departing from the spirit of the invention. All falling within the scope of the present invention.

The invention provides a method for judging the stability of multi-lathe parallel turning based on a differential quadrature method, which comprises the following steps:

(1) and establishing a dynamic equation of multi-lathe tool parallel turning, namely a multi-time-lag second-order differential equation. And determining relevant parameters in the equation according to the vibration test/modal test result.

(2) And carrying out state space transformation on the differential equation to obtain a state space equation.

(3) And (4) carrying out normalization processing on the state space equation, and converting the state space equation into a standard form convenient for solving by a differential integration method.

(4) Based on a differential quadrature method, discrete points on adjacent unit intervals [0,1] or [ -1,0] are used for representing derivative terms and time lag terms in the state space equation, and the multi-time lag equation is equivalently converted into a group of algebraic equation sets only containing cutting parameters.

(5) And judging the stability of the parallel turning process of the multiple turning tools based on the Floquet theorem, and drawing a stability map in the cutting parameter space.

More specifically, the following description will be given of a specific embodiment of the present invention with reference to a specific machining example, the parameters of which are given in experiment 2 in document 2. two turning tools are adopted for parallel turning, and the turning tools are mounted on different tool holders without vibration coupling therebetween, and a schematic diagram thereof is shown in fig. 1. Recording a turning tool with smaller cutting depth as a first turning tool, recording a turning tool with larger cutting depth as a second turning tool, wherein the diameter d of a processed workpiece is 35mm, the rotating speed omega of a main shaft is 1800rpm, and the cutting force coefficient K in the feeding direction_fIs K_f872MPa, modal test parameters are shown in table 1:

TABLE 1 Multi-cutter parallel cutting example Modal test parameters

Mode vibration mode	fn(Hz)	ζ(％)	k(N/m)
				1	2238.9	3.23	4.769*107
2	2372.3	4.51	1.166*108

Wherein f is_nDenotes the natural frequency, ζ denotes the relative damping ratio, and k denotes the stiffness coefficient.

Step (1), establishing a multi-turning tool parallel turning dynamics equation as shown in a formula (1):

$\begin{array}{c}\left[\begin{array}{c}{\stackrel{\&CenterDot;\&CenterDot;}{z}}_1\left(t\right)\\ {\stackrel{\&CenterDot;\&CenterDot;}{z}}_2\left(t\right)\end{array}\right]+\left[\begin{array}{cc}2{\mathrm{\&zeta;}}_1{\mathrm{\&omega;}}_{n1}& 0\\ 0& 2{\mathrm{\&zeta;}}_2{\mathrm{\&omega;}}_{n2}\end{array}\right]\left[\begin{array}{c}{\stackrel{\&CenterDot;}{z}}_1\left(t\right)\\ {\stackrel{\&CenterDot;}{z}}_2\left(t\right)\end{array}\right]+\left[\begin{array}{cc}{\mathrm{\&omega;}}_{n1}^2& 0\\ 0& {\mathrm{\&omega;}}_{n2}^2\end{array}\right]\left[\begin{array}{c}{z}_1\left(t\right)\\ z_2\left(t\right)\end{array}\right]=\\ K_f\left[\begin{array}{c}\frac{{\mathrm{\&omega;}}_{n1}^2}{k_1}{a}_1(-z_1(t)+z_2(t-\frac{\mathrm{\&tau;}}{2}\left)\right)\\ \frac{{\mathrm{\&omega;}}_{n2}^2}{k_2}{a}_1(-z_2(t)+z_1(t-\frac{\mathrm{\&tau;}}{2}\left)\right)+\frac{{\mathrm{\&omega;}}_{n2}^2}{k_2}(a_2-a_1)(-z_2(t)+z_2(t-\mathrm{\&tau;}\left)\right)\end{array}\right]\end{array}---\left(1\right)$

wherein subscript 1 represents a first turning tool, subscript 2 represents a second turning tool, z (t) represents a kinematically responsive displacement, ω_nWhich represents the natural circular frequency of the frequency,K_fdenotes a cutting force coefficient in a feed direction, τ denotes a workpiece rotation period, ζ denotes a relative damping ratio, k denotes a stiffness coefficient, a₁Indicating the depth of cut of the first turning tool, a₂Indicating the depth of cut of the second turning tool,the acceleration is represented by the acceleration of the vehicle,representing speed and t time.

The formula (1) is collated, and a multiple-time-lag second-order differential equation is shown as a formula (2):

$\begin{array}{c}\left[\begin{array}{c}{\stackrel{\&CenterDot;\&CenterDot;}{z}}_1\left(t\right)\\ {\stackrel{\&CenterDot;\&CenterDot;}{z}}_2\left(t\right)\end{array}\right]+\left[\begin{array}{cc}2{\mathrm{\&zeta;}}_1{\mathrm{\&omega;}}_{1n}& 0\\ 0& 2{\mathrm{\&zeta;}}_2{\mathrm{\&omega;}}_{n2}\end{array}\right]\left[\begin{array}{c}{\stackrel{\&CenterDot;}{z}}_1\left(t\right)\\ {\stackrel{\&CenterDot;}{z}}_2\left(t\right)\end{array}\right]+\left[\begin{array}{cc}{\mathrm{\&omega;}}_{n1}^2+K_f\frac{{\mathrm{\&omega;}}_{n1}^2}{k_1}{a}_1& 0\\ 0& {\mathrm{\&omega;}}_{n2}^2+K_f\frac{{\mathrm{\&omega;}}_{n2}^2}{k_2}{a}_2\end{array}\right]\left[\begin{array}{c}{z}_1\left(t\right)\\ z_2\left(t\right)\end{array}\right]=\\ \left[\begin{array}{cc}0& K_f\frac{{\mathrm{\&omega;}}_{n1}^2}{k_1}{a}_1\\ K_f\frac{{\mathrm{\&omega;}}_{n2}^2}{k_2}{a}_1& 0\end{array}\right]\left[\begin{array}{c}{z}_1(t-\frac{\mathrm{\&tau;}}{2})\\ z_2(t-\frac{\mathrm{\&tau;}}{2})\end{array}\right]+\left[\begin{array}{cc}0& 0\\ 0& K_f\frac{{\mathrm{\&omega;}}_{n2}^2}{k_2}(a_2-a_1)\end{array}\right]\left[\begin{array}{c}{z}_1(t-\mathrm{\&tau;})\\ z_2(t-\mathrm{\&tau;})\end{array}\right]\end{array}---\left(2\right)$

and (2) carrying out state space transformation on the formula (2) to obtain a state space equation shown as a formula (3):

$\begin{array}{c}\left[\begin{array}{c}{\stackrel{\&CenterDot;}{z}}_1\left(t\right)\\ {\stackrel{\&CenterDot;}{z}}_2\left(t\right)\\ {\stackrel{\&CenterDot;\&CenterDot;}{z}}_1\left(t\right)\\ {\stackrel{\&CenterDot;\&CenterDot;}{z}}_2\left(t\right)\end{array}\right]=\left[\begin{array}{cccc}0& 0& 1& 0\\ 0& 0& 0& 1\\ -{\mathrm{\&omega;}}_{n1}^2-K_f\frac{{\mathrm{\&omega;}}_{n1}^2}{k_1}{a}_1& 0& -2{\mathrm{\&zeta;}}_1{\mathrm{\&omega;}}_{n1}& 0\\ 0& -{\mathrm{\&omega;}}_{n2}^2-K_f\frac{{\mathrm{\&omega;}}_{n2}^2}{k_2}{a}_2& 0& -2{\mathrm{\&zeta;}}_2{\mathrm{\&omega;}}_{n2}\end{array}\right]\left[\begin{array}{c}{z}_1\left(t\right)\\ z_2\left(t\right)\\ {\stackrel{\&CenterDot;}{z}}_1\left(t\right)\\ {\stackrel{\&CenterDot;}{z}}_2\left(t\right)\end{array}\right]+\\ \left[\begin{array}{cccc}0& 0& 0& 0\\ 0& 0& 0& 0\\ 0& K_f\frac{{\mathrm{\&omega;}}_{n1}^2}{k_1}{a}_1& 0& 0\\ K_f\frac{{\mathrm{\&omega;}}_{n2}^2}{k_2}{a}_1& 0& 0& 0\end{array}\right]\left[\begin{array}{c}{z}_1(t-\frac{\mathrm{\&tau;}}{2})\\ z_2(t-\frac{\mathrm{\&tau;}}{2})\\ {\stackrel{\&CenterDot;}{z}}_1(t-\frac{\mathrm{\&tau;}}{2})\\ {\stackrel{\&CenterDot;}{z}}_2(t-\frac{\mathrm{\&tau;}}{2})\end{array}\right]+\left[\begin{array}{cccc}0& 0& 0& 0\\ 0& 0& 0& 0\\ 0& 0& 0& 0\\ 0& K_f\frac{{\mathrm{\&omega;}}_{n2}^2}{k_2}(a_2-a_1)& 0& 0\end{array}\right]\left[\begin{array}{c}{z}_1(t-\mathrm{\&tau;})\\ z_2(t-\mathrm{\&tau;})\\ {\stackrel{\&CenterDot;}{z}}_1(t-\mathrm{\&tau;})\\ {\stackrel{\&CenterDot;}{z}}_2(t-\mathrm{\&tau;})\end{array}\right]\end{array}---\left(3\right)$

order toEquation (3) becomes:

$\stackrel{\&CenterDot;}{Z}\left(t\right)=AZ\left(t\right)+B_{1}Z(t-\frac{\mathrm{\&tau;}}{2})+B_{2}Z(t-\mathrm{\&tau;})---\left(4\right)$

wherein,representing the state velocity term, A representing the state displacement term coefficient matrix, B₁Representing a matrix of time-lapse state shift terms, B₂Representing a time-lag state shift term coefficient matrix.

And (3) normalizing the state space equation if t is not ξ & tauThe normalized state space equation becomes the standard form, as shown in equation (5):

$\stackrel{\&CenterDot;}{Z}\left(\mathrm{\&xi;}\right)=\mathrm{\&tau;}AZ\left(\mathrm{\&xi;}\right)+{\mathrm{\&tau;B}}_{1}Z(\mathrm{\&xi;}-0.5)+{\mathrm{\&tau;B}}_{2}Z(\mathrm{\&xi;}-1)---\left(5\right)$

step (4), periodically dispersing the formula (5), converting the equivalent into a group of algebraic equations, and taking n +1 second class Chebyshev discrete points ξ in the interval_iAs shown in equation (6):

${\mathrm{\&xi;}}_i=\frac{1}{2}\&lsqb;1-cos\left(\frac{i\mathrm{\&pi;}}{n}\right)\&rsqb;,i=0,...,n---\left(6\right)$

then at each discrete point equation (5) needs to be satisfied, namely:

$\stackrel{\&CenterDot;}{Z}\left({\mathrm{\&xi;}}_i\right)=\mathrm{\&tau;}AZ\left({\mathrm{\&xi;}}_i\right)+{\mathrm{\&tau;B}}_{1}Z({\mathrm{\&xi;}}_i-0.5)+{\mathrm{\&tau;B}}_{2}Z({\mathrm{\&xi;}}_i-1)---\left(7\right)$

then, over the interval [0,1],

$\left[\begin{array}{c}\stackrel{\&CenterDot;}{Z}\left({\mathrm{\&xi;}}_0\right)\\ \stackrel{\&CenterDot;}{Z}\left({\mathrm{\&xi;}}_1\right)\\ .\\ .\\ .\\ \stackrel{\&CenterDot;}{Z}\left({\mathrm{\&xi;}}_n\right)\end{array}\right]=\mathrm{\&tau;}\&CenterDot;I_{n+1}\&CircleTimes;A\left[\begin{array}{c}Z\left({\mathrm{\&xi;}}_0\right)\\ Z\left({\mathrm{\&xi;}}_1\right)\\ .\\ .\\ .\\ Z\left({\mathrm{\&xi;}}_n\right)\end{array}\right]+\mathrm{\&tau;}\&CenterDot;I_{n+1}\&CircleTimes;B_1\left[\begin{array}{c}Z({\mathrm{\&xi;}}_0-0.5)\\ Z({\mathrm{\&xi;}}_1-0.5)\\ .\\ .\\ .\\ Z({\mathrm{\&xi;}}_n-0.5)\end{array}\right]+\mathrm{\&tau;}\&CenterDot;I_{n+1}\&CircleTimes;B_2\left[\begin{array}{c}Z({\mathrm{\&xi;}}_0-1)\\ Z({\mathrm{\&xi;}}_1-1)\\ .\\ .\\ .\\ Z({\mathrm{\&xi;}}_n-1)\end{array}\right]---\left(8\right)$

wherein,representing the Kronecker product, i.e.I_n+1Represents an (n +1) -dimensional square matrix.

Step (5), differential integration method

The velocity term at the left end of equation (8) is expressed in terms of displacement by differential integration. For convenience of presentation, noteIs composed off(ξ_i) Representing state displacement terms scalar representation notation first uses (n +1) points (ξ)₀,f(ξ₀)),(ξ₁,f(ξ₁)),…,(ξ_n,f(ξ_n) Lagrange interpolation results as shown in equation (9):

$f\left(\mathrm{\&xi;}\right)=\underset{j=0}{\overset{n}{\&Sigma;}}{l}_j\left(\mathrm{\&xi;}\right)f\left({\mathrm{\&xi;}}_j\right)---\left(9\right)$

wherein the interpolation basis function is:

$l_i\left(\mathrm{\&xi;}\right)=\underset{k\&NotEqual;i}{\underset{k=0}{\overset{n}{\&Pi;}}}\frac{(\mathrm{\&xi;}-{\mathrm{\&xi;}}_k)}{({\mathrm{\&xi;}}_i-{\mathrm{\&xi;}}_k)}---\left(10\right)$

derivation is performed on the lagrange interpolation function, and then each time point of the left-end speed term of the formula (8) is substituted, and the result is shown in the formula (11):

wherein, the expression of H is shown as formula (12):

$H_{ij}=\frac{{\mathrm{dl}}_j\left(\mathrm{\&xi;}\right)}{d\mathrm{\&xi;}}{|}_{\mathrm{\&xi;}={\mathrm{\&xi;}}_i}=\left\{\begin{array}{cc}\frac{\underset{k\&NotEqual;i,j}{\underset{k=0}{\overset{n}{\&Pi;}}}({\mathrm{\&xi;}}_i-{\mathrm{\&xi;}}_k)}{\underset{k\&NotEqual;j}{\underset{k=0}{\overset{n}{\&Pi;}}}({\mathrm{\&xi;}}_j-{\mathrm{\&xi;}}_k)},& j\&NotEqual;i\\ \underset{k\&NotEqual;i}{\underset{k=0}{\overset{n}{\&Sigma;}}}\frac{1}{{\mathrm{\&xi;}}_j-{\mathrm{\&xi;}}_k},& j=i\end{array}\right.---\left(12\right)$

the time lag shift term Z (ξ -0.5) is [ -1,0 [ -1]And [0,1]]The second Chebyshev point representation on the interval can be found so that ξ is found_i-0.5 < 0 and ξ_i+1-0.5 > 0(i ═ 0,1, …, n); with [ -1,0 [ ]]And [0,1]]And taking the second class Chebyshev point on the interval as an interpolation point to obtain a Lagrange interpolation function, and substituting the time point of the time-lag displacement term to obtain:

wherein,

$T_{ij}^d=\left\{\begin{array}{cc}\frac{\underset{k\&NotEqual;j}{\underset{k=0}{\overset{n}{\&Pi;}}}\&lsqb;({\mathrm{\&xi;}}_i-0.5)-({\mathrm{\&xi;}}_k-1)\&rsqb;}{\underset{k\&NotEqual;j}{\underset{k=0}{\overset{n}{\&Pi;}}}\&lsqb;({\mathrm{\&xi;}}_j-1)-({\mathrm{\&xi;}}_k-1)\&rsqb;},& ({\mathrm{\&xi;}}_i-0.5)\&NotElement;\{{\mathrm{\&xi;}}_k-1|k=0,1,\mathrm{...},n\}\\ 0,& ({\mathrm{\&xi;}}_i-0.5)\&Element;\{{\mathrm{\&xi;}}_k-1|k=0,1,\mathrm{...},j-1,j+1,\mathrm{...},n\}\\ 1,& {\mathrm{\&xi;}}_i-0.5={\mathrm{\&xi;}}_j-1\end{array}\right.---\left(14\right)$

wherein,

$T_{ij}^p=\left\{\begin{array}{cc}\frac{\underset{k\&NotEqual;j}{\underset{k=0}{\overset{n}{\&Pi;}}}\&lsqb;({\mathrm{\&xi;}}_i-0.5)-{\mathrm{\&xi;}}_k\&rsqb;}{\underset{k\&NotEqual;j}{\underset{k=0}{\overset{n}{\&Pi;}}}({\mathrm{\&xi;}}_i-{\mathrm{\&xi;}}_k)},& ({\mathrm{\&xi;}}_i-0.5)\&NotElement;\{{\mathrm{\&xi;}}_k-1|k=0,1,\mathrm{...},n\}\\ 0,& ({\mathrm{\&xi;}}_i-{0.5}_1)\&Element;\left\{{\mathrm{\&xi;}}_k\right|k=0,1,\mathrm{...},j-1,j+1,\mathrm{...},n\}\\ 1,& {\mathrm{\&xi;}}_i-0.5={\mathrm{\&xi;}}_j\end{array}\right.---\left(16\right)$

the dispersion of the time lag term Z (xi-1) in the interval [ -1,0] is quite simple and convenient, as shown in the formula (17):

is easy to know, T ═ I_n+1，

Wherein T represents

Step (6), constructing a Floquet transfer matrix

Matrix H, T^d,T^pAnd after T disappears the first line is respectively marked as

Wherein, T^dTo representT^pTo represent

By combining the above formulas, the state transition relationship between two adjacent intervals [ -1,0] and [0,1] can be obtained as follows:

$P\left[\begin{array}{c}Z\left({\mathrm{\&xi;}}_0\right)\\ Z\left({\mathrm{\&xi;}}_1\right)\\ .\\ .\\ .\\ Z\left({\mathrm{\&xi;}}_n\right)\end{array}\right]=Q\left[\begin{array}{c}Z({\mathrm{\&xi;}}_0-1)\\ Z({\mathrm{\&xi;}}_1-1)\\ .\\ .\\ .\\ Z({\mathrm{\&xi;}}_n-1)\end{array}\right]---\left(18\right)$

wherein,

the Floquet state transition matrix phi is shown as the formula:

wherein,representing the Penrose-Moor generalized inverse;

according to the Floquet theory, if the moduli of all the characteristic values of phi are less than 1, the system is stable; if the modulus of any characteristic value in phi is larger than 1, the system is unstable.

In the depth of cut parameter space (a) of the first and second turning tools₁,a₂) The stability map of the multi-lathe parallel turning is shown in fig. 2.

According to the document 1, if other experimental conditions are kept unchanged, the limit cutting depth of the stability machining is 3.4mm by using only the first lathe tool for the conventional turning, and the limit cutting depth of the stability machining is 12.6mm by using only the second lathe tool for the conventional turning. According to the attached figure 2, after the parallel turning of a plurality of lathes, the limit cutting depth of the stable processing is greatly improved, if (a) is adopted₁,a₂) And (8,26) cutting depth parameter combination is adopted to perform non-regeneration flutter stability machining, the cutting depth of the first turning tool is increased from 3.4mm to 8mm, and is increased by 135.3%, and the cutting depth of the second turning tool is increased from 12.6mm to 26mm, and is increased by 106.3%. Therefore, compared with the traditional single turning tool cutting, the multi-turning tool parallel turning greatly improves the processing efficiency, the improvement of the processing efficiency is not only reflected in the simultaneous turning of a plurality of turning tools, but also reflected in the stable area of the traditional single turning tool parallel turning.

The foregoing description of specific embodiments of the present invention has been presented. It is to be understood that the present invention is not limited to the specific embodiments described above, and that various changes and modifications may be made by one skilled in the art within the scope of the appended claims without departing from the spirit of the invention.

Claims (9)

1. A multi-lathe parallel turning stability judging method based on a differential quadrature method is characterized by comprising the following steps:

step 1: establishing a multi-lathe tool parallel turning kinetic equation, and sorting the multi-lathe tool parallel turning kinetic equation to obtain a multi-time-lag second-order differential equation;

step 2: carrying out state space transformation on the multi-time-lag second-order differential equation to obtain a state space equation;

and step 3: carrying out normalization processing on the state space equation to obtain a state space equation in a standard form;

and 4, step 4: carrying out periodic dispersion on a state space equation in a standard form, and equivalently converting the state space equation into a group of algebraic equations serving as a state space equation expression;

and 5: based on a Lagrange interpolation function, taking a second class Chebyshev point on a unit interval [0,1] as a discrete point, and expressing a derivative term in a state space equation expression by using a displacement term by using a differential integration method;

step 6: judging the interval where the time lag term in the state space equation expression is located based on a Lagrange interpolation function, and if the discrete point of the time lag term belongs to the interval [0,1], representing the time lag term by taking the second class Chebyshev point on [0,1] as the discrete point; if the discrete point of the time lag term belongs to the interval [ -1,0], the time lag term is represented by taking the second type Chebyshev point on [ -1,0] as the discrete point;

and 7: and constructing a state transition matrix between adjacent unit intervals of [0,1] and [ -1,0], and judging the stability of the original system according to the Floquet theory.

2. The method for determining the stability of the multi-lathe tool parallel turning based on the differential integration method according to claim 1, wherein in the step 1, the multi-lathe tool parallel turning kinetic equation is shown as formula (1):

$\begin{array}{c}\left[\begin{array}{c}{\stackrel{\&CenterDot;\&CenterDot;}{z}}_1\left(t\right)\\ {\stackrel{\&CenterDot;\&CenterDot;}{z}}_2\left(t\right)\end{array}\right]+\left[\begin{array}{cc}2{\mathrm{\&zeta;}}_1{\mathrm{\&omega;}}_{n1}& 0\\ 0& 2{\mathrm{\&zeta;}}_2{\mathrm{\&omega;}}_{n2}\end{array}\right]\left[\begin{array}{c}{\stackrel{\&CenterDot;}{z}}_1\left(t\right)\\ {\stackrel{\&CenterDot;}{z}}_2\left(t\right)\end{array}\right]+\left[\begin{array}{cc}{\mathrm{\&omega;}}_{n1}^2& 0\\ 0& {\mathrm{\&omega;}}_{n2}^2\end{array}\right]\left[\begin{array}{c}{z}_1\left(t\right)\\ z_2\left(t\right)\end{array}\right]=\\ K_f\left[\begin{array}{c}\frac{{\mathrm{\&omega;}}_{n1}^2}{k_1}{a}_1(-z_1(t)+z_2(t-\frac{\mathrm{\&tau;}}{2}\left)\right)\\ \frac{{\mathrm{\&omega;}}_{n2}^2}{k_2}{a}_1(-z_2(t)+z_1(t-\frac{\mathrm{\&tau;}}{2}\left)\right)+\frac{{\mathrm{\&omega;}}_{n2}^2}{k_2}(a_2-a_1)(-z_2(t)+z_2(t-\mathrm{\&tau;}\left)\right)\end{array}\right]\end{array}---\left(1\right)$

wherein subscript 1 represents a first turning tool, subscript 2 represents a second turning tool, z (t) represents a kinematically responsive displacement, ω_nRepresenting natural circular frequency, K_fDenotes a cutting force coefficient in a feed direction, τ denotes a workpiece rotation period, ζ denotes a relative damping ratio, k denotes a stiffness coefficient, a₁Indicating the depth of cut of the first turning tool, a₂Indicating the depth of cut of the second turning tool,it is indicative of the dynamic acceleration response,representing the dynamic velocity response, t represents time;

the multi-lag second-order differential equation is shown in formula (2):

$\begin{array}{c}\left[\begin{array}{c}{\stackrel{\&CenterDot;\&CenterDot;}{z}}_1\left(t\right)\\ {\stackrel{\&CenterDot;\&CenterDot;}{z}}_2\left(t\right)\end{array}\right]+\left[\begin{array}{cc}2{\mathrm{\&zeta;}}_1{\mathrm{\&omega;}}_{1n}& 0\\ 0& 2{\mathrm{\&zeta;}}_2{\mathrm{\&omega;}}_{n2}\end{array}\right]\left[\begin{array}{c}{\stackrel{\&CenterDot;}{z}}_1\left(t\right)\\ {\stackrel{\&CenterDot;}{z}}_2\left(t\right)\end{array}\right]+\left[\begin{array}{cc}{\mathrm{\&omega;}}_{n1}^2+K_f\frac{{\mathrm{\&omega;}}_{n1}^2}{k_1}{a}_1& 0\\ 0& {\mathrm{\&omega;}}_{n2}^2+K_f\frac{{\mathrm{\&omega;}}_{n2}^2}{k_2}{a}_2\end{array}\right]\left[\begin{array}{c}{z}_1\left(t\right)\\ z_2\left(t\right)\end{array}\right]=\\ \left[\begin{array}{cc}0& K_f\frac{{\mathrm{\&omega;}}_{n1}^2}{k_1}{a}_1\\ K_f\frac{{\mathrm{\&omega;}}_{n2}^2}{k_2}{a}_1& 0\end{array}\right]\left[\begin{array}{c}{z}_1(t-\frac{\mathrm{\&tau;}}{2})\\ z_2(t-\frac{\mathrm{\&tau;}}{2})\end{array}\right]+\left[\begin{array}{cc}0& 0\\ 0& K_f\frac{{\mathrm{\&omega;}}_{n2}^2}{k_2}(a_2-a_1)\end{array}\right]\left[\begin{array}{c}{z}_1(t-\mathrm{\&tau;})\\ z_2(t-\mathrm{\&tau;})\end{array}\right]\end{array}---\left(2\right).$

3. the differential quadrature method-based multi-lathe parallel turning stability determination method according to claim 2, wherein in step 2, the state space equation is as shown in formula (3):

$\begin{array}{c}\left[\begin{array}{c}{\stackrel{\&CenterDot;}{z}}_1\left(t\right)\\ {\stackrel{\&CenterDot;}{z}}_2\left(t\right)\\ {\stackrel{\&CenterDot;\&CenterDot;}{z}}_1\left(t\right)\\ {\stackrel{\&CenterDot;\&CenterDot;}{z}}_2\left(t\right)\end{array}\right]=\left[\begin{array}{cccc}0& 0& 1& 0\\ 0& 0& 0& 1\\ -{\mathrm{\&omega;}}_{n1}^2-K_f\frac{{\mathrm{\&omega;}}_{n1}^2}{k_1}{a}_1& 0& -2{\mathrm{\&zeta;}}_1{\mathrm{\&omega;}}_{n1}& 0\\ 0& -{\mathrm{\&omega;}}_{n2}^2-K_f\frac{{\mathrm{\&omega;}}_{n2}^2}{k_2}{a}_2& 0& -2{\mathrm{\&zeta;}}_2{\mathrm{\&omega;}}_{n2}\end{array}\right]\left[\begin{array}{c}{z}_1\left(t\right)\\ z_2\left(t\right)\\ {\stackrel{\&CenterDot;}{z}}_1\left(t\right)\\ {\stackrel{\&CenterDot;}{z}}_2\left(t\right)\end{array}\right]+\\ \left[\begin{array}{cccc}0& 0& 0& 0\\ 0& 0& 0& 0\\ 0& K_f\frac{{\mathrm{\&omega;}}_{n1}^2}{k_1}{a}_1& 0& 0\\ K_f\frac{{\mathrm{\&omega;}}_{n2}^2}{k_2}{a}_1& 0& 0& 0\end{array}\right]\left[\begin{array}{c}{z}_1(t-\frac{\mathrm{\&tau;}}{2})\\ z_2(t-\frac{\mathrm{\&tau;}}{2})\\ {\stackrel{\&CenterDot;}{z}}_1(t-\frac{\mathrm{\&tau;}}{2})\\ {\stackrel{\&CenterDot;}{z}}_2(t-\frac{\mathrm{\&tau;}}{2})\end{array}\right]+\left[\begin{array}{cccc}0& 0& 0& 0\\ 0& 0& 0& 0\\ 0& 0& 0& 0\\ 0& K_f\frac{{\mathrm{\&omega;}}_{n2}^2}{k_2}(a_2-a_1)& 0& 0\end{array}\right]\left[\begin{array}{c}{z}_1(t-\mathrm{\&tau;})\\ z_2(t-\mathrm{\&tau;})\\ {\stackrel{\&CenterDot;}{z}}_1(t-\mathrm{\&tau;})\\ {\stackrel{\&CenterDot;}{z}}_2(t-\mathrm{\&tau;})\end{array}\right]\end{array}---\left(3\right)$

order toEquation (3) becomes:

$\stackrel{\&CenterDot;}{Z}\left(t\right)=AZ\left(t\right)+B_{1}Z(t-\frac{\mathrm{\&tau;}}{2})+B_{2}Z(t-\mathrm{\&tau;})---\left(4\right)$

wherein,representing the velocity of the state, A representing the coefficient matrix of the state displacement term, B₁Representing a matrix of time-lapse state shift terms, B₂Representing a time-lag state shift term coefficient matrix.

4. The method for judging the stability of the multi-lathe parallel turning based on the differential quadrature method according to claim 3, wherein the step 3 specifically comprises the following steps:

let t equal ξ. tau, then

The normalized state space equation becomes the standard form, as shown in equation (5):

$\stackrel{\&CenterDot;}{Z}\left(\mathrm{\&xi;}\right)=\mathrm{\&tau;}AZ\left(\mathrm{\&xi;}\right)+{\mathrm{\&tau;B}}_{1}Z(\mathrm{\&xi;}-0.5)+{\mathrm{\&tau;B}}_{2}Z(\mathrm{\&xi;}-1)---\left(5\right).$

5. the method for judging the stability of the multi-lathe parallel turning based on the differential quadrature method according to claim 4, wherein the step 4 specifically comprises:

taking n +1 second class Chebyshev discrete points ξ in the interval_iAs shown in equation (6):

${\mathrm{\&xi;}}_i=\frac{1}{2}\&lsqb;1-cos\left(\frac{i\mathrm{\&pi;}}{n}\right)\&rsqb;,i=0,1,...,n---\left(6\right)$

then at each discrete point equation (5) needs to be satisfied, namely:

$\stackrel{\&CenterDot;}{Z}\left({\mathrm{\&xi;}}_i\right)=\mathrm{\&tau;}AZ\left({\mathrm{\&xi;}}_i\right)+{\mathrm{\&tau;B}}_{1}Z({\mathrm{\&xi;}}_i-0.5)+{\mathrm{\&tau;B}}_{2}Z({\mathrm{\&xi;}}_i-1)---\left(7\right)$

then, over the interval [0,1],

$\left[\begin{array}{c}\stackrel{\&CenterDot;}{Z}\left({\mathrm{\&xi;}}_0\right)\\ \stackrel{\&CenterDot;}{Z}\left({\mathrm{\&xi;}}_1\right)\\ .\\ .\\ .\\ \stackrel{\&CenterDot;}{Z}\left({\mathrm{\&xi;}}_n\right)\end{array}\right]=\mathrm{\&tau;}\&CenterDot;I_{n+1}\&CircleTimes;A\left[\begin{array}{c}Z\left({\mathrm{\&xi;}}_0\right)\\ Z\left({\mathrm{\&xi;}}_1\right)\\ .\\ .\\ .\\ Z\left({\mathrm{\&xi;}}_n\right)\end{array}\right]+\mathrm{\&tau;}\&CenterDot;I_{n+1}\&CircleTimes;B_1\left[\begin{array}{c}Z({\mathrm{\&xi;}}_0-0.5)\\ Z({\mathrm{\&xi;}}_1-0.5)\\ .\\ .\\ .\\ Z({\mathrm{\&xi;}}_n-0.5)\end{array}\right]+\mathrm{\&tau;}\&CenterDot;I_{n+1}\&CircleTimes;B_2\left[\begin{array}{c}Z({\mathrm{\&xi;}}_0-1)\\ Z({\mathrm{\&xi;}}_1-1)\\ .\\ .\\ .\\ Z({\mathrm{\&xi;}}_n-1)\end{array}\right]---\left(8\right)$

wherein,representing the Kronecker product, i.e.I_n+1Representing an n +1 dimensional square matrix.

6. The method for judging the stability of the multi-lathe parallel turning based on the differential quadrature method according to claim 5, wherein the step 5 specifically comprises:

expressing a velocity term at the left end of the formula (8) by using a displacement term by using a differential quadrature method; for convenience of presentation, noteIs composed off(ξ_i) A scalar expression symbol representing a state displacement term;

first using (n +1) points (ξ)₀,f(ξ₀)),(ξ₁,f(ξ₁)),…,(ξ_n,f(ξ_n) Lagrange interpolation is performed, the result is shown in equation (9):

$f\left(\mathrm{\&xi;}\right)=\underset{j=0}{\overset{n}{\&Sigma;}}{l}_j\left(\mathrm{\&xi;}\right)f\left({\mathrm{\&xi;}}_j\right)---\left(9\right)$

wherein the interpolation basis function is:

$l_i\left(\mathrm{\&xi;}\right)=\underset{k\&NotEqual;i}{\underset{k=0}{\overset{n}{\&Pi;}}}\frac{(\mathrm{\&xi;}-{\mathrm{\&xi;}}_k)}{({\mathrm{\&xi;}}_i-{\mathrm{\&xi;}}_k)}---\left(10\right)$

derivation is performed on the lagrange interpolation function, and then each time point of the left-end speed term of the formula (8) is substituted, and the result is shown in the formula (11):

wherein, the expression of H is shown as formula (12):

$H_{ij}=\frac{{\mathrm{dl}}_j\left(\mathrm{\&xi;}\right)}{d\mathrm{\&xi;}}{|}_{\mathrm{\&xi;}={\mathrm{\&xi;}}_i}=\left\{\begin{array}{cc}\frac{\underset{k\&NotEqual;i,j}{\underset{k=0}{\overset{n}{\&Pi;}}}({\mathrm{\&xi;}}_i-{\mathrm{\&xi;}}_k)}{\underset{k\&NotEqual;j}{\underset{k=0}{\overset{n}{\&Pi;}}}({\mathrm{\&xi;}}_j-{\mathrm{\&xi;}}_k)},& j\&NotEqual;i\\ \underset{k\&NotEqual;i}{\underset{k=0}{\overset{n}{\&Sigma;}}}\frac{1}{{\mathrm{\&xi;}}_j-{\mathrm{\&xi;}}_k},& j=i\end{array}\right.---\left(12\right).$

7. the method for judging the stability of the multi-lathe parallel turning based on the differential quadrature method according to claim 6, wherein the step 6 specifically comprises:

the time lag shift term Z (ξ -0.5) is [ -1,0 [ -1]And [0,1]]The second Chebyshev point representation on the interval, find i so that ξ_i-0.5 < 0 and ξ_i+1-0.5 > 0, i ═ 0,1, …, n; respectively with [ -1,0 [ ]]And [0,1]]And taking the second class of Chebyshev points on the interval as interpolation points to obtain a Lagrange interpolation function, and substituting time points of time-lag displacement terms to obtain:

wherein,

$T_{ij}^d=\left\{\begin{array}{cc}\frac{\underset{k\&NotEqual;j}{\underset{k=0}{\overset{n}{\&Pi;}}}\&lsqb;({\mathrm{\&xi;}}_i-0.5)-({\mathrm{\&xi;}}_k-1)\&rsqb;}{\underset{k\&NotEqual;j}{\underset{k=0}{\overset{n}{\&Pi;}}}\&lsqb;({\mathrm{\&xi;}}_j-1)-({\mathrm{\&xi;}}_k-1)\&rsqb;},& ({\mathrm{\&xi;}}_i-0.5)\&NotElement;\{{\mathrm{\&xi;}}_k-1|k=0,1,\mathrm{...},n\}\\ 0,& ({\mathrm{\&xi;}}_i-0.5)\&Element;\{{\mathrm{\&xi;}}_k-1|k=0,1,\mathrm{...},j-1,j+1,\mathrm{...},n\}\\ 1,& {\mathrm{\&xi;}}_i-0.5={\mathrm{\&xi;}}_j-1\end{array}\right.---\left(14\right)$

wherein,

$T_{ij}^p=\left\{\begin{array}{cc}\frac{\underset{k\&NotEqual;j}{\underset{k=0}{\overset{n}{\&Pi;}}}\&lsqb;({\mathrm{\&xi;}}_i-0.5)-{\mathrm{\&xi;}}_k\&rsqb;}{\underset{k\&NotEqual;j}{\underset{k=0}{\overset{n}{\&Pi;}}}({\mathrm{\&xi;}}_i-{\mathrm{\&xi;}}_k)},& ({\mathrm{\&xi;}}_i-0.5)\&NotElement;\{{\mathrm{\&xi;}}_k-1|k=0,1,\mathrm{...},n\}\\ 0,& ({\mathrm{\&xi;}}_i-0.5)\&Element;\left\{{\mathrm{\&xi;}}_k\right|k=0,1,\mathrm{...},j-1,j+1,\mathrm{...},n\}\\ 1,& {\mathrm{\&xi;}}_i-0.5={\mathrm{\&xi;}}_j\end{array}\right.---\left(16\right)$

the dispersion of the time lag term Z (ξ -1) over the interval [ -1,0], as shown in equation (17):

can obtain T ═ I_n+1，

Wherein T represents

8. The method for judging the stability of the multi-lathe parallel turning based on the differential quadrature method according to claim 7, wherein the step 7 specifically comprises:

constructing a Floquet transfer matrix:

matrix H, T^d,T^pAnd after T disappears the first line is respectively marked as

Wherein, T^dTo representT^pTo represent

By combining the above formulas, the state transition relationship between two adjacent intervals [ -1,0] and [0,1] can be obtained as follows:

$P\left[\begin{array}{c}Z\left({\mathrm{\&xi;}}_0\right)\\ Z\left({\mathrm{\&xi;}}_1\right)\\ .\\ .\\ .\\ Z\left({\mathrm{\&xi;}}_n\right)\end{array}\right]=Q\left[\begin{array}{c}Z({\mathrm{\&xi;}}_0-1)\\ Z({\mathrm{\&xi;}}_1-1)\\ .\\ .\\ .\\ Z({\mathrm{\&xi;}}_n-1)\end{array}\right]---\left(18\right)$

wherein,

in the formula (19) and the formula (20), 0_4×(4n+4)Zero matrix, I, with dimension 4 × (4n +4)_4×4An identity matrix of dimension 4 × 4, 0_4×4Zero matrix with dimension 4 × 4, 0_4n×4Zero matrix, I, with dimension 4n × 4_n×nAn identity matrix of dimension n × n, 0_4n×4nZero matrix, 0, with dimension 4n × 4n_4×(4n+4)A zero matrix with dimension 4 × (4n + 4);

the Floquet state transition matrix phi is shown as the formula:

wherein,representing the Penrose-Moor generalized inverse;

according to the Floquet theory, if the moduli of all the characteristic values of phi are less than 1, the system is stable; if the modulus of any characteristic value in phi is larger than 1, the system is unstable.

9. The method for judging the stability of the multi-lathe parallel turning based on the differential quadrature method as claimed in claim 1, further comprising the steps of:

and 8: and drawing a stability map of the system in a time-lag parameter space.