Motorized spindle dynamic thermal error modeling method based on time sequence analysis

Abstract

The invention provides a motorized spindle dynamic thermal error modeling method based on time sequence analysis, and belongs to the technical field of machining. Comprising the following steps: acquiring all temperature rise data and all thermal error data of the motorized spindle under different working conditions; performing iterative clustering on all the temperature rise data by using an elbow rule to obtain a plurality of final clusters; performing grey correlation analysis on all clustering temperature rise data and all thermal error data in each final clustering cluster, sorting all correlation degrees in each obtained correlation degree set from large to small, and constructing a feature data set according to a sorting result; selecting an initial dynamic model, and sequentially performing order determination, parameter identification and verification on the initial dynamic model by utilizing all the feature data and all the thermal error data to obtain a final thermal error dynamic model; and performing thermal error prediction and thermal error compensation on the motorized spindle by using the final thermal error dynamic model. According to the invention, the precision and robustness of thermal error prediction of the motorized spindle can be improved.

Description

Technical Field

[0001] This application relates to the technical field of machining, and particularly relates to a method for dynamically modeling the thermal error of an electric spindle based on time series analysis. Background Art

[0002] "Made in China 2025" includes numerically controlled machine tools in the "strategic areas that must be broken through quickly". Numerically controlled machine tools are not only the technical cornerstone of modern industry but also an important pillar of scientific and technological innovation. Numerically controlled machine tools play a crucial role in a country's industrial fields such as metal processing, aerospace, and steel. With the continuous development and progress of industrial technology, numerically controlled machine tools, as a symbol of modern mechanical manufacturing level, have received increasing attention, and the requirements for their machining accuracy have also been continuously improved. Therefore, it is of great significance to improve the machining accuracy of numerically controlled machine tools.

[0003] Numerically controlled machine tools are commonly used to machine threaded parts, thin-walled parts, gear parts, mirror parts, etc. During the machining process of numerically controlled machine tools, due to the friction of components such as motors and bearings inside the electric spindle, the temperature of the electric spindle rises, and then thermal deformation occurs. This thermal deformation will cause the deviation of the axis position of the electric spindle, thereby generating thermal error, that is, the thermal error of the electric spindle. This is the main source of thermal error in the machining process of numerically controlled machine tools, which will reduce the machining accuracy, reliability, and execution accuracy of numerically controlled machine tools.

[0004] Currently, generally, error suppression technology or error compensation technology is used to reduce the adverse effects of thermal error on numerically controlled machine tools. Error suppression technology refers to reducing the sources of thermal error by improving the manufacturing accuracy and assembly accuracy of various components of numerically controlled machine tools. For example, by adopting symmetric structure design, thermal harmony design, etc. to eliminate or reduce the influence of thermal error on numerically controlled machine tools. The cost required for error suppression technology is relatively high, the flexibility is poor, and it is easily interfered by external environmental factors. Therefore, it is rarely used in the production of numerically controlled machine tools. And error compensation technology has attracted much attention from the industrial community and scholars due to its economic and efficient characteristics. Error compensation technology mainly predicts or directly measures thermal error through a model, and then inputs the thermal error compensation value into the numerically controlled machine tool system, thereby improving the machining accuracy in actual production. The strategy of predicting and compensating with a thermal error model has also become the most cost-effective solution to solve the thermal error problem of numerically controlled machine tools.

[0005] However, the existing thermal error compensation technology generally adopts a static modeling method, and this method will have the following problems:

[0006] First, static modeling uses the current temperature or speed of the CNC machine tool as the input of the thermal error model, considering only the relationship between heat input and output within the current time period. However, in actual machining, CNC machine tools are easily affected by factors such as the environment, coolant, and tool wear, which can make static modeling less accurate in predicting long-term thermal errors.

[0007] Secondly, static modeling overemphasizes the static characteristics of thermal errors and ignores their dynamic characteristics. Therefore, when the actual working conditions of the CNC machine tool are different from the test conditions, the prediction performance of static modeling will be reduced.

[0008] Finally, when the temperature field of a CNC machine tool undergoes significant changes, static modeling makes it difficult to describe the mapping between temperature and thermal error. This is because thermal deformation at a given moment is affected not only by the current temperature but also by previous temperatures. Static modeling relies on information at the current moment or at isolated points in time, making it less stable and robust than dynamic modeling.

[0009] Therefore, it is necessary to propose a solution to improve one or more problems existing in the above-mentioned related technical solutions.

[0010] It should be noted that the information disclosed in the above background technology section is only used to enhance the understanding of the background of this application, and therefore may include information that does not constitute prior art known to ordinary technicians in this field. Summary of the Invention

[0011] The present application provides a method for modeling dynamic thermal errors of an electric spindle based on timing analysis, the method comprising the following steps:

[0012] By performing thermal characteristic testing on the electric spindle of the CNC machine tool, all temperature rise data and all thermal error data of the electric spindle under different working conditions are obtained;

[0013] Preprocessing all the temperature rise data, and iteratively clustering all the preprocessed temperature rise data using the elbow rule to obtain a temperature clustering result, wherein the temperature clustering result includes a plurality of final clusters having a cluster center and a plurality of clustered temperature rise data;

[0014] performing grey correlation analysis on all the clustered temperature rise data and all the thermal error data in each of the final clusters to obtain a plurality of correlation degree sets, sorting all the correlation degrees in each of the correlation degree sets from large to small, and using the clustered temperature rise data corresponding to the maximum correlation degree in each of the correlation degree sets as feature data to form a feature data set;

[0015] Select an initial dynamic model, and successively perform order determination, parameter identification, and verification on the initial dynamic model using all the feature data and all the thermal error data to obtain a final thermal error dynamic model;

[0016] Use the final thermal error dynamic model to predict and compensate for the thermal error of the motorized spindle.

[0017] Further, the step of obtaining all the temperature rise data and all the thermal error data of the motorized spindle under different working conditions by performing thermal characteristic detection on the motorized spindle of the CNC machine tool includes:

[0018] Set a plurality of temperature sensors on the motorized spindle, and each of the temperature sensors corresponds to a temperature monitoring point;

[0019] Perform the thermal characteristic detection on the motorized spindle, and use all the temperature sensors to respectively obtain all the temperature rise data and all the thermal error data of the motorized spindle under different working conditions.

[0020] Further, the step of preprocessing all the temperature rise data includes:

[0021] Perform cleaning processing on all the temperature rise data, remove all the missing temperature rise data and all the abnormal temperature rise data therein to obtain a normal temperature rise data set, and the normal temperature rise data set includes a plurality of normal temperature rise data;

[0022] Perform standardization processing on the normal temperature rise data set to obtain a standardized temperature rise data set, and the standardized temperature rise data set includes a plurality of standardized temperature rise data.

[0023] Further, the expression of the matrix form of the normal temperature rise data set is:

[0024]

[0025] wherein, G represents the normal temperature rise data set, G n [ represents all the normal temperature rise data of the nth temperature sensor in the normal temperature rise data set, G n [ (m) represents the mth normal temperature rise data of the nth temperature sensor in the normal temperature rise data set;

[0026] The expression of the standardized temperature rise data set is:

[0027]

[0028] wherein, represents the average value of all the normal temperature rise data of the nth temperature sensor, m represents the number of normal temperature rise data, j represents the jth normal temperature rise data, G n [(j) represents the j-th normal temperature rise data of the n-th temperature sensor, s n represents the standard deviation of all the normal temperature rise data of the n-th temperature sensor, M n represents all the standardized temperature rise data of the n-th temperature sensor, and N represents the number of temperature sensors.

[0029] Further, the step of iteratively clustering all the preprocessed temperature rise data by using the elbow method to obtain the temperature clustering result includes:

[0030] Randomly select a plurality of the standardized temperature rise data from the standardized temperature rise data set, and respectively use all the randomly selected standardized temperature rise data as the initial clustering centers;

[0031] Respectively calculate the sum of the squared distances between each of the standardized temperature rise data in the standardized temperature rise data set and all the initial clustering centers, and respectively assign each of the standardized temperature rise data to the initial clustering cluster represented by the initial clustering center corresponding to the minimum value of the sum of the squared distances;

[0032] Respectively calculate the mean value of all the standardized temperature rise data in each of the initial clustering clusters, and respectively use the mean value of all the standardized temperature rise data as the new clustering center;

[0033] Iteratively perform the steps of calculating the sum of the squared distances and updating the clustering centers until the preset number of iterations is reached or the updated clustering centers no longer change, to obtain a plurality of the final clustering clusters;

[0034] Wherein, each of the final clustering clusters respectively contains a plurality of the clustering temperature rise data.

[0035] Further, the expression for calculating the sum of the squared distances is:

[0036] d nx =‖M n -P x ‖ 2 (3)

[0037] Wherein, d nx represents the sum of the squared distances between all the standardized temperature rise data of the n-th temperature sensor and the x-th clustering center, M n represents all the standardized temperature rise data of the n-th temperature sensor, P x represents the x-th clustering center;

[0038] The expression for updating the clustering center is:

[0039]

[0040] Among them, P x ′ represents the new cluster center, and |P x | represents the number of cluster centers.

[0041] Furthermore, the expression for calculating the correlation degree is:

[0042]

[0043] Among them, φ(y, h n ) represents the correlation degree between all the thermal error data and the cluster temperature rise data of the nth temperature sensor respectively, y represents all the thermal error data, and h n represents all the cluster temperature rise data of the nth temperature sensor, y k represents the kth thermal error data, and h ni represents the ith cluster temperature rise data of the nth temperature sensor, N represents the number of temperature sensors, and γ represents the correlation coefficient;

[0044] The expression for the correlation coefficient is:

[0045]

[0046] Among them, γ(y k , h ni ) represents the correlation coefficient between the kth thermal error data and the ith cluster temperature rise data of the nth temperature sensor, min n min i |y k - h ni | represents the minimum absolute value of the difference between the kth thermal error data and the ith cluster temperature rise data of the nth temperature sensor, max n max i |y k - h ni | represents the maximum absolute value of the difference between the kth thermal error data and the ith cluster temperature rise data of the nth temperature sensor, η represents the discrimination coefficient, and |y k - h ni | represents the absolute value of the difference between the kth thermal error data and the ith cluster temperature rise data of the nth temperature sensor.

[0047] Furthermore, the steps of selecting the initial dynamic model and sequentially performing order determination, parameter identification, and verification on the initial dynamic model by using all the feature data and all the thermal error data to obtain the final thermal error dynamic model include:

[0048] Taking the feature time series corresponding to all the feature data as the input and the thermal error time series corresponding to all the thermal error data as the output, define the initial dynamic model;

[0049] Order the initial dynamic model using the Akaike information criterion;

[0050] Identify the parameters of the initial dynamic model after ordering using the least squares method;

[0051] Perform cross-validation and performance index verification on the initial dynamic model that has successively undergone the ordering and parameter identification to obtain the final thermal error dynamic model.

[0052] Furthermore, the expression of the initial dynamic model is:

[0053]

[0054] where u[t] represents the input term of the initial dynamic model, v[t] represents the output term of the initial dynamic model, e[t] represents the disturbance term of the initial dynamic model, A(q) represents the influence polynomial of the output term, B(q) represents the first influence polynomial of the input term, F(q) represents the second influence polynomial of the input term, C(q) represents the first influence polynomial of the disturbance term, D(q) represents the second influence polynomial of the disturbance term, t represents the t-th moment, n a represents the order of the influence polynomial of the output term, represents the undetermined coefficient of the n a -th order, n b represents the order of the first influence polynomial of the input term, represents the undetermined coefficient of the n b -th order, n f represents the order of the second influence polynomial of the input term, represents the undetermined coefficient of the n f -th order, n c represents the order of the first influence polynomial of the disturbance term, represents the undetermined coefficient of the n c -th order, n d represents the order of the second influence polynomial of the disturbance term, represents the undetermined coefficient of the n d -th order, and q represents the backward shift operator;

[0055] When C(q) = 1, D(q) = 1, and F(q) = 1, the initial dynamic model is an ARX model;

[0056] When B(q) = 1, C(q) = 1, and F(q) = 1, the initial dynamic model is an AR model;

[0057] When B(q) = 1, D(q) = 1, and F(q) = 1, the initial dynamic model is an ARMA model.

[0058] In an exemplary embodiment of the present application, the expression for order determination is:

[0059] AIC = 2w - 2ln(L) (8)

[0060] Where AIC represents the information complexity of the initial dynamic model, w represents the number of parameters of the initial dynamic model, and L represents the maximum likelihood function value of the initial dynamic model;

[0061] The expression for parameter identification is:

[0062]

[0063] Where, represents the identification estimate value of the parameter matrix, t represents the t-th moment, represents the time series data of the input and output of the initial dynamic model at the t-th moment, E(t) represents the actual thermal error data measured at the t-th moment, T represents the transpose, n is the n-th temperature sensor, and N represents the number of temperature sensors.

[0064] Beneficial effects:

[0065] The present application provides a method for modeling the dynamic thermal error of an electric spindle based on time series analysis, which has at least the following beneficial effects:

[0066] (1) By detecting the thermal characteristics of the electric spindle and using the elbow method and gray correlation analysis method to construct a feature dataset, the optimal thermal error sensitive points are selected, effectively solving the collinearity problem between temperature sensitive points;

[0067] (2) By using all feature data and all thermal error data to perform order determination, parameter identification, and verification on the initial dynamic model in sequence, the final thermal error dynamic model is obtained, solving the pseudo-lag problem caused by thermal deformation of the electric spindle and improving the accuracy and robustness of electric spindle thermal error prediction. Description of the Drawings

[0068] The drawings here are incorporated into the specification and form a part of this specification, showing embodiments consistent with the present application, and are used together with the specification to explain the principles of the present application. Obviously, the drawings in the following description are only some embodiments of the present application. For those of ordinary skill in the art, other drawings can be obtained based on these drawings without creative efforts.

[0069] Figure 1Schematic diagram showing the steps of a dynamic thermal error modeling method for an electric spindle based on time series analysis in an exemplary embodiment of the present application;

[0070] Figure 2 Schematic diagram showing the process speed spectrum of the electric spindle of a longitudinal cutting lathe under different working conditions in an exemplary embodiment of the present application;

[0071] Figure 3 Schematic diagram showing the principle of the elbow method in an exemplary embodiment of the present application;

[0072] Figure 4 Schematic diagram showing the optimal number of clusters of the electric spindle under different working conditions in an exemplary embodiment of the present application;

[0073] Figure 5 Schematic diagram showing the structure of the ARX model in an exemplary embodiment of the present application;

[0074] Figure 6 Curve graph showing the thermal error prediction of the electric spindle under different working conditions in an exemplary embodiment of the present application;

[0075] Figure 7 Curve graph showing the thermal error residuals of the electric spindle under different working conditions in an exemplary embodiment of the present application;

[0076] Figure 8 Schematic diagram showing the thermal error prediction effect of the electric spindle under different working conditions in an exemplary embodiment of the present application. Detailed implementation manners

[0077] Now, example embodiments will be described more fully with reference to the accompanying drawings. However, the example embodiments can be implemented in various forms and should not be construed as limited to the examples set forth herein; rather, these embodiments are provided so that this application will be more complete and comprehensive, and will fully convey the concept of the example embodiments to those skilled in the art. The described features, structures, or characteristics can be combined in any suitable manner in one or more embodiments.

[0078] In addition, the accompanying drawings are only schematic illustrations of the present application and are not necessarily drawn to scale. The same reference numerals in the drawings denote the same or similar parts, and thus repeated descriptions thereof will be omitted. Some of the block diagrams shown in the drawings are functional entities and do not necessarily correspond to physically or logically independent entities. These functional entities can be implemented in software form, or in one or more hardware modules or integrated circuits, or in different networks and / or processor devices and / or microcontroller devices.

[0079] For this reason, the present example embodiment provides a dynamic thermal error modeling method for an electric spindle based on time series analysis, as Figure 1As shown, the method may include the following steps:

[0080] Step S101: By detecting the thermal characteristics of the electric spindle of the CNC machine tool, all temperature rise data and all thermal error data of the electric spindle under different working conditions are obtained.

[0081] Step S102: Preprocess all the temperature rise data, and use the elbow method to perform iterative clustering on the preprocessed all temperature rise data to obtain a temperature clustering result. The temperature clustering result includes several final clustering clusters composed of one clustering center and multiple clustering temperature rise data.

[0082] Step S103: Respectively perform grey relational analysis on all the clustering temperature rise data in each final clustering cluster and all the thermal error data to obtain multiple sets of correlation degrees. Sort all the correlation degrees in each set of correlation degrees from large to small, and respectively use the clustering temperature rise data corresponding to the maximum correlation degree in each set of correlation degrees as characteristic data to form a characteristic data set.

[0083] Step S104: Select an initial dynamic model, and use all the characteristic data and all the thermal error data to perform order determination, parameter identification, and verification on the initial dynamic model in sequence to obtain a final thermal error dynamic model.

[0084] Step S105: Use the final thermal error dynamic model to predict and compensate the thermal error of the electric spindle.

[0085] The embodiment of the present application proposes a dynamic thermal error modeling method for an electric spindle based on time series analysis, which has at least the following beneficial effects:

[0086] (1) In the present application, by detecting the thermal characteristics of the electric spindle and using the elbow method and grey relational analysis method to construct a characteristic data set, the optimal thermal error sensitive points are selected from it, thus effectively solving the collinearity problem between temperature sensitive points;

[0087] (2) In the present application, by using all the characteristic data and all the thermal error data to perform order determination, parameter identification, and verification on the initial dynamic model in sequence to obtain a final thermal error dynamic model, the pseudo-lag problem caused by the thermal deformation of the electric spindle is solved, and the accuracy and robustness of the thermal error prediction of the electric spindle are improved.

[0088] Next, a dynamic thermal error modeling method for an electric spindle based on time series analysis proposed in this exemplary embodiment will be described in more detail.

[0089] In step S101 of this embodiment, by detecting the thermal characteristics of the electric spindle of the CNC machine tool, all temperature rise data and all thermal error data of the electric spindle under different working conditions are obtained. Step S101 of this embodiment may include the following sub-steps:

[0090] Sub-step S1011: multiple temperature sensors are set on the electric spindle, and each temperature sensor is equivalent to a temperature monitoring point.

[0091] Furthermore, in this embodiment, it is preferred to use PT100 temperature sensors and contact displacement sensors to measure all temperature rise data and all thermal error data of each temperature monitoring point of the electric spindle of the longitudinal lathe, and these data are transmitted to the supporting host computer software through the HIOKI8423 data acquisition instrument. Figure 2 As shown in the figure, based on the machining characteristics of the sliding head lathe, multiple speed profiles for thermal error detection have been designed. These profiles fully consider key factors such as the lathe tool cutting frequency and workpiece rotation rate, and are constructed in conjunction with the actual delivery process, providing a scientific and comprehensive basis for setting working conditions for thermal error detection.

[0092] Sub-step S1012: performing thermal characteristic detection on the electric spindle, and using all temperature sensors to respectively obtain all temperature rise data and all thermal error data of the electric spindle under different working conditions.

[0093] In this embodiment, by performing thermal characteristic detection on the electric spindle, all temperature rise data and all thermal error data of the electric spindle under different working conditions can be obtained, thereby providing a data basis for subsequent optimization of the electric spindle monitoring points and establishment of the final thermal error dynamic model.

[0094] In step S102 of this embodiment, all temperature rise data are preprocessed and iteratively clustered using the elbow rule to obtain temperature clustering results. The temperature clustering results include a number of final clusters each having a cluster center and multiple clustered temperature rise data. Step S102 of this embodiment may include the following sub-steps:

[0095] Sub-step S1021: cleaning all temperature rise data, removing all missing temperature rise data and all abnormal temperature rise data, and obtaining a normal temperature rise data set, which includes multiple normal temperature rise data.

[0096] Furthermore, for the convenience of representation, the normal temperature rise data set is expressed in matrix form, and its expression is:

[0097]

[0098] Among them, G represents the normal temperature rise data set, G n represents all normal temperature rise data of the nth temperature sensor in the normal temperature rise data set, G n (m) represents the mth normal temperature rise data of the nth temperature sensor in the normal temperature rise data set.

[0099] Sub-step S1022: Since the values of the original normal temperature rise data fluctuate greatly, in order to eliminate the differences in various data types and reduce the calculation error of cluster analysis, it is necessary to standardize the normal temperature rise data set to obtain a standardized temperature rise data set. The standardized temperature rise data set contains multiple standardized temperature rise data.

[0100] Furthermore, the expression of the standardized temperature rise data set is:

[0101]

[0102] Wherein, represents the average value of all normal temperature rise data of the nth temperature sensor, m represents the number of normal temperature rise data, j represents the jth normal temperature rise data, G n (j) represents the jth normal temperature rise data of the nth temperature sensor, s n represents the standard deviation of all normal temperature rise data of the nth temperature sensor, M n represents all the standardized temperature rise data of the nth temperature sensor, and N represents the number of temperature sensors.

[0103] Sub-step S1023: Randomly select multiple standardized temperature rise data from the standardized temperature rise data set, and use all the randomly selected standardized temperature rise data as the initial cluster centers respectively.

[0104] Sub-step S1024: Calculate the sum of the squares of the distances between each standardized temperature rise data in the standardized temperature rise data set and all the initial cluster centers respectively, and assign each standardized temperature rise data to the initial cluster with the initial cluster center corresponding to the minimum sum of the squares of the distances as the representative.

[0105] Furthermore, the expression of the sum of the squares of the distances is:

[0106] d nx =‖M n -P x ‖ 2 (3)

[0107] Wherein, d nx represents the sum of the squares of the distances between all the standardized temperature rise data of the nth temperature sensor and the xth cluster center respectively, M n represents all the standardized temperature rise data of the nth temperature sensor, and P x represents the xth cluster center.

[0108] Sub-step S1025: Calculate the mean value of all the standardized temperature rise data in each initial cluster respectively, and use the mean value of all the standardized temperature rise data as the new cluster center respectively.

[0109] Sub-step S1026: Iteratively perform the steps of calculating the sum of squared distances and updating the cluster centers until the preset number of iterations is reached or the updated cluster centers no longer change, obtaining multiple final cluster clusters;

[0110] Among them, each final cluster cluster contains multiple cluster temperature rise data respectively.

[0111] Furthermore, the expression for updating the cluster center is:

[0112]

[0113] Among them, P x ′ represents the new cluster center, and |P x | represents the number of cluster centers.

[0114] In this embodiment, the Elbow Method (EM) is a simple and intuitive method that determines the number of classifications of a data set by observing the inflection point of the objective function curve. Its principle is as Figure 3 shown. The calculation principle of the Elbow Method is based on the objective function J k , and this objective function measures the sum of squared distances from the data points in each data category to the corresponding cluster center. If the data points in a data category are denser, the objective function of the data category is smaller; on the contrary, if the data points in a data category are more dispersed, the objective function of the data category is larger. Its calculation formula is expressed as follows:

[0115]

[0116] Among them, J k represents the objective function, k represents the kth classification number, and x represents the set of classification numbers.

[0117] The number of final cluster clusters of the temperature monitoring points of the motorized spindle can be determined by the Elbow Method.

[0118] In step S103 of this embodiment, perform gray relational analysis on all the cluster temperature rise data in each final cluster cluster and all the thermal error data respectively to obtain multiple sets of correlation degrees. Sort all the correlation degrees in each set of correlation degrees from large to small, and use the cluster temperature rise data corresponding to the maximum correlation degree in each set of correlation degrees as feature data to form a feature data set.

[0119] Furthermore, the expression for calculating the correlation degree is:

[0120]

[0121] Among them, φ(y,h n) represents the correlation degree of all thermal error data with all clustering temperature rise data of the nth temperature sensor, y represents all thermal error data, h n represents all clustering temperature rise data of the nth temperature sensor, y k represents the kth thermal error data, h ni represents the ith clustering temperature rise data of the nth temperature sensor, N represents the number of temperature sensors, and γ represents the correlation coefficient.

[0122] Furthermore, the expression of the correlation coefficient is:

[0123]

[0124] Among them, γ(y k ,h ni ) represents the correlation coefficient between the kth thermal error data and the ith clustering temperature rise data of the nth temperature sensor, min n min i |y k -h ni | represents the minimum absolute value of the difference between the kth thermal error data and the ith clustering temperature rise data of the nth temperature sensor, max n max i |y k -h ni | represents the maximum absolute value of the difference between the kth thermal error data and the ith clustering temperature rise data of the nth temperature sensor, η represents the resolution coefficient, and |y k -h ni | represents the absolute value of the difference between the kth thermal error data and the ith clustering temperature rise data of the nth temperature sensor.

[0125] In step S104 of this embodiment, an initial dynamic model is selected, and the initial dynamic model is successively order-determined, parameter-identified, and verified using all feature data and all thermal error data to obtain a final thermal error dynamic model. Step S104 of this embodiment may include the following sub-steps:

[0126] Sub-step S1041: Define an initial dynamic model with the feature time series corresponding to all feature data as the input and the thermal error time series corresponding to all thermal error data as the output.

[0127] Furthermore, the expression of the initial dynamic model is:

[0128]

[0129] Among them, u[t] represents the input term of the initial dynamic model, v[t] represents the output term of the initial dynamic model, e[t] represents the disturbance term of the initial dynamic model, A(q) represents the influence polynomial of the output term, B(q) represents the first influence polynomial of the input term, F(q) represents the second influence polynomial of the input term, C(q) represents the first influence polynomial of the disturbance term, D(q) represents the second influence polynomial of the disturbance term, t represents the t-th moment, n a represents the order of the influence polynomial of the output term, represents the n a order undetermined coefficient, n b represents the order of the first influence polynomial of the input term, represents the n b order undetermined coefficient, n f represents the order of the second influence polynomial of the input term, represents the n f order undetermined coefficient, n c represents the order of the first influence polynomial of the disturbance term, represents the n c order undetermined coefficient, n d represents the order of the second influence polynomial of the disturbance term, represents the n d order undetermined coefficient, q represents the backward shift operator.

[0130] As Figure 5 shown, when C(q) = 1, and D(q) = 1, and F(q) = 1, the initial dynamic model is an ARX model (Auto-Regressive with Extra Inputs Model);

[0131] When B(q) = 1, and C(q) = 1, and F(q) = 1, the initial dynamic model is an AR model (Auto-Regressive Model);

[0132] When B(q) = 1, and D(q) = 1, and F(q) = 1, the initial dynamic model is an ARMA model (Auto-Regressive Model).

[0133] Furthermore, the existing FIR model maps past inputs to the current output, but the FIR model ignores the dynamic information of the output in the lag system, while the ARX model maps both past input data and output data to the current output, and its specific form is as follows:

[0134]

[0135] In the formula, n a represents the order of the output term of the ARX model, and n b represents the order of the input term of the ARX model, represents the undetermined coefficient of the nth a order in the ARX model, represents the undetermined coefficient of the nth b order in the ARX model, E[t - n a represents the nth a undetermined coefficient corresponding to the E[t - n a th output data, u[t - n b represents the u[t - n b th input data corresponding to the nth b undetermined coefficient, and e[t] represents the disturbance term.

[0136] To make the above formula more concise, the backward shift operator q is introduced, and we get:

[0137]

[0138] In the formula, q represents the backward shift operator, u[t], u[t + 1], u[t - 1] are input data, and the so-called input data refers to the obtained temperature rise data and thermal error data.

[0139] Therefore, the difference equation form can be re-expressed as:

[0140]

[0141] The model described by this formula is called a source autoregressive model, which is a time series model and has become the preferred choice for predicting dynamic linear models due to its advantages of simple calculation and small computational load. Here, AR refers to the autoregressive part A(q)E[t], that is, the influence of the previous value of the output on the current value, and X refers to the additional input B(q)u[t], which is also called the exogenous variable. As Figure 6 shown in the structure of the ARX model, where u[t] represents the input term, E[t] represents the output term, and e[t] represents the disturbance term, that is, noise. Both A(q) and B(q) are polynomials.

[0142] For most system identification problems, the model structure of the system to be identified is usually unknown. Therefore, only the structural form of the mathematical model can be inferred based on the prior knowledge of the system to be identified first, and then its coefficients can be estimated. In addition, since the actually collected data is easily affected by various measurement noises, it is not sufficient to model only based on the minimum value of the fitting error of the measured data. It is also necessary to ensure that the model can effectively predict the test data under similar conditions and reduce the order of the model as much as possible while maintaining similar performance.

[0143] Sub-step S1042: Determine the order of the initial dynamic model using the Akaike Information Criterion.

[0144] The Akaike Information Criterion (AIC) is a method proposed by Japanese statistician Hirotugu Akaike for selecting the model order. This method is proposed based on the information-theoretic interpretation of the maximum likelihood method and is essentially to determine the order of the model using the basic information content.

[0145] Furthermore, the expression for determining the order is:

[0146] AIC = 2w - 2ln(L) (8)

[0147] Where, AIC represents the information complexity of the initial dynamic model, w represents the number of parameters of the initial dynamic model, and L represents the maximum likelihood function value of the initial dynamic model.

[0148] Sub-step S1043: Identify the parameters of the initial dynamic model with the determined order using the least squares method.

[0149] Furthermore, as can be seen from formula (8), by collecting the temperature rise data and thermal error data of the motorized spindle, the least squares method is used to identify the two undetermined coefficients of the ARX model. For convenience, the undetermined parameter matrix θ of the ARX model and the information matrix of the collected input data and output data Are denoted in the following vector form:

[0150]

[0151] Furthermore, formula (7) is expressed as the product form of the parameter matrix θ and the information matrix To obtain:

[0152]

[0153] Furthermore, using the least squares method to perform parameter identification on it, the identification estimate value of the parameter matrix θ can be obtained The specific expression for parameter identification is:

[0154]

[0155] Among them, represents the identification and estimation value of the parameter matrix, t represents the t-th moment, represents the time series data of the input and output of the initial dynamic model at the t-th moment, E(t) represents the actual thermal error data measured at the t-th moment, T represents the transpose, n is the n-th temperature sensor, and N represents the number of temperature sensors.

[0156] Sub-step S1044: Perform cross-validation and performance index verification on the initial dynamic model after order determination and parameter identification in sequence to obtain the final thermal error dynamic model. Cross-validation and performance index verification belong to conventional technologies in this field and will not be elaborated here.

[0157] In order to verify the excellent effect of a dynamic thermal error modeling method for an electric spindle based on time series analysis proposed in this application, the following simulation experiments were conducted.

[0158] In this simulation experiment, since the electric spindle is greatly affected by the ambient temperature, the temperature monitoring point (T4) at the machine tool bed is selected as one of the thermal temperature sensitive points. Combining the installation positions of the temperature sensors, mean clustering based on the elbow method is used to classify the other temperature monitoring points under different working conditions (working condition one, working condition two, and working condition three), and the number of clusters is as Figure 4 shown, Figure 4 In it, 4a shows the number of final clustering clusters obtained by the electric spindle under working condition one; 4b shows the number of final clustering clusters obtained by the electric spindle under working condition two; 4c shows the number of final clustering clusters obtained by the electric spindle under working condition three.

[0159] First, the temperature clustering results of the electric spindle under different working conditions obtained are as shown in Table 1 below;

[0160] Table 1 Temperature clustering results of the electric spindle under different working conditions

[0161]

[0162] Secondly, the correlation degrees between the temperature rise data and the thermal error data are as shown in Table 2 below:

[0163] Table 2 Correlation degrees between the temperature rise data and the thermal error data

[0164]

[0165]

[0166] Then, the thermal error sensitive points are screened, as shown in Table 3 below,

[0167] Table 3 Screening results of thermal error sensitive points

[0168]

[0169] As can be seen from Table 2 and Table 3, the screening results of the motorized spindle are the same under Working Condition 1 and Working Condition 3. There are differences in the second type of screening results between Working Condition 2 and the other two working conditions. However, both Temperature Sensors T7 and T8 are used to monitor the temperature rise data of the stator and rotor, and the correlation degrees of the two temperature monitoring points are not very different. Therefore, considering comprehensively the measurement positions of the temperature sensors and the magnitudes of their correlation degrees, the temperature rise data of T1, T4, and T7 are selected as the optimal inputs for the final thermal error dynamic model. Among them, T1 is the temperature rise data at the coupling in the X direction, T4 is the temperature rise data at the machine tool bed, and T7 is the temperature rise data of the stator and rotor.

[0170] In this simulation experiment, the Finite Impulse Response (FIR) model can use the current and previous temperature rise data to calculate thermal deformation. The Multivariable Linear Regression (MLR) model is widely used in the thermal deformation prediction of machine tools due to its simple structure and convenient identification.

[0171] In order to verify the stability and accuracy of the final thermal error dynamic model constructed in this application, the temperature rise data of the motorized spindle collected under Working Condition 1 are used as the training set, and the temperature rise data of the motorized spindle collected under Working Condition 2 and Working Condition 3 are used as the test set. A final thermal error dynamic model based on identification theory is established, and the proposed dynamic model is compared with the FIR model and the MLR model.

[0172] As Figure 6 shown, Figure 6 a is the predicted value of the thermal error under Working Condition 2, Figure 6 b is the predicted value of the thermal error under Working Condition 3; it can be seen from Figure 6 this that the prediction effect of the ARX model scheme is significantly better than that of the static FIR model and MLR model.

[0173] As Figure 7 shown, Figure 7 a is the thermal error residual under Working Condition 2, Figure 7 b is the thermal error residual under Working Condition 3. From the thermal error prediction results of the motorized spindle, it can be seen that the thermal error residual of the ARX thermal error model is within 3μm, while the thermal error residuals of the FIR thermal error model and the MLR thermal error model are within 5μm and 7.2μm respectively.

[0174] Furthermore, the Mean Absolute Error (MAE), Root Mean Square Error (RMSE), and coefficient of determination R2 To measure the accuracy and stability of the predictions of these three models.

[0175] As Figure 8 shown, Figure 8 it shows the cases of the motorized spindle under different operating conditions. Figure 8 a is the case under operating condition 2, Figure 8 b is the case under operating condition 3. The ARX thermal error model, FIR thermal error model, and MLR thermal error model respectively calculate the MAE, REMSE, and R of the thermal error prediction residuals of the motorized spindle. 2 . By comparing with the FIR thermal error model and MLR thermal error model, it is found that for the ARX thermal error model, the RMSE residuals under operating condition 2 are respectively reduced by 1.67 μm and 3.45 μm; the MAE is also respectively reduced by 1.56 μm and 3.22 μm; and R 2 is respectively increased by 20% and 62%. For the ARX thermal error model under operating condition 3, the RMSE residuals are respectively reduced by 1.8 μm and 2.95 μm, and the MAE is also respectively reduced by 1.59 μm and 2.61 μm; and R 2 is respectively increased by 16% and 34%. Thus, it can be seen that the prediction results of the ARX thermal error model are better than those of the MLR thermal error model and FIR thermal error model, and it has better prediction accuracy and robustness.

[0176] In addition, the terms "first" and "second" are only used for descriptive purposes and cannot be construed as indicating or implying relative importance or implicitly specifying the quantity of the indicated technical features. Thus, the features defined with "first" and "second" may explicitly or implicitly include one or more of such features. In the description of the embodiments of the present application, "a plurality" means two or more, unless otherwise specifically defined.

[0177] In the description of this specification, the description with reference to terms such as "one embodiment", "some embodiments", "example", "specific example", or "some examples" means that the specific features, structures, materials, or characteristics described in connection with the embodiment or example are included in at least one embodiment or example of the present application. In this specification, the schematic descriptions of the above terms do not necessarily refer to the same embodiment or example. Moreover, the specific features, structures, materials, or characteristics described can be combined in a suitable manner in any one or more embodiments or examples. In addition, those skilled in the art can combine and combine the different embodiments or examples described in this specification.

[0178] As described above, it is only the specific implementation manner of the present application, but the protection scope of the present application is not limited thereto. Any person skilled in the art within the technical scope disclosed by the present application can easily think of various equivalent modifications or substitutions, and these modifications or substitutions should all be covered within the protection scope of the present application.

[0179] Those skilled in the art will readily conceive of other embodiments of the present application after considering the specification and practicing the invention disclosed herein. The present application is intended to cover any variations, uses, or adaptations of the present application, which follow the general principles of the present application and include the common general knowledge or conventional technical means in the technical field not disclosed by the present application.

Claims

1. A dynamic thermal error modeling method for an electric spindle based on time series analysis, characterized in that The method includes the following steps: By performing thermal characteristic detection on the electric spindle of a numerically controlled machine tool, all temperature rise data and all thermal error data of the electric spindle under different working conditions are obtained; Preprocess all the temperature rise data, and perform iterative clustering on all the preprocessed temperature rise data using the elbow method to obtain a temperature clustering result, where the temperature clustering result includes several final clustering clusters each composed of one clustering center and multiple clustered temperature rise data; Respectively perform grey relational analysis on all the clustered temperature rise data in each final clustering cluster and all the thermal error data to obtain multiple sets of correlation degrees. Sort all the correlation degrees in each set of correlation degrees from largest to smallest, and respectively use the clustered temperature rise data corresponding to the largest correlation degree in each set of correlation degrees as characteristic data to form a characteristic data set; Select an initial dynamic model, and perform order determination, parameter identification, and verification on the initial dynamic model in sequence using all the characteristic data and all the thermal error data to obtain a final thermal error dynamic model; Use the final thermal error dynamic model to perform thermal error prediction and thermal error compensation on the electric spindle.

2. The method for modeling the dynamic thermal error of an electric spindle based on timing analysis according to claim 1, characterized in that, The step of obtaining all the temperature rise data and all the thermal error data of the electric spindle under different working conditions by performing thermal characteristic detection on the electric spindle of a numerically controlled machine tool includes: Set multiple temperature sensors on the electric spindle, and each temperature sensor respectively corresponds to a temperature monitoring point; Perform the thermal characteristic detection on the electric spindle, and use all the temperature sensors to respectively obtain all the temperature rise data and all the thermal error data of the electric spindle under different working conditions.

3. The method for dynamically modeling the thermal error of an electric spindle based on timing analysis according to claim 1, wherein, The step of preprocessing all the temperature rise data includes: Perform cleaning processing on all the temperature rise data to remove all the missing temperature rise data and all the abnormal temperature rise data therein, and obtain a normal temperature rise data set, where the normal temperature rise data set contains multiple normal temperature rise data; Perform standardization processing on the normal temperature rise data set to obtain a standardized temperature rise data set, where the standardized temperature rise data set contains multiple standardized temperature rise data.

4. The method for modeling the dynamic thermal error of an electric spindle based on timing analysis according to claim 3, wherein The expression of the normal temperature rise data set in matrix form is: Among them, G represents the normal temperature rise dataset, and G n represents all the normal temperature rise data of the nth temperature sensor in the normal temperature rise dataset, and G n (m) represents the mth normal temperature rise data of the nth temperature sensor in the normal temperature rise dataset; The expression of the standardized temperature rise data set is: Among them, represents the average value of all normal temperature rise data of the nth temperature sensor, m represents the number of normal temperature rise data, j represents the jth normal temperature rise data, and G n (j) represents the jth normal temperature rise data of the nth temperature sensor, and s n represents the standard deviation of all normal temperature rise data of the nth temperature sensor, and M n represents all the standardized temperature rise data of the nth temperature sensor, and N represents the number of temperature sensors.

5. The method for modeling the dynamic thermal error of an electric spindle based on time series analysis according to claim 3, wherein, The step of performing iterative clustering on all the preprocessed temperature rise data using the elbow method to obtain a temperature clustering result includes: Randomly select multiple standardized temperature rise data from the standardized temperature rise data set, and respectively use all the randomly selected standardized temperature rise data as initial clustering centers; Respectively calculate the sum of the squares of the distances between each standardized temperature rise data in the standardized temperature rise data set and all the initial clustering centers, and respectively assign each standardized temperature rise data to the initial clustering cluster represented by the initial clustering center corresponding to the minimum value of the sum of the squares of the distances; Respectively calculate the mean value of all the standardized temperature rise data in each initial clustering cluster, and respectively use the mean value of all the standardized temperature rise data as a new clustering center; Iteratively perform the steps of calculating the sum of squared distances and updating the cluster centers until the preset number of iterations is reached or the updated cluster centers no longer change, and obtain multiple final cluster clusters; Among them, each final cluster cluster contains multiple pieces of the cluster temperature rise data respectively.

6. The method for dynamically modeling the thermal error of an electric spindle based on timing analysis according to claim 5, wherein The expression for calculating the sum of squared distances is: d nx = ‖M n - P x ‖ 2 (3) Among them, d nx represents the sum of the squares of the distances between all the standardized temperature rise data of the nth temperature sensor and the xth cluster center, M n represents all the standardized temperature rise data of the nth temperature sensor, P x represents the xth cluster center; The expression for updating the cluster centers is: Among them, P' x represents the new cluster center, and |P x | represents the number of cluster centers.

7. The dynamic thermal error modeling method of the motorized spindle based on timing analysis according to claim 1, wherein The expression for calculating the correlation degree is: Among them, φ(y, h n ) represents the correlation degree between all the thermal error data and all the clustering temperature rise data of the nth temperature sensor. y represents all the thermal error data, and h n represents all the clustering temperature rise data of the nth temperature sensor. y k represents the kth thermal error data, and h ni represents the ith clustering temperature rise data of the nth temperature sensor. N represents the number of temperature sensors, and γ represents the correlation coefficient; The expression for the correlation coefficient is: where γ(y k ,h ni ) represents the correlation coefficient between the k-th thermal error data and the i-th clustering temperature rise data of the n-th temperature sensor, min n min i |y k -h ni | represents the minimum absolute value of the difference between the k-th thermal error data and the i-th clustering temperature rise data of the n-th temperature sensor, max n max i |y k -h ni | represents the maximum absolute value of the difference between the k-th thermal error data and the i-th clustering temperature rise data of the n-th temperature sensor, η represents the discrimination coefficient, and |y k -h ni | represents the absolute value of the difference between the k-th thermal error data and the i-th clustering temperature rise data of the n-th temperature sensor.

8. The method for dynamically modeling the thermal error of an electric spindle based on timing analysis according to claim 1, characterized in that The steps of selecting the initial dynamic model and sequentially performing model order determination, parameter identification, and verification on the initial dynamic model using all the feature data and all the thermal error data to obtain the final thermal error dynamic model include: Define the initial dynamic model by taking the feature time series corresponding to all the feature data as the input and the thermal error time series corresponding to all the thermal error data as the output. Perform the model order determination on the initial dynamic model using the Akaike information criterion. Perform the parameter identification on the initial dynamic model after the model order determination using the least squares method. Perform cross-validation and performance index verification on the initial dynamic model that has sequentially undergone the model order determination and the parameter identification to obtain the final thermal error dynamic model.

9. The method for modeling the dynamic thermal error of an electric spindle based on timing analysis according to claim 8, wherein The expression for the initial dynamic model is: where \(u[t]\) represents the input term of the initial dynamic model, \(v[t]\) represents the output term of the initial dynamic model, \(e[t]\) represents the disturbance term of the initial dynamic model, \(A(q)\) represents the influence polynomial of the output term, \(B(q)\) represents the first influence polynomial of the input term, \(F(q)\) represents the second influence polynomial of the input term, \(C(q)\) represents the first influence polynomial of the disturbance term, \(D(q)\) represents the second influence polynomial of the disturbance term, \(t\) represents the \(t\)-th moment, \(n\) a represents the order of the influence polynomial of the output term, represents the \(n\) a order undetermined coefficient, \(n\) b represents the order of the first influence polynomial of the input term, represents the \(n\) b order undetermined coefficient, \(n\) f represents the order of the second influence polynomial of the input term, represents the \(n\) f order undetermined coefficient, \(n\) c represents the order of the first influence polynomial of the disturbance term, represents the \(n\) c order undetermined coefficient, \(n\) d represents the order of the second influence polynomial of the disturbance term, represents the \(n\) d order undetermined coefficient, \(q\) represents the backward shift operator; When C(q)=1, and D(q)=1, and F(q)=1, the initial dynamic model is an ARX model; When B(q)=1, and C(q)=1, and F(q)=1, the initial dynamic model is an AR model; When B(q)=1, and D(q)=1, and F(q)=1, the initial dynamic model is an ARMA model.

10. The method for modeling the dynamic thermal error of an electric spindle based on timing analysis according to claim 8, characterized in that, The expression for performing the model order determination is: AIC = 2w - 2ln(L) (8) Among them, AIC represents the information complexity of the initial dynamic model, w represents the number of parameters of the initial dynamic model, and L represents the maximum likelihood function value of the initial dynamic model; The expression for performing the parameter identification is: Among them, represents the identification and estimation value of the parameter matrix, t represents the t-th moment, represents the time series data of the input and output of the initial dynamic model at the t-th moment, E(t) represents the actual thermal error data measured at the t-th moment, T represents the transpose, n is the n-th temperature sensor, and N represents the number of temperature sensors.

Images