Multi-physics field coupling simulation method for motorized spindle

Abstract

The invention belongs to the technical field of simulation analysis, and particularly relates to an electric spindle multi-physics field coupling simulation method. The motorized spindle multi-physics field coupling simulation method comprises the following steps that S100, a multi-physics field coupling model is established through finite element software, wherein the steps of geometric model importing, material setting, rotor dynamics module setting, solid heat transfer module setting and grid division are included; and S200, solving through the finite element software: carrying out multi-physics field coupling on the bearing nonlinear restoring force model, the motor unbalanced magnetic pulling force model and the main shaft thermal model with the motorized main shaft rotor model, and solving to obtain the dynamic performance. By using the simulation model for motorized spindle dynamic characteristic research in the technical scheme, the dynamic characteristic research requirement of the motorized spindle under the actual complex working condition can be met.

Description

Technical Field

[0001] The present invention belongs to the technical field of simulation analysis, and particularly relates to a multi-physical field coupling simulation method for an electric spindle. Background Art

[0002] Existing simulation models for studying the dynamic characteristics of electric spindles mostly include: bearing-rotor-thermal coupling models, bearing-rotor-motor (electromagnetic) coupling models, and motor (electromagnetic)-rotor-thermal coupling models. Under complex working conditions, the electric spindle has strong non-linear characteristics during operation, and the existing simulation models for studying the dynamic characteristics of electric spindles cannot meet the research needs of the dynamic characteristics of the rotor under actual complex working conditions.

[0003] Therefore, there is an urgent need to propose a multi-physical field coupling simulation method for an electric spindle that is closer to actual complex working conditions to solve the above problems. Summary of the Invention

[0004] The purpose of the present invention is to provide a multi-physical field coupling simulation method for an electric spindle that can realize the research of dynamic characteristics under complex working conditions. This purpose is achieved through the following technical solutions:

[0005] The first aspect of the present invention proposes a multi-physical field coupling simulation method for an electric spindle, including the following steps:

[0006] S100. Establish an electric spindle rotor model through finite element software, including: importing a geometric model, setting materials, setting a rotor dynamics module, setting a solid heat transfer module, and mesh generation. The geometric model includes a shaft structure model;

[0007] Setting the rotor dynamics module includes applying a non-linear restoring force and an unbalanced magnetic pull force to the shaft structure model;

[0008] Setting the solid heat transfer module includes applying the surface temperature T motor of the shaft caused by the motor and the surface temperature T bearing of the shaft caused by the bearing to the surface of the shaft structure model;

[0009] Solve the non-linear restoring force through a bearing non-linear restoring force model, solve the unbalanced magnetic pull force through a motor unbalanced magnetic pull force model, and solve the surface temperature T motor of the shaft caused by the motor and the surface temperature T bearing of the shaft caused by the bearing through a spindle thermal model;

[0010] S200. Solve through the finite element software, including: coupling the bearing non-linear restoring force model, the motor unbalanced magnetic pull force model, and the spindle thermal model with the electric spindle rotor model to solve for the dynamic performance.

[0011] The electro-spindle simulation model in this technical solution is a multi-physical field coupling simulation model, that is: a coupled dynamic model integrating the heat of bearings - motors - rotating shafts - spindles. This model can more accurately reflect the spindle's time-domain rotational accuracy, critical speed, and order distribution under actual complex working conditions, providing a theoretical basis for the subsequent optimization design of high-performance spindles.

[0012] In addition, the multi-physical field coupling simulation method of the electro-spindle of the present invention may also have the following additional technical features:

[0013] In some embodiments of the present invention, the geometric model further includes a tool head structure model, a tool holder structure model, and a motor rotor structure model. The tool head structure model and the rotating shaft structure model are respectively connected to both ends of the tool holder structure model, and the motor rotor structure model is sleeved on the outer periphery of the rotating shaft structure model.

[0014] In some embodiments of the present invention, the rotating shaft structure model has an axially penetrating cavity.

[0015] In some embodiments of the present invention, in S100, applying the unbalanced magnetic pull force to the rotating shaft structure model includes the following steps:

[0016] S111: Divide the motor rotor structure model into multiple axial segments along its axis on average;

[0017] S112: Establish a coordinate system with the center point of the end face far from the motor on each axial segment as the origin, and obtain the position information of the selected nodes in real time through the probe built into the finite element software;

[0018] S113: Use the position information as the input variable of the motor unbalanced magnetic pull force model;

[0019] S114: Apply each unbalanced magnetic pull force to the corresponding axial segment.

[0020] In some embodiments of the present invention, in S100, solving the bearing nonlinear restoring force model includes the following steps:

[0021] S121: Set the bearing parameters and convergence accuracy. The bearing parameters include the contact angle α between the ball and the inner and outer raceways 0 , the ball diameter D; obtain the relative radial displacement δ of the inner and outer rings of the bearing in the X-axis direction x , the relative radial displacement δ of the inner and outer rings of the bearing in the Y-axis direction y , the relative axial displacement δ of the inner and outer rings of the bearing z , the relative angular displacement θ of the inner and outer rings of the bearing around the X-axis x , the relative angular displacement θ of the inner and outer rings of the bearing around the Y-axis directiony , the X-axis and the Y-axis are perpendicular;

[0022] S122. Calculate the axial distance A between the curvature center O of the outer raceway groove of a single ball and the final position O2 of the curvature center of the inner raceway groove 1j , and calculate the radial distance A between the curvature center O of the outer raceway groove of the single ball and the final position O2 of the curvature center of the inner raceway groove 2j ;

[0023]

[0024] where BD is the distance between the curvature centers of the inner and outer raceway grooves; B = f i +f o -1, f i , f o are respectively the groove curvature radius coefficients of the inner and outer raceways; Δ 1j and Δ 2j satisfy the following relationship:

[0025]

[0026] In the formula, ψ j is the azimuth angle of the ball; R i is the radius of the locus of the curvature center of the inner raceway and satisfies the following relationship:

[0027]

[0028] In the formula, d m is the pitch diameter;

[0029] Constraint X 1j , X 2j The initial value of and substitute it into the following formula:

[0030] (A 1j -X 1j ) 2 +(A 2j -X 2j ) 2 -[(f i -0.5)D+δ ij ) 2 =0

[0031]

[0032] In the formula, X 1j , X 2j are respectively the projections of the distance between the final position of the rolling ball center and the curvature center of the outer raceway in the axial and radial directions;

[0033] According to Hertz contact theory, the contact deformation amount δ between the ball and the inner ring raceway ijAnd the contact deformation δ between the ball and the outer ring raceway oj It can be calculated by the following formula:

[0034]

[0035] Where Q represents the normal pressure of contact deformation, E', R is calculated by the following formula:

[0036]

[0037] Where, E1, E2 are the elastic moduli of the two contact bodies; ν1, ν2 are the Poisson's ratios of the two contact bodies;

[0038] R x With R y are the curvature radii of the two contact bodies on the two main planes (the ball contacts the inner raceway):

[0039]

[0040] Contact between the ball and the outer raceway:

[0041]

[0042] Where, γ * =Dcosα 0 / d m ;

[0043] Since rolling bearings are lubricated, the effect of lubricating oil film thickness should also be considered:

[0044]

[0045] Where, L ij , L oj h is the distance between the center of curvature of the inner and outer raceway grooves and the center of the ball; ic , h oc are the center oil film thickness between the ball and the inner and outer ring raceways, respectively, which can be calculated by the following formula:

[0046]

[0047] in, e is a natural constant; other unknown quantities are determined by the following formula:

[0048]

[0049] According to the positional relationship between the center of curvature of the inner and outer raceway grooves and the center of the ball, the following geometric relationship exists:

[0050]

[0051] Where, α ij , α oj are the contact angles between the ball and the inner and outer raceways, respectively, where j is the ball number, satisfying 1≤j≤Z, and Z is the total number of balls;

[0052] S123, establish a single ball subjected to centrifugal force F cj and gyroscopic torque M gj The force relationship between the inner and outer raceways under the action of and solve:

[0053]

[0054] Among them, λ ij ,λ oj are the friction coefficients between the ball and the inner and outer raceways, α ij , α oj are the contact angles between the ball and the inner and outer raceways, D is the ball diameter, Q ij , Q oj Indicates the positive pressure of the contact deformation between the ball and the inner and outer raceways; centrifugal force F cj and gyroscopic torque M gj Determined by the following formula:

[0055]

[0056] Among them, ω R is the angular velocity of the ball, ω m is the angular velocity of the ball, J is the moment of inertia of the ball, ω is the speed of the inner ring of the bearing, β is the angle between the rotation axis and the revolution axis, and m is the mass of the ball;

[0057] According to Jones's raceway control hypothesis, that is, the bearing is closer to outer raceway control when running at high speed. The following relationship can be established:

[0058]

[0059] Where γ′=D / d m , tanβ=sinα oj / (cosα oj +γ′);

[0060] The inner ring of the bearing is simultaneously subjected to the forces and moments of the shaft and the balls. By analyzing the forces acting on it, the nonlinear restoring force of a single ball can be calculated as follows:

[0061]

[0062] Where r i =f i D;

[0063] S124. Repeat steps S2 - S3 until all the balls are calculated;

[0064] S125. Sum them up to obtain the non - linear restoring force.

[0065] In some embodiments of the present invention, in S100, the solution of the motor unbalanced magnetic pull force model includes the following steps:

[0066] S131. Set the motor parameters, where the motor parameters include the initial air - gap length δ and the number of stator slots Z;

[0067] S132. Calculate the eccentric air - gap length δ(α, t):

[0068] δ(α,t)≈δ0 - rcos(α - γ)

[0069] S133. Calculate the air - gap permeance Λ(α, t):

[0070]

[0071] where Λ n is as follows:

[0072]

[0073] S134. Calculate the motor magnetomotive force F(α, t);

[0074]

[0075] where F s is the fundamental - wave magnetomotive force of the stator winding, ω0 is the electrical frequency; I max is the maximum value of the current; N is the number of turns in series per phase winding; p is the number of pole pairs of the motor; is the power - factor angle; k w is the stator - winding coefficient, and the expression is as follows:

[0076]

[0077] In the formula, q is the number of slots per pole per phase; α is the slot - pitch angle related to the number of poles respectively, and the expression is as follows:

[0078]

[0079] where z is the number of stator slots;

[0080]

[0081] In the formula, is the winding pitch ratio;

[0082] S135. Calculate the air - gap magnetic density B(α, t):

[0083] B(α,t)=F j (α,t)Λ(a,t)

[0084] S136. Calculate Maxwell stress σ:

[0085]

[0086] S137. Calculate the unbalanced magnetic pull:

[0087]

[0088] Among them, f1, f2, f3, and f4 satisfy the following relationship:

[0089]

[0090] In some embodiments of the present invention, solving the spindle thermal model includes solving a bearing thermal model, a motor thermal model, and a heat transfer model.

[0091] In some embodiments of the present invention, solving the bearing thermal model includes the following steps:

[0092] Calculate the frictional heat P between the ball and the inner and outer raceways spin :

[0093]

[0094] In the formula, μ si is the coefficient of rotational friction between the ball and the inner raceway. For ball bearings, ω spin is the angular velocity of the ball's spin; Q jtot is the contact recovery force;

[0095]

[0096] Calculate the frictional heat P generated by the balls and cage ca :

[0097]

[0098] Where D b Ball diameter; D m The bearing pitch circle is equal to the average value of the inner and outer diameters of the bearing D m =(D o +D i ) / 2; α0 initial contact angle; m cage is the cage mass; μ sc is the friction factor between the ball and the cage; ω cage is the cage angular velocity, ω cage =ωD i / (D i +D o );

[0099] Calculate the total heat generation power P of the bearing b :

[0100] P b =n b (P spin +P ca ).

[0101] In some embodiments of the present invention, solving the motor thermal model includes the following steps:

[0102] Calculation of copper loss P of asynchronous motor c :

[0103] P c =P c1 +P c2

[0104] Stator copper loss P c1 and rotor copper loss P c2 Calculated by the following formula:

[0105]

[0106]

[0107] Where, I1 is the stator winding current; R1 is the winding resistance; I2 is the rotor winding current; R2 is the rotor bar resistance;

[0108] Calculate the iron loss P of the asynchronous motor Fe :

[0109] P Fe =P Fet +P Fej

[0110] Stator (rotor) yoke iron loss P Fej Calculated by the following formula:

[0111] p Fej =K a p hej G j

[0112] Where K a is the empirical coefficient, when the capacity P N <100kV□A, K a =1.59, otherwise 1.3; G j is the weight of the yoke; p hej is the yoke loss coefficient, which is expressed as follows:

[0113]

[0114] Where P 10 / 50 To be B j =1T, f = 50Hz, the loss per unit weight of silicon steel sheet, in W / kg; B j0 is the magnetic induction intensity of the unloaded yoke;

[0115] Stator (rotor) tooth iron loss P Fet Calculated by the following formula:

[0116] p Fet =K a p het G t

[0117] Where G t is the weight of the yoke; K a is the empirical coefficient, which is 1.8 for induction motors; p het is the tooth loss coefficient;

[0118]

[0119] Where P 10 / 50 To be B j =1T, f = 50Hz, the loss per unit weight of silicon steel sheet, in W / kg; B t0 is the magnetic induction intensity of the unloaded tooth;

[0120] Calculate the rotor and air friction loss P w :

[0121]

[0122] In the formula, f r D ro is the diameter of the rotor outer ring; L r is the axial length of the rotor; η air is the dynamic viscosity of air; h g is the average air gap length;

[0123] Calculation of stray losses P of asynchronous motors s :

[0124]

[0125] In the formula, copper conductor Take 0.005, P out is the output power;

[0126] Calculate the total heating power P:

[0127] P=P cu +P Fe +Pw +P s .

[0128] In some embodiments of the present invention, solving the heat transfer model comprises the following steps:

[0129] Calculate the convective heat transfer P from the stator (rotor) core to the air gap 1(2) :

[0130]

[0131] Where σ SB =5.67×10 -8 is the Stefan-Boltzmann constant; ε thr is the relative emissivity, usually taken as 0.85; T2 is the temperature of the heat outflow surface; T1 is the temperature of the heat absorption surface; S is the surface area of the heat sink;

[0132] The surface temperature rise T of the stator (rotor) core is calculated by the following formula:

[0133]

[0134] Where P Fe is the stator (rotor) iron loss; O is the stator (rotor) heat dissipation area, k is the heat transfer coefficient per unit surface area to the external air;

[0135] Calculate the convective heat transfer P3 from the stator (rotor) core to the air gap:

[0136] P3=α th |T1-T2|S

[0137] Heat transfer coefficient α from stator (rotor) core to air gap th Calculated by the following formula:

[0138]

[0139] Where λ is the thermal conductivity of air. At room temperature (20°C), the thermal conductivity of air is 0.0267W / mK. At 0°C, the thermal conductivity of air is 0.0251W / mK. At 100°C, the thermal conductivity of air is 0.0321W / mK.

[0140] The Nusselt number Nu is calculated by the following formula:

[0141]

[0142] Where, Ta≈Ta m ;

[0143] Taylor number Ta is calculated by the following formula:

[0144]

[0145] where ρ is the mass density of the fluid; Ω is the angular velocity of the rotor; r m is the average radius of the stator and rotor; δ is the air-gap length; μ is the dynamic viscosity of the fluid;

[0146] Calculate the convective heat transfer amount P4 from the end of the rotor core to the air:

[0147] P4 = α th |T1 - T2|S

[0148] The natural convection heat transfer coefficient α th is calculated by the following formula:

[0149]

[0150] where T1 is the temperature of the heat outflow surface; T2 is the temperature of the heat absorption surface; D is the outer diameter of the rotor;

[0151] Calculate the convective heat transfer amount P5 from the stator core to the coolant:

[0152]

[0153] The heat transfer coefficient α from the stator core to the coolant th is calculated by the following formula:

[0154]

[0155] where v is the velocity of the coolant; l is the length of the motor body;

[0156] Calculate the convective heat transfer amount P from the stator core to the air 6(7) :

[0157] P 6(7) = α th |T1 - T2|S

[0158] The natural convection heat transfer coefficient α th is calculated by the following formula:

[0159]

[0160] where T1 is the temperature of the radiation surface; T2 is the temperature of the absorption surface; D is the outer diameter of the stator;

[0161] Calculate the thermal radiation amount from the stator core to the air

[0162]

[0163] where σ SB = 5.67×10-8 is the Stefan-Boltzmann constant; ε thr is the relative emissivity, usually taken as 0.85; T1 is the temperature of the radiation surface; T2 is the temperature of the absorption surface; S is the surface area of the radiator;

[0164] Calculate the heat dissipation P of the bearing bs :

[0165]

[0166] Calculate the remaining heat W:

[0167]

[0168] Calculate the surface temperature T of the rotating shaft under the action of the motor motor :

[0169]

[0170] In the formula, O motor is the surface area of the rotating shaft under the action of the heat generated by the motor;

[0171] Calculate the surface temperature T of the rotating shaft under the action of the bearing bearing :

[0172]

[0173] In the formula, O bearing is the surface area of the rotating shaft under the action of the heat generated by the bearing. Description of the drawings

[0174] In the drawings:

[0175] Figure 1 Schematically shows the structural schematic diagram of the geometric model according to the embodiment of the present invention;

[0176] Figure 2 Schematically shows the cross-sectional view of the rotating shaft structure model according to the embodiment of the present invention;

[0177] Figure 3 Schematically shows the solution flow chart of the multi-physical field coupling model of the motorized spindle according to the embodiment of the present invention;

[0178] Figure 4 Schematically shows the position relationship diagram of the curvature centers of the inner and outer raceways of the bearing and the center of the ball according to the embodiment of the present invention;

[0179] Figure 5 Schematically shows the force diagram of the bearing balls according to the embodiment of the present invention;

[0180] Figure 6Schematically shows a dual-cylinder model of the stator and rotor of a motor considering eccentricity according to an embodiment of the present invention;

[0181] Figure 7 Schematically shows a flow chart for solving the bearing thermal model according to an embodiment of the present invention;

[0182] Figure 8 Schematically shows a flow chart for solving the motor thermal model according to an embodiment of the present invention;

[0183] Figure 9 Schematically shows a graph of the non-linear restoring force of the front bearing varying with time according to an embodiment of the present invention;

[0184] Figure 10 Schematically shows a graph of the non-linear restoring force of the rear bearing varying with time according to an embodiment of the present invention;

[0185] Figure 11 Schematically shows a graph of the unbalanced magnetic pull force varying with time according to an embodiment of the present invention;

[0186] Figure 12 Schematically shows a graph of the temperature of the front bearing varying with time according to an embodiment of the present invention;

[0187] Figure 13 Schematically shows a graph of the temperature of the rear bearing varying with time according to an embodiment of the present invention;

[0188] Figure 14 Schematically shows a graph of the temperature of the motor varying with time according to an embodiment of the present invention;

[0189] Figure 15 Schematically shows a graph of the radial displacement of the rotating shaft varying with time according to an embodiment of the present invention;

[0190] Figure 16 Schematically shows a Campbell diagram of the dynamic performance analysis results of Example 1;

[0191] Figure 17 Schematically shows a locus diagram of the dynamic performance analysis results of Example 1;

[0192] Figure 18 Schematically shows a spectrogram of the dynamic performance analysis results of Example 1.

[0193] The reference numerals in the drawings are represented as follows:

[0194] 100, tool head structure model; 200, tool holder structure model; 300, rotating shaft structure model; 301, cavity; 400, motor rotor structure model. Detailed implementation mode

[0195] Figure 1 The structural schematic diagram of the geometric model according to the embodiment of the present invention is schematically shown.

[0196] Figure 2 The structural schematic diagram of the rotating shaft structure model according to the embodiment of the present invention is schematically shown.

[0197] Figure 3 The solution flow chart of the multi-physical field coupling model of the motorized spindle according to the embodiment of the present invention is schematically shown. As Figures 1 to 3 shown, the present invention proposes a multi-physical field coupling simulation method for a motorized spindle, including the following steps:

[0198] S100. Establish a motorized spindle rotor model through finite element software, including: importing a geometric model, setting materials, setting a rotor dynamics module, setting a solid heat transfer module, and meshing. The geometric model includes a rotating shaft structure model 300;

[0199] Set the rotor dynamics module, including applying a nonlinear restoring force (Nonlinear Restoring Force, NRF) and an unbalanced magnetic pulling force (Unbalanced Magnetic Pulling Force, UMP) to the rotating shaft structure model 300;

[0200] Set the solid heat transfer module, including applying the surface temperature T of the motor acting on the rotating shaft motor and the surface temperature T of the bearing acting on the rotating shaft bearing to the surface of the rotating shaft structure model 300;

[0201] Solve the nonlinear restoring force through the bearing nonlinear restoring force model, solve the unbalanced magnetic pulling force through the motor unbalanced magnetic pulling force model, and solve the surface temperature T of the motor acting on the rotating shaft through the spindle thermal model motor and the surface temperature T of the bearing acting on the rotating shaft bearing ;

[0202] S200. Solve through finite element software, including: coupling the bearing nonlinear restoring force model, the motor unbalanced magnetic pulling force model, and the spindle thermal model with the motorized spindle rotor model to solve for the dynamic performance.

[0203] It can be understood that in the finite element software, the rotating shaft structure model 300 is equivalent to the rotating shaft in the actual working condition. As Figure 3As shown, the nonlinear restoring force is calculated through the bearing nonlinear restoring force model and applied to the front and rear rotating shafts. At the same time, after eccentricity occurs during the rotation of the rotating shaft, it acts on the bearing nodes, thereby affecting the bearing nonlinear restoring force model. Similarly, the unbalanced magnetic pull force model is calculated through the motor unbalanced magnetic pull force model and applied to the rotating shaft. At the same time, after eccentricity occurs during the rotation of the rotating shaft, it acts on the motor unbalanced magnetic pull force model. As Figure 4 and Figure 5 shown, through the spindle thermal model, the motor temperature T motor and the bearing temperature T bearing are applied to the surface of the rotating shaft. The electric spindle simulation model in this technical solution is a multi-physics field coupling simulation model, that is: a coupled dynamics model integrating bearings - motors - rotating shafts - spindle heat is integrated at the same time. This model can more accurately reflect the spindle time-domain rotation accuracy, critical speed, and order distribution, etc. of the spindle under actual complex working conditions, providing a theoretical basis for the subsequent optimization design of high-performance spindles.

[0204] Furthermore, referring to Figure 1 , the geometric model further includes a tool head structure model 100, a tool holder structure model 200, and a motor rotor structure model 400. The tool head structure model 100 and the rotating shaft structure model 300 are respectively connected to both ends of the tool holder structure model 200, and the motor rotor structure model 400 is sleeved on the outer periphery of the rotating shaft structure model 300.

[0205] In the prior art, the geometric model is only a rotating shaft model without other additional structures, which cannot meet the research requirements of the dynamic characteristics of the rotating shaft under actual complex working conditions. The geometric model in this technical solution includes a tool holder, a tool head, and a motor rotor at the same time, and can accurately reflect the dynamic performance of the spindle during the actual complex working process. Optionally, the tool head can be a standard rod or other cutting tools.

[0206] Furthermore, referring to Figure 2 , the rotating shaft structure model 300 has an axially penetrating cavity 301.

[0207] In the prior art, the geometric model is usually a solid rotating shaft model, but in actual situations, the rotating shaft is a hollow structure. Therefore, this technical solution sets an axially penetrating cavity 301 in the rotating shaft structure model 300, which is closer to the actual working conditions, so as to obtain more accurate analysis results.

[0208] Furthermore, continuing to refer to Figure 1 , in S100, applying the unbalanced magnetic pull force to the rotating shaft structure model 300 includes the following steps:

[0209] S111. Divide the motor rotor structure model 400 into multiple axial segments along its axis on average;

[0210] S112, establishing a coordinate system with the center point of the end face of each shaft segment away from the motor as the origin, and obtaining position information of the selected node in real time using a probe built into the finite element software;

[0211] S113, using the position information as an input variable of a motor unbalanced magnetic pull model;

[0212] S114. Apply each unbalanced magnetic pull to the corresponding shaft segment.

[0213] Among them, the probe is a built-in function of the finite element software, which can be used to obtain the position information of the selected nodes of the finite element model in real time.

[0214] During actual operation, the radial eccentricity of the rotating shaft will vary due to different axial positions. However, in the prior art, the unbalanced magnetic pull modeling method only applies a total unbalanced magnetic pull force to the overall structure of the motor rotor, making it difficult to reflect the actual unbalanced magnetic pull distribution at different axial positions along the rotating shaft. In the present technical solution, by dividing the motor rotor portion into multiple parts along the axial direction and applying the unbalanced magnetic pull force to the corresponding motor rotor shaft segments respectively, the difference in unbalanced magnetic pull force along the axial direction of the rotating shaft under different eccentricities can be reflected, thereby accurately reflecting the dynamic performance of the rotating shaft in actual complex working conditions. For example, in the present technical solution, the motor rotor portion is divided into four parts, and the length of each part along the axial direction of the geometric model is equal. In other embodiments, it can be divided into three, five, six, etc. parts as needed, and the length of each part can be equal or different.

[0215] Furthermore, in this technical solution, the finite element software used in the model building and solution process is COMSOL Multiphysics.

[0216] COMSOL Multiphysics is an industry-leading multiphysics simulation platform, providing the ability to simulate a single physics field and flexibly couple multiple physics fields. This allows for precise analysis of equipment, processes, and procedures across various engineering fields. The software's built-in Model Builder includes a complete modeling workflow, facilitating all simulation steps, from geometric modeling, material parameter and physics field setup, to solution and result processing.

[0217] Furthermore, if Figure 3 As shown, in S100, solving the bearing nonlinear restoring force model includes the following steps:

[0218] S121, set bearing parameters and convergence accuracy, bearing parameters include the contact angle α between the ball and the inner and outer raceways 0 , ball diameter D; obtain the relative radial displacement δ of the inner and outer rings of the bearing in the X-axis direction x , the relative radial displacement of the inner and outer rings of the bearing in the Y-axis direction δ y, the relative axial displacement δ between the inner and outer rings of the bearing z , the relative angular displacement θ of the inner and outer rings of the bearing about the X-axis x , the relative angular displacement θ of the inner and outer rings of the bearing about the Y-axis y , the X-axis and the Y-axis are perpendicular;

[0219] S122. As Figure 4 shown, calculate the axial distance A between the curvature center O of the outer raceway groove of a single ball and the final position O2 of the curvature center of the inner raceway groove 1j , calculate the radial distance A between the curvature center O of the outer raceway groove of the single ball and the final position O2 of the curvature center of the inner raceway groove 2j ;

[0220]

[0221] where BD is the distance between the curvature centers of the inner and outer raceway grooves; B = f i +f o , f i , f o are the groove curvature radius coefficients of the inner and outer raceways respectively; Δ 1j and Δ 2j satisfy the following relationship:

[0222]

[0223] In the formula, ψ j is the azimuth angle of the ball; R i is the radius of the locus of the curvature center of the inner raceway and satisfies the following relationship:

[0224]

[0225] In the formula, d m is the pitch diameter;

[0226] Constrain X 1j , X 2j with its initial value and substitute it into the following formula:

[0227]

[0228] In the formula, X 1j , X 2j are the projections of the distance between the final position of the ball center and the curvature center of the outer raceway in the axial and radial directions respectively;

[0229] According to Hertz contact theory, the contact deformation δ ij of the ball with the inner ring raceway and the contact deformation δ oj of the ball with the outer ring raceway can be calculated by the following formula:

[0230]

[0231] In the formula, Q represents the normal pressure of the contact deformation, E', R are calculated by the following formulas respectively:

[0232]

[0233] In the formula, E1 and E2 are the elastic moduli of the two contacting bodies; ν1 and ν2 are the Poisson's ratios of the two contacting bodies;

[0234] R x and R y are the radii of curvature of the two contacting bodies on the two principal planes (contact between the ball and the inner raceway):

[0235]

[0236] For the contact between the ball and the outer raceway:

[0237]

[0238] In the formula, γ * = Dcosα 0 / d m ;

[0239] Since the rolling bearing is lubricated, the influence of the lubricating oil film thickness should also be considered:

[0240]

[0241] In the formula, L ij , L oj are the distances between the curvature centers of the inner and outer raceway grooves and the center of the ball respectively; h ic , h oc are the central oil film thicknesses between the ball and the inner and outer raceway grooves respectively, and can be calculated by the following formula:

[0242]

[0243] Among them, e is the natural constant; other unknowns are determined by the following formula:

[0244]

[0245] As Figure 4 shown, there are also the following geometric relationships according to the positional relationship between the curvature centers of the inner and outer raceway grooves and the center of the ball:

[0246]

[0247] In the formula, α ij , α ojare the contact angles between the ball and the inner and outer raceways respectively, where j is the ball number, satisfying 1 ≤ j ≤ Z, and Z is the total number of balls;

[0248] S123. As Figure 5 shown, establish the force relationship between a single ball and the inner and outer raceways under the action of centrifugal force F cj and gyroscopic moment M gj and solve it:

[0249]

[0250] where λ ij , λ oj are the friction coefficients between the ball and the inner and outer raceways respectively, α ij , α oj are the contact angles between the ball and the inner and outer raceways respectively, D is the ball diameter, Q ij , Q oj represent the normal pressures of the contact deformation between the ball and the inner and outer raceways; the centrifugal force F cj and the gyroscopic moment M gj are determined by the following formula:

[0251]

[0252] where ω R is the angular velocity of the ball's rotation about its own axis, ω m is the angular velocity of the ball's revolution, J is the moment of inertia of the ball, ω is the rotational speed of the inner ring of the bearing, β is the angle between the self-rotation axis and the revolution axis, and m is the mass of the ball;

[0253] According to Jones' raceway control hypothesis, that is: when the bearing is operating at high speed, it is closer to the outer raceway control, the following relational expression can be established:

[0254]

[0255] where γ′ = D / d m , tanβ = sinα oj / (cosα oj +γ′);

[0256] The inner ring of the bearing is simultaneously subjected to the forces and torques of the rotating shaft and the balls. By analyzing its forces, the calculation formula for the non-linear restoring force of a single ball can be obtained:

[0257]

[0258] In the formula, r i = f i D;

[0259] S124. Repeat steps S2 - S3 until all the balls are calculated;

[0260] S125. Sum and obtain the nonlinear restoring force.

[0261] Furthermore, if Figure 6 As shown, the motor unbalanced magnetic pull model adopts the motor stator and rotor double cylinder model. In S100, solving the motor unbalanced magnetic pull model includes the following steps:

[0262] S131, setting motor parameters, including initial air gap length δ and number of stator slots Z;

[0263] S132. Calculate the eccentric air gap length δ(α, t):

[0264] δ(α,t)≈δ0-rcos(α-γ)(17)

[0265] S133. Calculate the air gap permeability Λ(α, t):

[0266]

[0267] Among them, Λ n as follows:

[0268]

[0269] S134, calculating the motor magnetic potential F(α, t);

[0270]

[0271] Among them, F s is the fundamental magnetic potential of the stator winding, ω0 is the electrical frequency; I max is the maximum current; N is the number of series turns of each phase winding; p is the number of motor pole pairs; is the power factor angle; k w is the stator winding coefficient, and its expression is as follows:

[0272] k w =k d k p (twenty one)

[0273]

[0274] Where q is the number of slots per pole and per phase; α is the slot pitch angle related to the number of poles, and the expression is as follows:

[0275]

[0276] Where z is the number of stator slots;

[0277]

[0278] In the formula, is the winding pitch ratio;

[0279] S135. Calculate the air gap magnetic flux density B(α, t):

[0280] B(α,t)=F j (α,t)Λ(a,t) (25)

[0281] S136. Calculate Maxwell stress σ:

[0282]

[0283] S137. Calculate the unbalanced magnetic pull:

[0284]

[0285] Among them, f1, f2, f3, and f4 satisfy the following relationship:

[0286]

[0287] Furthermore, solving the spindle thermal model includes solving the bearing thermal model, the motor thermal model and the heat transfer model.

[0288] Furthermore, if Figure 7 As shown in Figure 2, the solution of the bearing thermal model includes the following steps:

[0289] Calculate the frictional heat P between the ball and the inner and outer raceways spin :

[0290]

[0291] In the formula, μ si is the coefficient of rotational friction between the ball and the inner raceway. For ball bearings, ω spin is the angular velocity of the ball's spin; Q jtot is the contact recovery force;

[0292]

[0293] Calculate the frictional heat P generated by the balls and cage ca :

[0294]

[0295] Where D b is the ball diameter; D m The bearing pitch circle is equal to the average value of the inner and outer diameters of the bearing D m =(D o +D i ) / 2; α0 is the initial contact angle; m cage is the cage mass; μ scis the friction factor between the ball and the cage; ω cage is the cage angular velocity, ω cage =ωD i / (D i +D o );

[0296] Calculate the total heat generation power P of the bearing b :

[0297] P b =n b (P spin +P ca ) (32)

[0298] Furthermore, if Figure 8 As shown in Figure 2, solving the motor thermal model includes the following steps:

[0299] Calculation of copper loss P of asynchronous motor c :

[0300] P c =P c1 +P c2 (33)

[0301] Stator copper loss P c1 and rotor copper loss P c2 Calculated by the following formula:

[0302]

[0303] Where, I1 is the stator winding current; R1 is the winding resistance; I2 is the rotor winding current; R2 is the rotor bar resistance;

[0304] Calculate the iron loss P of the asynchronous motor Fe :

[0305] P Fe =P Fet +P Fej (36)

[0306] Stator (rotor) yoke iron loss P Fej Calculated by the following formula:

[0307] p Fej =K a p hej G j (37)

[0308] Where K a is the empirical coefficient, when the capacity When K a =1.59, otherwise 1.3; G j is the weight of the yoke; p hejis the yoke loss coefficient, which is expressed as follows:

[0309]

[0310] Where P 10 / 50 To be B j =1T, f = 50Hz, the loss per unit weight of silicon steel sheet, in W / kg; B j0 is the magnetic induction intensity of the unloaded yoke;

[0311] Stator (rotor) tooth iron loss P Fet Calculated by the following formula:

[0312] p Fet =K a p het G t (39)

[0313] Where G t Yoke weight; K a is the empirical coefficient, which is 1.8 for induction motors; p het Tooth loss coefficient;

[0314]

[0315] Where P 10 / 50 To be B j =1T, f = 50Hz, the loss per unit weight of silicon steel sheet, in W / kg, B t0 is the magnetic induction intensity of the unloaded tooth;

[0316] Calculate the rotor and air friction loss P w :

[0317]

[0318] In the formula, f r D ro is the diameter of the rotor outer ring; L r is the axial length of the rotor; η air is the dynamic viscosity of air; h g is the average air gap length;

[0319] Calculation of stray losses P of asynchronous motors s :

[0320]

[0321] In the formula, copper conductor Take 0.005, P out is the output power;

[0322] Calculate the total heating power P:

[0323] P=P cu +P Fe +P w +P s (43)

[0324] Furthermore, the solution of the heat transfer model includes the following steps:

[0325] Calculate the convective heat transfer P from the stator (rotor) core to the air gap 1(2) :

[0326]

[0327] In the formula, is the Stefan-Boltzmann constant; ε thr is the relative emissivity, usually taken as 0.85; T2 is the temperature of the heat outflow surface; T1 is the temperature of the heat absorption surface; S is the surface area of the heat sink;

[0328] The surface temperature rise T of the stator (rotor) core is calculated by the following formula:

[0329]

[0330] Where P Fe is the stator (rotor) iron loss; O is the stator (rotor) heat dissipation area, k is the heat transfer coefficient per unit surface area to the external air;

[0331] Calculate the convective heat transfer P3 from the stator (rotor) core to the air gap:

[0332] P3=α th |T1-T2|S (46)

[0333] Heat transfer coefficient α from stator (rotor) core to air gap th Calculated by the following formula:

[0334]

[0335] Where λ is the thermal conductivity of air. At room temperature (20°C), the thermal conductivity of air is 0.0267W / mK. At 0°C, the thermal conductivity of air is 0.0251W / mK. At 100°C, the thermal conductivity of air is 0.0321W / mK.

[0336] The Nusselt number Nu is calculated by the following formula:

[0337]

[0338] Where, Ta≈Ta m ;

[0339] Taylor number Ta is calculated by the following formula:

[0340]

[0341] where ρ is the mass density of the fluid; Ω is the angular velocity of the rotor; r m is the average radius of the stator and the rotor; δ is the air-gap length; μ is the dynamic viscosity of the fluid;

[0342] Calculate the convective heat transfer amount P4 from the end of the rotor core to the air:

[0343] P4 = α th |T1 - T2|S(50)

[0344] The natural convection heat transfer coefficient α th is calculated by the following formula:

[0345]

[0346] where T1 is the temperature of the heat outflow surface; T2 is the temperature of the heat absorption surface; D is the outer diameter of the rotor;

[0347] Calculate the convective heat transfer amount P5 from the stator core to the coolant:

[0348]

[0349] The heat transfer coefficient α from the stator core to the coolant th is calculated by the following formula:

[0350]

[0351] where v is the velocity of the coolant; l is the length of the motor body;

[0352] Calculate the convective heat transfer amount P from the stator core to the air 6(7) :

[0353] P 6(7) = α th |T1 - T2|S(54)

[0354] The natural convection heat transfer coefficient α th is calculated by the following formula:

[0355]

[0356] where T1 is the temperature of the radiation surface; T2 is the temperature of the absorption surface; D is the outer diameter of the stator;

[0357] Calculate the thermal radiation amount from the stator core to the air

[0358] φ th = ε thr σ SB(T1 4 -T2 4 )S (56)

[0359] In the formula, is the Stefan-Boltzmann constant; ε thr is the relative emissivity, usually taken as 0.85; T1 is the radiating surface temperature; T2 is the absorbing surface temperature; S is the surface area of the radiator;

[0360] Calculation of bearing heat dissipation P bs :

[0361]

[0362] Calculate the remaining heat W:

[0363]

[0364] Calculate the surface temperature T of the motor shaft motor :

[0365]

[0366] Where, O motor The surface area of the shaft affected by the heat generated by the motor;

[0367] Calculate the surface temperature T of the bearing shaft bearing :

[0368]

[0369] Where, O bearing The surface area of the shaft to which the heat generated by the bearing acts.

[0370] Example 1

[0371] In this embodiment, the shaft material is structural steel, and the rotation speed is 15000 rpm. The simulation results show the nonlinear restoring force of the front and rear bearings, the unbalanced magnetic pull, and the front and rear bearing temperatures T bearing , motor temperature T motor And the radial component of the axis of rotation is Figures 9 to 15 The dynamic performance analysis results are shown in Figures 16 to 18 shown.

[0372] The above are merely preferred embodiments of the present invention, but the scope of protection of the present invention is not limited thereto. Any changes or substitutions that can be easily conceived by a person skilled in the art within the technical scope disclosed in the present invention should be included in the scope of protection of the present invention. Therefore, the scope of protection of the present invention should be based on the scope of protection of the claims.

Claims

1. A multi-physical field coupling simulation method for an electric spindle, characterized in that, It includes the following steps: S100. Establish an electric spindle rotor model through finite element software, including: importing a geometric model, setting materials, setting a rotor dynamics module, setting a solid heat transfer module, and mesh generation. The geometric model includes a shaft structure model (300); Set the rotor dynamics module, including applying a nonlinear restoring force and an unbalanced magnetic pull force to the shaft structure model (300); Set the solid heat transfer module, including acting on the surface of the shaft structure model (300) the surface temperature T of the shaft acted on by the motor motor and the surface temperature T of the shaft acted on by the bearing bearing ​ Solve the non-linear restoring force through the bearing non-linear restoring force model, solve the unbalanced magnetic pull force through the motor unbalanced magnetic pull force model, and solve the surface temperature T of the rotating shaft affected by the motor through the spindle thermal model motor and the surface temperature T of the rotating shaft affected by the bearing bearing ; S200. Solve through the finite element software, including: performing multi-physics field coupling of the bearing nonlinear restoring force model, the motor unbalanced magnetic pull force model, and the spindle thermal model with the electric spindle rotor model to solve for the dynamic performance.

2. The multi-physical field coupling simulation method of the electric spindle according to claim 1, wherein The geometric model further includes a tool head structure model (100), a tool holder structure model (200), and a motor rotor structure model (400). The tool head structure model (100) and the shaft structure model (300) are respectively connected to both ends of the tool holder structure model (200), and the motor rotor structure model (400) is sleeved on the outer periphery of the shaft structure model (300).

3. The multi-physical field coupling simulation method of the motorized spindle according to claim 1, characterized in that The shaft structure model (300) has an axially penetrating cavity (301).

4. The multi-physical field coupling simulation method of the electric spindle according to claim 1, wherein, In S100, applying the unbalanced magnetic pull force to the shaft structure model (300) includes the following steps: S111. Divide the motor rotor structure model (400) into multiple axial segments along its axis on average; S112. Establish a coordinate system with the center point of the end face of each axial segment away from the motor as the origin, and obtain the position information of selected nodes in real time through the probe built in the finite element software; S113. Use the position information as the input variable of the motor unbalanced magnetic pull force model; S114. Apply each unbalanced magnetic pull force to the corresponding axial segment.

5. The multi-physical field coupling simulation method of the electric spindle according to claim 1, wherein In S100, the solution of the bearing nonlinear restoring force model includes the following steps: S121. Set the bearing parameters and convergence accuracy. The bearing parameters include the contact angle α between the ball and the inner and outer raceways 0 , the ball diameter D; obtain the relative radial displacement δ of the inner and outer rings of the bearing in the X-axis direction x , the relative radial displacement δ of the inner and outer rings of the bearing in the Y-axis direction y , the relative axial displacement δ of the inner and outer rings of the bearing z , the relative angular displacement θ of the inner and outer rings of the bearing about the X-axis x , the relative angular displacement θ of the inner and outer rings of the bearing about the Y-axis direction y , the X-axis and the Y-axis are perpendicular; S122. Calculate the axial distance A between the curvature center O of the outer raceway groove of a single ball and the final position O2 of the curvature center of the inner raceway groove. 1j , calculate the radial distance A between the curvature center O of the outer raceway groove of the single ball and the final position O2 of the curvature center of the inner raceway groove. 2j ; Among them, BD is the distance between the curvature centers of the inner and outer raceway grooves; B = f i + f o - 1, f i , f o are the groove curvature radius coefficients of the inner and outer raceways respectively; Δ 1j and Δ 2j satisfy the following relationship: where ψ j is the azimuth angle of the ball; R i is the radius of the locus of the center of curvature of the inner raceway, and satisfies the following relationship: where d m is the pitch diameter; Constraint X 1j , X 2j Take the initial value of and substitute it into the following formula: (A 1j -X 1j ) 2 +(A 2j -X 2j ) 2 -[(f i -0.5)D+δ ij ) 2 = 0 where X 1j and X 2j are respectively the projections of the distance between the final position of the center of the rolling ball and the center of curvature of the outer raceway in the axial and radial directions; According to Hertz contact theory, the contact deformation δ between the ball and the inner raceway ij and the contact deformation δ between the ball and the outer raceway oj can be calculated by the following formula: where Q represents the normal pressure of contact deformation, E', R are calculated respectively by the following formulas: In the formula, E1 and E2 are the elastic moduli of the two contacting bodies; ν1 and ν2 are the Poisson's ratios of the two contacting bodies; R x With R y are the radii of curvature of the two contacting bodies on two main planes (contact between the ball and the inner raceway): For the contact between the ball and the outer raceway: where γ * = D cos α 0 / d m ; Since the rolling bearing is lubricated, the influence of the lubricating oil film thickness should also be considered: where L ij , L oj are the distances from the curvature centers of the inner and outer raceway grooves to the center of the ball respectively; h ic , h oc are the central oil film thicknesses between the ball and the inner and outer raceway grooves respectively, and can be calculated by the following formula: Among them, e is the natural constant; other unknowns are determined by the following formula: There are also the following geometric relationships according to the position relationship between the curvature centers of the inner and outer raceway grooves and the center of the ball: where α ij , α oj are the contact angles between the balls and the inner and outer raceways respectively, where j is the ball number, satisfying 1 ≤ j ≤ Z, and Z is the total number of balls; S123. Establish the force relationship between a single ball and the inner and outer raceways under the action of centrifugal force F cj and gyroscopic moment M gj and solve it: where, λ ij , λ oj are the friction coefficients between the balls and the inner and outer raceways respectively, α ij , α oj are the contact angles between the balls and the inner and outer raceways respectively, D is the ball diameter, Q ij , Q oj represent the normal pressures of the contact deformation between the balls and the inner and outer raceways; the centrifugal force F cj and the gyroscopic moment M gj are determined by the following formula: where ω R is the angular velocity of the ball's self-rotation, ω m is the angular velocity of the ball's revolution, J is the moment of inertia of the ball, ω is the rotational speed of the inner ring of the bearing, β is the angle between the self-rotation axis and the revolution axis, and m is the mass of the ball; According to Jones' raceway control hypothesis, that is: the bearing is closer to the outer raceway control during high-speed operation. The following relationship can be established: where γ′ = D / d m , tanβ = sinα oj / (cosα oj + γ′); The inner ring of the bearing is simultaneously subjected to the forces and torques of the shaft and the ball. By performing a force analysis on it, the calculation formula for the nonlinear restoring force of a single ball can be obtained: where r i = f i D; S124. Repeat steps S2 - S3 until all the balls are calculated; S125. Sum to obtain the nonlinear restoring force.

6. The multi-physical field coupling simulation method of the electric spindle according to claim 1, characterized in that In S100, the solution of the motor unbalanced magnetic pull force model includes the following steps: S131. Set the motor parameters, and the motor parameters include the initial air gap length δ and the number of stator slots Z; S132. Calculate the eccentric air gap length δ(α, t): δ(α,t)≈δ0 - rcos(α - γ) S133. Calculate the air gap permeance Λ(α, t): Among them, Λ n are as follows: S134. Calculate the motor magnetomotive force F(α, t); Among them, F s is the fundamental magnetic potential of the stator winding, ω0 is the electrical frequency; I max is the maximum value of the current; N is the number of turns in series of each phase winding; p is the number of pole pairs of the motor; is the power factor angle; k w is the stator winding coefficient, and the expression is as follows: k w = k d k p Wherein, q is the number of slots per pole per phase; α is the slot pitch angle related to the number of poles, and the expression is as follows: Wherein, z is the number of stator slots; In the formula, is the winding pitch ratio; S135. Calculate the air-gap magnetic density B(α, t): B(α,t) = F j (α,t)Λ(a,t) S136. Calculate the Maxwell stress σ: S137. Calculate the unbalanced magnetic pull: Wherein, f1, f2, f3, f4 satisfy the following relationship:

7. The multi-physical field coupling simulation method of the electric spindle according to claim 1, characterized in that The solution of the spindle thermal model includes solving the bearing thermal model, the motor thermal model and the heat transfer model.

8. The multi-physical field coupling simulation method of the electric spindle according to claim 7, wherein The solution of the bearing thermal model includes the following steps: Calculate the heat generated by friction between the ball and the inner and outer raceways P spin : where μ si is the rotational friction coefficient between the ball and the inner raceway. For ball bearings, ω spin is the spin angular velocity of the ball; Q jtot is the contact restoring force; Calculate the heat generated by the friction between the ball and the cage P ca : Where D b is the ball diameter; D m is the bearing pitch circle, whose value is equal to the average value of the inner and outer diameters of the bearing D m =(D o +D i ) / 2; α0 is the initial contact angle; m cage is the cage mass; μ sc is the friction coefficient between the ball and the cage; ω cage is the cage angular velocity, ω cage =ωD i / (D i +D o ); Calculate the total heat generation power P of the bearing b : P b = n b (P spin + P ca ).

9. The multi-physical field coupling simulation method of the motorized spindle according to claim 7, characterized in that The solution of the motor thermal model includes the following steps: Calculate the copper loss P of the asynchronous motor c : P c = P c1 + P c2 Stator copper loss \(P\) c1 and rotor copper loss \(P\) c2 are calculated by the following formula: In the formula, I1 is the stator winding current; R1 is the winding resistance; I2 is the rotor winding current; R2 is the rotor bar resistance; Calculate the iron loss P of the asynchronous motor Fe : P Fe = P Fet + P Fej Stator (rotor) yoke iron loss P Fej Calculated by the following formula: p Fej = K a p hej G j where K a is an empirical coefficient. When the capacity is, K a = 1.59, otherwise it is 1.3; G j is the weight of the yoke; p hej is the yoke loss coefficient, and the expression is as follows: where P 10 / 50 is the loss per unit weight of the silicon steel sheet when B j = 1T, f = 50Hz, with the unit of W / kg; B j0 is the magnetic induction intensity of the no-load yoke; Stator (rotor) tooth iron loss P Fet Calculated by the following formula: p Fet = K a p het G t Where G t is the weight of the yoke; K a is an empirical coefficient, taking 1.8 for induction motors; p het is the tooth loss coefficient; where P 10 / 50 is the loss per unit weight of the silicon steel sheet when B j = 1T and f = 50Hz, with the unit of W / kg; B t0 is the magnetic induction intensity of the no-load teeth; Calculate the friction loss P between the rotor and the air w : where f r is the rotational frequency; D ro is the outer diameter of the rotor; L r is the axial length of the rotor; η air is the aerodynamic viscosity; h g is the average air gap length; Calculate the stray loss P of the asynchronous motor s : wherein, the copper bar is taken as 0.005, and P out is the output power; Calculate the total heat generation power P: P = P cu + P Fe + P w + P s 。 10. The multi-physical field coupling simulation method of the electric spindle according to claim 7, characterized in that The solution of the heat transfer model includes the following steps: Calculate the convective heat transfer P of the stator core to the air gap 1(2) : In the formula, is the Stefan-Boltzmann constant; ε thr is the relative emissivity, usually taken as 0.85; T2 is the temperature of the heat outflow surface; T1 is the temperature of the heat absorption surface; S is the surface area of the heat dissipating body; The temperature rise T on the surface of the stator (rotor) core is calculated by the following formula: where P Fe is the stator (rotor) iron loss; O is the stator (rotor) heat dissipation area, and k is the heat transfer coefficient of the unit surface area to the external air; Calculate the convective heat transfer amount P3 from the stator (rotor) core to the air gap: P3 = α th |T1 - T2|S The heat transfer coefficient α from the stator (rotor) core to the air gap th is calculated by the following formula: In the formula, λ is the thermal conductivity of air; at room temperature (20°C), the thermal conductivity of air is 0.0267 W / m-K; at 0°C, the thermal conductivity of air is 0.0251 W / m-K; at 100°C, the thermal conductivity of air is 0.0321 W / m-K; The Nusselt number Nu is calculated by the following formula: where Ta≈Ta m ; The Taylor number Ta is calculated by the following formula: where ρ is the mass density of the fluid; Ω is the angular velocity of the rotor; r m is the average radius of the stator and the rotor; δ is the air gap length; μ is the dynamic viscosity of the fluid; Calculate the convective heat transfer amount P4 from the end of the rotor core to the air: P4 = α th |T1 - T2|S Natural convection heat transfer coefficient α th It is calculated by the following formula: In the formula, T1 is the temperature of the heat outflow surface; T2 is the temperature of the heat absorption surface; D is the outer diameter of the rotor; Calculate the convective heat transfer amount P5 from the stator core to the coolant: The heat transfer coefficient α of the stator core to the coolant th is calculated by the following formula: In the formula, v is the coolant velocity; l is the length of the motor body; Calculate the convective heat transfer rate P from the stator core to the air 6(7) : P 6(7) = α th |T1 - T2|S Natural convection heat transfer coefficient α th It is calculated by the following formula: In the formula, T1 is the temperature of the radiation surface; T2 is the temperature of the absorption surface; D is the outer diameter of the stator; Calculate the amount of heat radiation from the stator core to the air In the formula, is the Stefan-Boltzmann constant; ε thr is the relative emissivity, usually taken as 0.85; T1 is the temperature of the radiation surface; T2 is the temperature of the absorption surface; S is the surface area of the radiator; Calculating the heat dissipation P of the bearing bs : Calculate the remaining heat W: Calculate the surface temperature T of the motor's acting rotating shaft motor : Wherein, O motor is the surface area of the rotating shaft affected by the heat generation temperature of the motor; Calculate the surface temperature T of the rotating shaft under the action of the bearing bearing : Where, O bearing is the surface area of the rotating shaft acted upon by the heat generation temperature of the bearing.

Images