阅读说明：本技术 控制力矩陀螺输出力矩的预示方法 (Indication method for controlling moment gyro output moment ) 是由 韩勤锴 褚福磊 于 2020-06-15 设计创作，主要内容包括：本发明公开了一种控制力矩陀螺输出力矩的预示方法,包括：采用理论分析和数值仿真相结合的方式,考虑控制力矩陀螺CMG(control moment gyros)高速转子微振动的输出力矩模型,明确高速转子微振动源与输出力矩的耦合传递特性,定量揭示关键设计参数对整机输出力矩的影响规律,最后通过试验验证,理论计算干扰频率的误差小于5％,并为整机输出力矩干扰成分的优化提供理论依据。(The invention discloses a method for indicating the output torque of a control moment gyro, which comprises the following steps: by adopting a mode of combining theoretical analysis and numerical simulation, an output torque model of the micro-vibration of a control torque gyro CMG (control moment gyros) high-speed rotor is considered, the coupling transfer characteristic of a micro-vibration source and the output torque of the high-speed rotor is determined, the influence rule of key design parameters on the output torque of the whole machine is quantitatively disclosed, and finally, the error of theoretically calculated interference frequency is less than 5% through test verification, and a theoretical basis is provided for the optimization of the interference component of the output torque of the whole machine.)

1. A method for predicting the output torque of a control moment gyro is characterized by comprising the following steps:

step S1, dividing a control moment gyroscope CMG high-speed rotor system into five subsystems, including a high-speed rotor, a bearing assembly, a flexible support, a vibration isolation device and a low-speed frame, establishing a complete machine system micro-vibration analysis model aiming at the high-speed rotor, the bearing assembly, the flexible support and the vibration isolation device, and solving the complete machine system micro-vibration analysis model through a numerical method to obtain a flywheel rotor system micro-vibration response;

step S2, establishing a low-speed frame coordinate system, calculating gyro moment vectors around three coordinate axes of the low-speed frame caused by rotation and revolution of the high-speed rotor, and processing the gyro moment vectors to obtain a mathematical expression of the CMG output moment of the control moment gyro;

and step S3, presetting the rotating speed of the low-speed frame on the basis of the micro-vibration response of the flywheel rotor system, and substituting the rotating speed into the mathematical expression of the control moment gyro CMG output moment to obtain the theoretical predicted value of the control moment gyro CMG output moment.

2. The method of predicting control moment gyro output moment of claim 1, wherein the low speed frame includes a frame main shaft, a frame bearing and a frame support structure, wherein the frame main shaft is supported by four-point contact rolling ball bearings and connected to the high speed rotor through a mounting bracket in the flexible bearing, and the frame support structure is connected to a satellite body, and the frame main shaft is rotated at a low speed by a motor to output a control moment to control the attitude of the satellite.

3. The method for predicting a control moment gyro output moment as claimed in claim 1, wherein said step S1 further includes:

respectively deducing a high-speed rotor disturbance vibration equation set, a bearing assembly nonlinear supporting force vector, a flexible support motion differential equation set and a vibration isolation device motion differential equation set, and obtaining a complete machine system micro-vibration analysis model through sorting and grouping:

wherein the content of the first and second substances,representing the system degree of freedom vector, M, K₁,C₁G is the mass, stiffness, damping and gyro matrix of the system, respectively, F_e,F_gRepresenting excitation of unbalanced masses and excitation of self-gravity, F_bRepresenting the non-linear disturbance force vector caused by the left and right bearings.

4. The method of predicting an output torque of a control moment gyro of claim 3, wherein the high-speed rotor considers an exciting force f caused by a static unbalance in an actual mounting process of the control moment gyro CMG_u1，f_u2Excitation torque T caused by dynamic unbalance_u1，T_u2In addition, the included angle between the rotation axis and the middle line of the bearing is consideredUnbalanced torque T of the engine_s1，T_s2The calculation formula is as follows:

T_s1＝-(I_dfw-I_pfw)βΩ²cosΩt

T_s2＝(I_dfw-I_pfw)βΩ²sinΩt

wherein, I_dfwIs the diametral moment of inertia of the high-speed rotor, I_pfwIs the polar moment of inertia of the high speed rotor.

5. The predictive method of controlling moment gyro output moment of claim 1 wherein the gyro moment vector is:

G_y＝T_gfJ_fω_f

wherein G is_yIs gyro moment vector, H_fIs the moment of inertia vector, T, of the high-speed rotor_gfRepresenting the angular velocity matrix, ω_fIs a transient angular velocity vector.

Technical Field

The invention relates to the technical field of spacecraft structure dynamics, in particular to a prediction method for controlling the output torque of a moment gyroscope.

Background

A Control Moment Gyro (CMG) is an ideal executing mechanism for realizing rapid attitude maneuver of spacecrafts such as agile satellites. How to optimize the CMG output torque and reduce the disturbance force and the disturbance torque output by the CMG as much as possible are key problems to be solved urgently by a high-performance spacecraft. At present, scholars at home and abroad have developed theoretical and experimental research work aiming at space CMG output torque interference factors, most of research focuses on analyzing the micro-vibration mechanism of a high-speed rotor system, but some problems still exist: (1) the CMG high-speed rotor micro-vibration model is inaccurate; (2) lack of an output torque model that accounts for rotor micro-vibrations; (3) the coupling transmission characteristics of the rotor micro-vibration and the output torque are not clear.

Disclosure of Invention

The present invention is directed to solving, at least to some extent, one of the technical problems in the related art.

Therefore, the invention aims to provide a prediction method for controlling the output torque of the moment gyroscope, which is helpful for clearing the mechanical characteristics of CMG disturbance vibration and provides a theoretical basis for optimizing the disturbance component of the output torque of the whole gyroscope.

The invention aims to provide a prediction method for controlling the output torque of a moment gyroscope.

In order to achieve the above object, an embodiment of the present invention provides a method for predicting an output torque of a control moment gyro, including the following steps: step S1, dividing a control moment gyroscope CMG high-speed rotor system into five subsystems, including a high-speed rotor, a bearing assembly, a flexible support, a vibration isolation device and a low-speed frame, establishing a complete machine system micro-vibration analysis model aiming at the high-speed rotor, the bearing assembly, the flexible support and the vibration isolation device, and solving the complete machine system micro-vibration analysis model through a numerical method to obtain a flywheel rotor system micro-vibration response; step S2, establishing a low-speed frame coordinate system, calculating gyro moment vectors around three coordinate axes of the low-speed frame caused by rotation and revolution of the high-speed rotor, and processing the gyro moment vectors to obtain a mathematical expression of the CMG output moment of the control moment gyro; and step S3, presetting the rotating speed of the low-speed frame on the basis of the micro-vibration response of the flywheel rotor system, and substituting the rotating speed into the mathematical expression of the control moment gyro CMG output moment to obtain the theoretical predicted value of the control moment gyro CMG output moment.

The prediction method for controlling the output torque of the torque gyroscope in the embodiment of the invention adopts a mode of combining theoretical analysis and numerical simulation, firstly provides an output torque model considering the micro-vibration of the CMG high-speed rotor, secondly defines the coupling transfer characteristic of the micro-vibration source of the high-speed rotor and the output torque, quantitatively discloses the influence rule of key design parameters on the output torque of the whole gyroscope, and finally corrects a relevant theoretical model through a test result and provides a theoretical basis for the optimization of the interference component of the output torque of the whole gyroscope.

In addition, the prediction method for controlling the output torque of the moment gyro according to the above embodiment of the present invention may further have the following additional technical features:

further, in an embodiment of the present invention, the low-speed frame includes a frame main shaft, a frame bearing and a frame support structure, wherein the frame main shaft is supported by a four-point contact ball rolling bearing and is connected to the high-speed rotor through a mounting bracket in the flexible bearing, the frame support structure is connected to the satellite body, and the frame main shaft rotates at a low speed under the driving of a motor to output a control torque so as to control the satellite attitude.

Further, in an embodiment of the present invention, the step S1 further includes:

respectively deducing a high-speed rotor disturbance vibration equation set, a bearing assembly nonlinear supporting force vector, a flexible support motion differential equation set and a vibration isolation device motion differential equation set, and obtaining a complete machine system micro-vibration analysis model through sorting and grouping:

wherein the content of the first and second substances,

representing the system degree of freedom vector, M, K₁,C₁G is the mass, stiffness, damping and gyro matrix of the system, respectively, F_e,F_gRepresenting excitation of unbalanced masses and excitation of self-gravity, F_bRepresenting the non-linear disturbance force vector caused by the left and right bearings.

Further, in one embodiment of the present invention, the high speed rotor accounts for the excitation force f caused by static imbalance during actual installation_u1，f_u2Excitation torque T caused by dynamic unbalance_u1，T_u2In addition, the unbalanced torque T caused by the included angle between the rotation axis and the middle line of the bearing is considered_s1，T_s2The calculation formula is as follows:

T_s1＝-(I_dfw-I_pfw)βΩ²cosΩt

T_s2＝(I_dfw-I_pfw)βΩ²sinΩt

wherein, I_dfwIs the diametral moment of inertia of the high-speed rotor, I_pfwIs the polar moment of inertia of the high speed rotor.

Further, in an embodiment of the present invention, the gyro moment vector is:

G_y＝T_gfJ_fω_f

wherein G is_yIs gyro moment vector, J_fIs the moment of inertia vector, T, of the high-speed rotor_gfRepresenting the angular velocity matrix, ω_fIs a transient angular velocity vector.

Additional aspects and advantages of the invention will be set forth in part in the description which follows and, in part, will be obvious from the description, or may be learned by practice of the invention.

Drawings

The foregoing and/or additional aspects and advantages of the present invention will become apparent and readily appreciated from the following description of the embodiments, taken in conjunction with the accompanying drawings of which:

FIG. 1 is a flow chart of a predictive method of controlling moment gyro output torque in accordance with one embodiment of the present invention;

FIG. 2 is a schematic structural diagram of a small control moment gyroscope CMG according to an embodiment of the invention;

FIG. 3 is a schematic view of a whole-machine microvibration analysis model according to an embodiment of the present invention;

FIG. 4 is a schematic view of an exemplary rotor analysis freedom of the present invention;

FIG. 5 is a schematic view of a flywheel rotor system having static and dynamic unbalanced masses according to an embodiment of the present invention;

FIG. 6 is a schematic diagram of an angular contact ball bearing analysis coordinate system according to an embodiment of the present invention;

FIG. 7 is a schematic view showing waviness of inner and outer races and rolling elements of a bearing according to an embodiment of the present invention;

FIG. 8 is a schematic view of the loaded fore-aft raceway groove center of curvature and ball center of the present invention;

FIG. 9 is a schematic diagram of rolling element forces and moments under load according to one embodiment of the present invention;

FIG. 10 is a schematic view of a bearing assembly analytical coordinate system in accordance with an embodiment of the present invention;

FIG. 11 is a flowchart illustrating a method for solving a whole-machine micro-vibration analysis model according to an embodiment of the present invention;

FIG. 12 is a schematic view of a high speed rotor angular velocity and frame coordinate system of one embodiment of the present invention;

FIG. 13 is a schematic diagram of a CMG experimental system.

Detailed Description

Reference will now be made in detail to embodiments of the present invention, examples of which are illustrated in the accompanying drawings, wherein like or similar reference numerals refer to the same or similar elements or elements having the same or similar function throughout. The embodiments described below with reference to the drawings are illustrative and intended to be illustrative of the invention and are not to be construed as limiting the invention.

A predictive method of controlling an output torque of a moment gyro proposed according to an embodiment of the present invention is described below with reference to the accompanying drawings.

Fig. 1 is a flowchart of a predictive method of controlling a moment gyro output moment according to an embodiment of the present invention.

As shown in fig. 1, the method for predicting the output torque of the control moment gyro comprises the following steps:

in step S1, the control moment gyro CMG high-speed rotor system is divided into five subsystems including a high-speed rotor, a bearing assembly, a flexible support, a vibration isolation device and a low-speed frame, a complete machine system micro-vibration analysis model is established for the high-speed rotor, the bearing assembly, the flexible support and the vibration isolation device, and the complete machine system micro-vibration analysis model is solved by a numerical method to obtain a flywheel rotor system micro-vibration response.

Specifically, as shown in fig. 2, the structure of the small CMG system mainly includes a high-speed rotor, a bearing assembly, a deck housing, a mounting bracket, a low-speed frame, and a support structure thereof, wherein the flexible support includes the deck housing and the mounting bracket. A pair of high-precision angular contact ball bearings are installed face to face, the inner ring is fixed on the main shaft, and certain axial preload is applied through the inner loading sleeve and the outer loading sleeve. The rotor is fixed with the bearing outer ring through the bearing mounting shell, and the motor is driven to rotate around the main shaft at a high speed in a working state. One end of the main shaft is connected with the mounting bracket, and the other end is supported on the deck plate shell. The low-speed frame main shaft four-point contact ball bearing supports and is connected with the high-speed rotor through the mounting bracket. The frame supporting structure is connected with the satellite body, and a lower frame main shaft is driven by a motor to rotate at a low speed so as to output a control torque for controlling the satellite attitude. In the CMG operating condition, the high speed rotor speed generally exceeds 5000rpm, while the low speed frame speed does not exceed 10rpm at most. The rotation speeds of the two are very different, the coupling of the two is weaker from the analysis of the vibration angle, and the low-speed frame mainly plays a role in supporting the high-speed rotor. In the embodiment of the invention, the low-speed frame is regarded as a rigid body, and the rotor system is regarded as a complex rigid-flexible coupling system, and the influence of the vibration response on the output torque is analyzed.

As shown in fig. 3, during vibration analysis modeling, a complete machine nonlinear spring-mass centralized parameter analysis model is established for four subsystems, namely a rotor, a bearing assembly, a deck plate shell and a mounting bracket, respectively, respective disturbance vibration analysis equations are derived for the subsystems, and finally, a rotor system vibration analysis model is obtained through equation set. The established system can be described as:

(1) the high-speed rotor is a complex spoke type structure, structural vibration of the high-speed rotor is ignored, the rotor is regarded as concentrated mass, five-degree-of-freedom motion of the space of the rotor is considered, translation in three directions and swing in two directions (torsional vibration around an axis is ignored), and an analysis coordinate system X is established_fw-Y_fw-Z_fw。

(2) The unbalanced disturbing vibration force in the working state is transmitted to the main shaft through the bearing supporting force, and then is transmitted to the satellite body through the mounting bracket and the low-speed frame. While transmitting the unbalanced disturbing force, the internal exciting force can be generated due to the Hertz contact, the elastohydrodynamic lubrication, the waviness and other microscopic factors in the bearing. In the embodiment of the invention, the non-linear vibration disturbing force in five directions of the angular contact bearing is established (namely, the vibration disturbing force in three directions and the vibration disturbing force in two directions except the torque around the axisVibration moment) to describe the disturbance excitation of the two aspects. For the left and right bearings, the disturbance force is f_bLi,f_bRiWherein i is 1,2,3,4, 5.

(3) Respectively using lumped masses m_s,m_fAnd simulating the inertia of the fixed main shaft and the mounting bracket. Like the rotor, in X_fw-Y_fw-Z_fwMass m is concentrated under the coordinate system_sHas translation in three directions and swinging in two directions, and nonlinear disturbance force f between the rotor and the rotor through a bearing_bLi,f_bRiForce coupling is performed. To simplify the analysis, the masses m are concentrated_fOnly the translational freedom degrees along three directions and the concentrated mass m are considered_sCoupled by a spring with the coupling stiffness of k_x1,k_y1,k_z1. Lumped mass n_fThe flexible support is connected with the satellite body, the flexible support can be simplified into three-direction linear springs, and the stiffness of the springs is k_fx,k_fy,k_fzAnd (4) showing. In addition, the upper end of the stationary spindle is usually also supported on the machine housing, with a support stiffness k in both directions_x2,k_y2。

Firstly, a high-speed rotor disturbance equation set is established as follows:

(as shown in fig. 4), the high-speed rotor is a complex spoke-type structure, the structural vibration of the rotor is ignored, the rotor is regarded as a concentrated mass, and an analysis coordinate system X is established_fw-Y_fw-Z_fw. Considering five degrees of freedom motion in rotor space, including translation (u) in three directions_fw,v_fw,w_fw) And around X_fw,Y_fwOscillation of the shaft

According to the balance relation of the forces, a high-speed rotor disturbance equation can be obtained:

in the formula, m_fw,I_dfw,I_pfwMass, diameter moment of inertia and polar moment of inertia, respectively, f_bLi,f_bRiAnd (i ═ 1,2, …,5) respectively represents the nonlinear disturbance vibration force and moment of the left and right angular contact bearings, and a concrete solution is given in the nonlinear disturbance vibration force of the bearings. f. of_u1,f_u2,T_u1,T_u2,T_s1,T_s2The unbalanced mass, the excitation force and the moment caused by the installation deflection.

As shown in fig. 5, the actual spokes and the hub are assembled by machining, and thus mass imbalance inevitably occurs. Two types of unbalanced excitation are considered in the embodiments of the present invention: static imbalance and dynamic imbalance.

For static imbalance, time-varying force excitation is induced, along X_fw,Y_fwComponent of the axis f_u1,f_u2Specifically, it is represented as:

wherein, U_s＝m_sr_msA static out-of-balance mass is indicated,

an initial phase angle representing the amount of static unbalance. To pairIn dynamic unbalance, time-varying moment excitation will be caused. It is likewise directed along X_fw,Y_fwShaft decomposition into T_u1,T_u2The method specifically comprises the following steps:

wherein, U_d＝m_dr_{md f}A representation of the mass of the dynamic unbalance,an initial phase angle representing the amount of dynamic unbalance.

In practical installation situations, a deflection, i.e. an angle β between the axis of rotation and the bearing centerline, is unavoidable. By theoretical analysis, it is believed that skewing likewise causes unbalanced torque T_s1,T_s2The method specifically comprises the following steps:

T_s1＝-(I_dfw-I_pfw)βΩ²cosΩt

T_s2＝(I_dfw-I_pfw)βΩ²sinΩt (4)

next, a bearing assembly nonlinear supporting force vector, namely a bearing nonlinear vibration disturbance force, is obtained, and the embodiment of the invention determines the contact deformation and the contact angle between the rolling element and the inner and outer raceways through load analysis, so as to obtain the bearing nonlinear vibration disturbance force. For angular contact ball bearings, it is often necessary to apply a preload during operation to maintain effective contact between the rolling elements and the raceways and thereby improve bearing slew accuracy and bearing stiffness. Based on the Hertz contact theory and the elastohydrodynamic lubrication theory, the bearing load under the preloading effect is analyzed, and then the nonlinear disturbance vibration force of the bearing is determined, wherein the method specifically comprises the following steps:

(1) analyzing bearing load under preload

As shown in FIG. 6, the outer ring rotates at a constant speed with Ω as the rotation speed, and FIG. 6 shows a bearing global coordinate system (X-Y-Z) and the secondLocal coordinate system (x) of j rolling elements_j-y_j-z_j) Wherein, the center of X-Y-Z is fixed at the axle center of the bearing, X_j-y_j-z_jThen by ω about the Z axis_cjRotation (omega)_cjIndicating the revolution speed of the rolling elements). And (3) considering the torsion of the bearing outer ring around the Z axis, the vectors of the degrees of freedom of the other five directions of the outer ring are as follows: u ═ u, v, w, θ^u,θ^v]^T. Accordingly, the preload vectors in the five directions are: f. of_p＝[f_u,f_v,f_w,m_u,m_v]^TWherein f is_u,f_v,f_wForce loads in three directions, respectively, m_u,m_vThen the torque applied to the X, Y axes.

At a certain time, the position angle of the jth rolling elementWherein phi is₀Initial position angle, N_bThe number of rolling elements. Without loss of generality, let phi₀0. At this moment, the load makes the outer ring move, and then the groove curvature center of the outer ring at the jth rolling element translates along the y and z axes and the swing displacement vector around the x axis is respectively:_j＝[ξ_j,η_j,θ_j]^T. Under the condition of small deformation, the deformation of the steel pipe is reduced,_jthe outer ring displacement u can be obtained through coordinate transformation:

_j＝T_ju+Δw_j(5)

in the formula T_jFor transforming the matrix, can be expressed as

Wherein R is_o＝0.5d_e-(r_go-0.5d_b)cosβ₀The distance from the bearing rotation axis to the curvature center of the outer ring raceway is shown. d_eIndicates the bearing pitch circle diameter, r_goDenotes the outer ring groove curvature radius, d_bDenotes the ball diameter, beta₀Representing the bearing nominal contact angle.

Δ w in formula (5)_jShowing additional displacement vectors due to the waviness of the inner and outer rings and the rolling elements of the bearing. As shown in FIG. 7, p_ij,p_ojRepresents the waviness of the inner and outer races in contact with the rolling elements j in the circumferential direction of the raceways, and q_ij,q_ojThe waviness in the axial direction is indicated. w is a_ij,w_ojRespectively showing the waviness of the rolling element j when contacting the inner raceway and the outer raceway. These waviness are periodic functions of time, which can be described as follows using harmonic functions:

p_ij＝∑A_incos(n_inω_cjt+2π(j-1)/N_b+ψ_in)

p_oj＝∑A_outcos(n_out(Ω-ω_cj)t+2π(j-1)/N_b+ψ_out)

in the formula, A_in,B_in,A_out,B_out,C_n，n_in,n_out,n_baAnd psi_in,

ψ_out,Respectively representing the amplitude, order and phase, omega, of the waviness of the various parts of the bearing_sjIs the rolling body spin speed. Pure rollerUnder the dynamic condition, considering the structural form that the inner ring is fixed and the outer ring rotates, the revolution speed of the rolling elements can be expressed as follows:

the rotation speed is as follows:

the geometrical relationships shown in fig. 6 and 7 are combined to see that:

as shown in FIG. 8, when the load is not applied, the rolling element is only contacted with the inner raceway and the outer raceway, but is not deformed, the center of mass of the rolling element is collinear with the curvature centers of the inner raceway and the outer raceway, and the contact angle is a nominal contact angle beta shown by a dotted line in FIG. 8₀. Once the outer ring bears load, the outer ring is firstly moved, and the rolling body is in contact deformation with the inner and outer raceways, so that the center of mass of the rolling body is moved along with the deformation, and the three are balanced again, as shown by a solid line in fig. 8. Since the inner ring is fixed, the curvature center of the inner ring groove is always kept unchanged, and the displacement of the curvature center of the outer ring groove is (xi) from the formula (5)_j,η_j). The position of the center of the ball is to be quantified in (v)_jy,v_jz) And (4) showing. Due to the loading action, the inner and outer contact angles are unequal, and the contact angle of the outer raceway is reduced to beta_ojWhile the inner raceway contact angle increases to β_ij。

As shown in FIG. 8, before deformation, the center of curvature of the inner ring groove is spaced from the center of the ball by a distance L_ijAfter deformation the distance becomes l_ij. For the distance between the center of curvature of the outer ring groove and the center of the ball, similarly, from L before deformation_ojChange to l_oj. If the curvature radius of the inner and outer raceway grooves is r_gi,r_goConsidering the waviness of the various parts of the bearing, there are:

L_ij＝r_gi+p_ij-h_ci-(0.5d_b+w_ij) (11)

L_oj＝r_go+p_oj-h_co-(0.5d_b+w_oj) (12)

wherein h is_ci,h_coThe thicknesses of the lubricating oil films between the rolling bodies and the inner and outer raceways are respectively shown. According to Hamrock and Dowson's research results, the thickness of the lubricating film can be calculated by the following empirical formula

h_c＝2.69R_xU^0.67G^0.53W^-0.067(1-0.61e^-0.73κ) (13)

In the formula (I), the compound is shown in the specification,representing dimensionless speed, wherein₀Which represents the viscosity of the lubricating oil at atmospheric pressure and at a temperature of 20 ℃,

the equivalent radius of the rolling direction of the rolling body is shown, the outer raceway contact is shown as "+", the inner raceway contact is shown as "-",

is an equivalent rotational speed, and d_ro,d_riThe diameters of the bearing outer raceway and the bearing inner raceway are respectively. G ═ E' c_ηpDenotes a dimensionless elastic modulus, wherein c_ηpIs a viscosity coefficient, E' ═ E/(1-v)²) The effective modulus of elasticity (E is the modulus of elasticity of the material and v is the Poisson's ratio) is shown.

Representing a dimensionless load, wherein Q_jIs the contact load between the rolling bodies and the raceways.

From the geometry shown in fig. 8, it is possible to obtain:

equations (14-17) are geometric equations to be satisfied in the load analysis. In addition to the geometrical conditions, load balancing conditions have to be met, respectively for the rolling elements and the bearing outer ring. The stress situation of the jth rolling element is shown in fig. 9. In the figure F_cjAnd M_gjThe centrifugal force and gyro moment to which the rolling elements rotate are shown

And M_gj＝I_bω_sjω_cjsinα_jWherein m is_b,I_bRespectively rolling body mass and moment of inertia, alpha_jIs the included angle between the rotation axis and the z axis. Under pure rolling condition, can obtain

The contact force between the rolling body and the inner and outer raceways is respectively Q_ij,Q_ojAnd (4) showing. According to the Hertz contact theory, it is possible to obtain:wherein, K_i,K_oAnd_ij,_ojthe contact rigidity coefficient and the contact deformation of the rolling body and the inner and outer raceways are respectively. As can be seen from the foregoing geometric analysis,_ij＝l_ij-L_ijand_oj＝l_oj-L_oj. When in use_ij,_oj>0, χ_ij＝1,χ_oj1 is ═ 1; when in use_ij,_ojWhen the concentration is less than or equal to 0, χ_ij＝0,χ_oj0. For the contact stiffness coefficient, it is determined by the geometric dimension of the contact area, the elasticity of the material, etc., there areWherein

Is the equivalent principal radius of curvature, κ is the ellipticity, ζ and ∈ are the first and second class ellipse integrals, respectively, related to k. According to Harris' special work, the equivalent curvature radius and the ellipticity corresponding to the rolling body and the inner raceway and the rolling body and the outer raceway can be respectively calculated, and then the contact rigidity coefficient K of the inner raceway and the outer raceway can be obtained_i,K_o。

As shown in FIG. 9, in the force analysis, out-of-plane friction is ignored, assuming the gyro moment is exactly balanced by the moment created by the friction in the y-z plane. Lambda [ alpha ]_ij,λ_ojThe coefficient of friction of the rolling elements in contact with the inner and outer raceways is shown as lambda_ij,＝λ_oj1. According to the stress relationship of the rolling body given in fig. 9, the following force balance equation of the rolling body can be obtained:

equations (18) and (19) are a pair of nonlinear algebraic equations. For each rolling element, a similar force balance equation exists.

In the load analysis, the load balance of the outer ring needs to be considered in addition to the stress balance of the rolling body. The magnitude of the load applied to the outer ring by the jth rolling body is Q_ojAndthe direction is opposite to that shown in fig. 9. These loads are converted to the center of curvature of the groove in the outer ring at the jth rolling element, with

The action of all the rolling bodies on the outer ring results in a resultant force:

wherein f is_b＝[f_bu,f_bv,f_bw,m_bu,m_bv]^TThe vector of the resultant force of all rolling bodies acting on the outer ring, namely the nonlinear bearing force vector of the bearing is shown. For the bearing outer ring, the following force balance equation exists:

f_p-f_b＝0 (22)

simultaneous equations (18), (19) and (22) are used to solve the nonlinear algebraic equation system iteratively by the Newton-Raphson method. Obviously, the dimensions of the system of equations are: 2N_b+5. Further, according to the formula (14-17), contact deformation between the rolling elements and the inner and outer raceways can be obtained_ij,_ojAnd contact angle beta_ij,β_oj. At the same time, the deformation displacement vector u of the outer ring under given preload can be obtained_o。

(2) Determining bearing nonlinear disturbance vibration force

As shown in FIG. 10, X_fw-Y_fw-Z_fwCoordinate system representing an analytical model of the bearing assembly, and X_bL-Y_bL-Z_bLAnd X_bR-Y_bR-Z_bRRespectively representing the left and right bearing analysis coordinate systems. The distance between the left and right support bearings and the center of mass is l₁,l₂. From the analysis of the high-speed rotor disturbance equation set, the X_fw-Y_fw-Z_fwThe displacement vector in the coordinate system is:the principal axis displacement vector is:

because the left and right bearing inner rings are fixedly connected with the main shaft, the bearing inner ring displacement can be obtained through the main shaft displacement through coordinate conversion.

As can be seen from FIG. 10, the pre-deformation u caused by the pre-load is taken into account by the coordinate transformation_oAnd the relative deformation vectors of the left bearing outer ring and the right bearing outer ring can be respectively obtained as follows:

wherein the content of the first and second substances,

andrespectively, the displacement vectors of the left and right bearing outer races, and T_LAnd T_RRepresenting the corresponding transformation matrix, which is obtained from the definition of the three coordinate systems in fig. 10:

and after the vibration displacement at the rotating body is converted into the displacement of the bearing outer ring, the nonlinear disturbance vibration force of the bearing is accurately solved through the nonlinear iterative equation.

The following describes in detail the calculation process of the bearing nonlinear support force by taking the left bearing as an example.

The displacement vector of the outer ring of the left bearing is known

The disturbance vibration displacement vector at the jth rolling body of the outer ring of the left bearing can be obtained in the formula (5)_jL. Similarly, the force balance equation satisfied by each rolling element of the left bearing can also be obtained, as shown in equations (18) and (19). The joint type (18) and (19) adopts a Newton-Raphson method to iteratively solve the nonlinear algebraic equation system (the dimension is 2N)_b) And solving to obtain the displacement of each rolling element in a random coordinate system, and further obtaining the load vector applied by the jth rolling element to the bearing outer ring according to the equations (14-16) and (20) as follows: [ Q ]_yjL,Q_zjL,Q_θxjL]^T. Therefore, under the global coordinate system, the resultant force of all the rolling bodies acting on the outer ring (i.e. the disturbance force of the left bearing on the fixed main shaft) can be expressed as:

similarly, the same disturbance force vector can be obtained for the right bearing as follows:

due to the Hertz contact nonlinearity, the lubricating oil film, and the unidirectional contact between the rolling elements and the raceways (only compression and no tension), there is a nonlinear relationship between bearing perturbation force and displacement.

And finally, respectively establishing a motion differential equation set of the flexible support and a motion differential equation set of the vibration isolation device.

Simplifying a fixed spindle to a concentrated mass m_sWith five degrees of freedom, including three-directional translation u_s,v_s,w_sAnd swinging in two directionsAt X_fw-Y_fw-Z_fwIn a coordinate system, the fixed main shaft and the rotor realize the coupling of five-direction freedom degree motion through the nonlinear vibration disturbing force of the bearing. One end of the fixed main shaft passes through the coupling rigidity k_x1,k_y1,k_z1Connected to the mass of the mounting bracket, the other end of which is supported by the deck plate housing and has a stiffness k_x2,k_y2. According to the force balance relationship, the motion differential equation of the fixed spindle can be obtained as follows:

wherein k is_sxz,k_syzRespectively representing the stiffness of the coupling of the transverse and axial movements of the spindle, I_sdRepresenting the principal axis transverse moment of inertia.

Lumped mass m for vibration isolation device_fBesides being coupled with the fixed main shaft through rigidity, the flexible support is connected with the satellite body through the flexible support. The connection rigidity of the flexible support in three directions is k_fx,k_fy,k_fz. Also according to the balance relation of the forces, the motion differential equation can be obtained as

The flexibility and mass inertia of the fixed main shaft, the vibration isolation device and the cabin plate support are respectively considered, respective motion differential equations are established, and the fixed main shaft and the rotor exert nonlinear force coupling action through the rolling bearing.

Further, the perturbation and vibration analysis equations of the subsystems are arranged and assembled, so that a complete machine system micro-vibration analysis model can be obtained:

wherein the content of the first and second substances,

representing the system degree of freedom vector, M, K₁,C₁G is the mass, stiffness, damping and gyro matrix of the system, respectively, F_e,F_gRepresenting excitation of unbalanced masses and excitation of self-gravity, F_bRepresenting the non-linear disturbance force vector caused by the left and right bearings. The matrix and the vector can be obtained by grouping the coefficient matrix and the external excitation vector of each subsystem respectively. The specific expression is as follows

M＝diag([m_fw,m_fw,m_fw,I_dfw,I_dfw,m_s,m_s,m_s,I_sd,I_sd,m_f,m_f,m_f]) (31)

G＝diag([ΩG_fw,0]) (34)

Because the coupling of the disturbance vibration force between the rotor and the main shaft concentrated mass through the nonlinear bearing is considered, the rigidity matrix K₁Submatrix K in_fw＝0，The remainder of the submatrix is represented as follows:

K_f＝diag([(k_x1+k_fx) (k_y1+k_fy) (k_z1+k_fz)]) (38)

the sub-matrices of the gyro matrix G may be represented as:

due to the fact that the nonlinear disturbance vibration force of the rolling bearing is considered, the differential equation set is a nonlinear coupling equation set and needs to be solved through a numerical method. Solving procedure As shown in FIG. 11, the pre-deformation of the outer ring at a given preload is first solved_ojAnd contact angle beta_oj(ii) a Respectively obtaining the nonlinear disturbance vibration force of the left bearing and the right bearing by iteratively solving a nonlinear algebraic equation; substituting the system differential equation into a complete machine system micro-vibration differential equation set formula (30), and calculating the system vibration response at the current moment by adopting a differential equation solver (such as ODE45) in MATLAB. At each calculation time step, nonlinear iterative analysis is required to solve the accurate nonlinear disturbance vibration force of the bearing. If the calculation time does not exceed the set time, repeating the solving process until the calculation time is equal to the set calculation time, and outputting a calculation result to obtain the micro-vibration response of the flywheel rotor system.

In step S2, a low-speed frame coordinate system is established, gyro moment vectors around three coordinate axes of the low-speed frame caused by rotation and revolution of the high-speed rotor are calculated, and the gyro moment vectors are processed to obtain a mathematical expression of the control moment gyro CMG output moment.

Specifically, by analyzing and solving a complete machine system micro-vibration analysis model, vibration displacement and speed response in five freedom directions of the system, namely translation displacement response (u) can be obtained_fw,v_fw,w_fw) Translation speed responseAnd angular displacement response

Angular velocity response

As shown in fig. 12, based on the coordinate transformation theory, the transient angular velocity vector can be obtained as

In the formula

Represents a winding w_fwOf rotational angular velocity, obviously

Represents a coordinate transformation matrix:

the formula (41-42) is substituted into the formula (40) and can be obtained by finishing

From the above equation, it can be seen that when the high speed rotor vibration is not considered, i.e.

Transient angular velocity component ω_fx＝0,ω_fy＝0,ω_fzΩ. Once the vibration induced yaw motion is taken into account, it not only results in a Z-turn_fwThe angular velocity of the shaft varies, and is inX_fw,Y_fwThe shaft induces an additional angular velocity component. This phenomenon will have a significant effect on the CMG output torque.

Further, FIG. 12 also shows a low speed frame coordinate system X_fg-Y_fg-Z_fg. Motor driving frame wound around X_fgRotation at an angular velocity of omega_cg. So that the high-speed rotor coordinate system revolves relatively to the low-speed frame coordinate system, and further outputs the control torque. In practice, the main shaft of the low-speed frame has inevitable installation deflection, so that the revolution speed of the high-speed rotor coordinate system not only winds around X_fgAxes, but there may also be components in the other two coordinate axes. If the mounting deflection is assumed to be X_fg-Z_fgIn the plane and with X_fgAngle of axis being beta_gzThen, the revolution angular velocity component of the high-speed rotor is:

similarly, if the mounting is skewed at X_fg-Y_fgIn the plane and with X_fgAngle of axis being beta_gyThen, the revolution angular velocity component of the high-speed rotor is:

according to the definition of the gyro moment, gyro moment vectors around three coordinate axes of the low-speed frame caused by the rotation and revolution of the high-speed rotor can be obtained through the following relations:

G_y＝T_gfJ_fω_f(46)

in the formula J_f＝diag([I_dfw,I_dfw,I_pfw]) Representing the moment of inertia vector, T, of the rotor_gfRepresenting angular velocity matrix

Considering that the output torque test in the experiment is a fixed direction test, the low-speed frame isContinuously rotating, so that the gyro moment G calculated by the formula (46)_yThe gyro moment G which can be actually tested can be obtained through coordinate transformation_ytAre compared to obtain

Therefore, by substituting the expression (43-45) into the expression (46,48), an expression of the CMG output torque considering the high-speed rotor vibration can be obtained.

For the presence of mounting deflection beta_gzThe case (2) is as follows:

similarly, for the presence of a mounting skew β_gyThe case (2) is as follows:

in step S3, on the basis of the micro-vibration response of the flywheel rotor system, the rotation speed of the low-speed frame is preset and substituted into the mathematical expression of the control moment gyro CMG output moment to obtain the theoretical prediction value of the control moment gyro CMG output moment.

That is, through the above analysis and derivation, the mathematical expressions of the CMG output torque can be obtained, as shown in equations (49-54), respectively. On the basis of obtaining the vibration response of the whole rotor, the theoretical predicted value of the CMG output torque can be obtained through an equation (49-54) by giving the rotating speed of the low-speed frame.

Further, the invention provides a detailed method for testing the dynamic test of the CMG output torque model and comparing, analyzing and verifying the CMG output torque model with the predicted value of a theoretical model, which comprises the following steps:

(1) construction test device

As shown in fig. 13, the vibration testing system includes: the CMG is supported by a three-axis air bearing platform and is placed on a Kistler9253B type force measuring platform; the signal testing system mainly comprises a charge amplifier, an OR35 data acquisition instrument, a testing computer and the like; the coordinate system of the frame is used to define the coordinate system, i.e. X, in the test_fgIn the direction of the axis of rotation of the frame (also perpendicular to the direction of the force-measuring table), Y_fg,Z_fgWhich are respectively the direction perpendicular to the axis of the frame, i.e. the horizontal width direction and the horizontal length direction of the force measuring table. During test, the CMG rotor system rotates around the axis of the frame, and the direction of the rotor rotation vector and Y are assumed initially_fgAnd (4) coinciding, and measuring the dynamic moment around three coordinate axes of the frame in real time by the force measuring platform.

(2) Model parameter analysis method

The parameters of a CMG high-speed system according to an embodiment of the present invention are listed in table 1, and include, support structure, imbalance excitation parameters, and the like. In the embodiment of the invention, the waviness amplitude of the inner and outer raceways along the circumferential direction and the axial direction is set to be 2 microns, the order is 12 orders (namely the number of rolling bodies), and the phase is set to be zero. As can be seen from equations (8) and (9), the characteristic frequencies due to the waviness of the inner and outer races are different. For the inner raceway waviness, the characteristic frequency is:

for the waviness of the outer race, its characteristic frequency is

TABLE 1 CMG high speed rotor System parameters

(3) Model comparison and verification method

When the frame rotating speed to be tested is a certain value in the test, the three shafts output torque and draw a corresponding FFT spectrogram. In order to compare with the experimental result, the theoretical model also needs to calculate a time domain oscillogram of the triaxial output torque and a corresponding FFT result when the corresponding frame rotating speed is calculated, and compares the experimental test with the theoretically predicted high-speed rotor vibration frequency value.

Further, the method for determining the key parameters influencing the CMG output torque according to the established CMG output torque model in the embodiment of the invention specifically comprises the following steps:

the influence of the static unbalance of the rotor is first analyzed. Calculate different U_sIn value, the output torque of the CMG and its frequency spectrum. Further analyzing the influence of dynamic unbalance and analyzing different U_dIn value, the output torque of the CMG and its frequency spectrum.

And analyzing the influence of the bearing waviness, increasing the amplitude of the bearing waviness to a certain value, and calculating a time domain variation curve of the output torque and a corresponding FFT frequency spectrum.

And analyzing the influence of the installation deflection angle, wherein the index for representing the installation deflection is the deflection angle beta. As can be seen from equation (4), the presence of this angle will cause additional unbalanced torque. And when the calculated beta is a certain value, outputting a time domain change curve of the moment and a corresponding FFT frequency spectrum.

And analyzing the influence of the CMG gravity on the output torque to simulate the difference between the actual space environment and the ground simulation environment.

The effect of the preload on the rotor bearings was analyzed to analyze the mechanism of preload contribution to the output torque.

Further, the embodiment of the invention can provide a method for optimizing the disturbance component of the output torque, which comprises two parts, namely disturbance component tracing and disturbance component optimizing methods.

The method is characterized in that the design parameters of the high-speed rotor system, the transmission characteristics between micro-vibration and output torque are determined, the high-frequency interference component of the output torque is traced, and the magnitude of the interference component of the output torque needs to be represented. For a small CMG, where high speed vibration is not considered, and the installation is considered good, and the frame speed is stable, equation (49-51) can be degraded to the ideal output torque, which can be expressed as:

the Y-axis and Z-axis output torques are relatively similar, but are 90 degrees out of phase. At the moment, the root mean square value of the interference component of the Y-axis output torque is used for representing the size of the interference component, namely:

in the formula, RMS (·) represents a root mean square value of the interference component.

According to the interference component index delta G defined by the formula (60), calculating the influence rule of unbalance, bearing waviness, installation error (deflection angle), consideration of gravity or not, bearing axial preload and support flexibility (main shaft rigidity and support rigidity) on the interference component of the output torque. Or aiming at different interference components, combining links such as materials, processing, assembly and the like, an optimization strategy for inhibiting the output torque interference components is provided.

And judging factors with obvious influence according to the interference component tracing and influence analysis results.

In the above examples, the factors judged to have significant influence are dynamic unbalance, bearing waviness, mounting errors and support flexibility, while static unbalance, gravity and bearing preload have little influence. Thus, from a design perspective, the optimization strategy is as follows:

from the perspective of output torque optimization, the static unbalance amount, bearing preload, and gravity parameters may not necessarily be the key control parameters. The amount of static unbalance is still strongly influenced from the micro-vibration angle. Therefore, from the engineering practice point of view, careful consideration should be given;

because the dynamic unbalance, the bearing waviness and the installation error obviously affect the interference components of the output torque and are positively correlated, the sizes of the dynamic unbalance, the bearing waviness and the installation error are strictly controlled in actual design and monitoring, and the influence of the dynamic unbalance, the bearing waviness and the installation error on the interference components is reduced as much as possible;

the influence of the support flexibility on the disturbance component of the output torque is complex and is an important factor for optimizing analysis. In practical design, the output torque disturbance component of the designed CMG is minimized by corresponding optimization design to ensure that the support flexibility is within a reasonable range.

To sum up, the method for predicting the output torque of the control moment gyroscope provided by the embodiment of the invention adopts a mode of combining theoretical analysis and numerical simulation, firstly provides an output torque model considering the micro-vibration of the CMG high-speed rotor, secondly determines the coupling transfer characteristic of the micro-vibration source of the high-speed rotor and the output torque, quantitatively discloses the influence rule of key design parameters on the output torque of the whole gyroscope, and finally corrects a relevant theoretical model through a test result and provides a theoretical basis for the optimization of the interference component of the output torque of the whole gyroscope.

Furthermore, the terms "first", "second" and "first" are used for descriptive purposes only and are not to be construed as indicating or implying relative importance or implicitly indicating the number of technical features indicated. Thus, a feature defined as "first" or "second" may explicitly or implicitly include at least one such feature. In the description of the present invention, "a plurality" means at least two, e.g., two, three, etc., unless specifically limited otherwise.

In the description herein, references to the description of the term "one embodiment," "some embodiments," "an example," "a specific example," or "some examples," etc., mean that a particular feature, structure, material, or characteristic described in connection with the embodiment or example is included in at least one embodiment or example of the invention. In this specification, the schematic representations of the terms used above are not necessarily intended to refer to the same embodiment or example. Furthermore, the particular features, structures, materials, or characteristics described may be combined in any suitable manner in any one or more embodiments or examples. Furthermore, various embodiments or examples and features of different embodiments or examples described in this specification can be combined and combined by one skilled in the art without contradiction.

Although embodiments of the present invention have been shown and described above, it is understood that the above embodiments are exemplary and should not be construed as limiting the present invention, and that variations, modifications, substitutions and alterations can be made to the above embodiments by those of ordinary skill in the art within the scope of the present invention.