CN106960093A - A kind of transmission system method for numerical simulation for considering gear and bearing Non-linear coupling - Google Patents

Abstract

The present invention relates to a kind of transmission system method for numerical simulation for considering gear and bearing Non-linear coupling, it is characterised in that this method comprises the following steps：1) the body unit FEM model of transmission system is set up；2) the transmission system reduced model comprising non-linear bearing unit and the equivalent engaging element of gear is set up；3) the nonlinear static mechanics of transmission system reduced model is solved and linear equivalence bearing rigidity is calculated；4) the transmission system finite element contact analysis model for including linear equivalence bearing rigidity is set up；5) Gear Contact state and the equilibrium iteration of non-linear bearing rigidity are calculated.The present invention combines finite element contact analysis model and advantage of the reduced model in transmission system numerical simulation comprising non-linear bearing unit and the equivalent engaging element of gear, the coupled relation of two kinds of models is set up using the equivalent meshing parameter of gear and equivalent bearing rigidity, by equilibrium iteration, the transmission system numerical simulation for considering gear and bearing Non-linear coupling is accurately realized.

Description

A kind of transmission system method for numerical simulation for considering gear and bearing Non-linear coupling

Technical field

It is particularly a kind of to consider gear and bearing Non-linear coupling the present invention relates to a kind of transmission system method for numerical simulation Transmission system method for numerical simulation, belong to technical field of mechanical transmission.

Background technology

Gear and bearing are the important components of the machine driven systems such as drive axle, gearbox, the engagement of gear, Nonlinear characteristic is presented with the change of system load operating mode in contact condition between bearing roller and raceway.Meanwhile, gear and There is also Non-linear coupling between bearing：On the one hand, the change of gear engagement states can be influenceed between drive shaft system and bearing Load distribution, so as to cause bearing rigidity to change；On the other hand, the change of bearing rigidity can influence gear in correspondence load Magnitude of misalignment under operating mode, and then influence the engagement of gear.Therefore, when being modeled and analyzing to transmission system, Ying Zhun Really consider the Non-linear coupling influence between gear and bearing.Existing research to transmission system when being modeled analysis, generally Using following two method for numerical simulation：

1) finite element contact computational methods：It can accurately realize that finite element is contacted based on the commercial finite element software such as ABAQUS Calculate, but transmission system generally comprises multiple rolling bearings, there is contact relation between the roller and raceway of each bearing, and Finite element contact calculates higher to grid required precision, is limited by convergence and calculation scale, it is difficult in contact computation model In consider the contact with raceway of Gear Contact, each bearing roller simultaneously, so existing research is generally only contacted to gear Definition, ignores the influence of bearing nonlinear characteristic, or individually to bearing progress contact calculating, so can not real embodiment bearing Non-linear coupling between gear.

2) the transmission system modeling and analysis methods based on non-linear bearing unit：The non-linear bearing unit of analytical form (Harris T A,Kotzalas M N.Essential concepts of bearing technology.5th ed.CRC Press, 2006.) nonlinear stiffness characteristic of bearing can be fast and effeciently simulated, but in order to realize the biography for including multiple bearings The Nonlinear Numerical of dynamic system model is solved, and the meshing relation of gear generally uses the equivalent engaging element mould of space form of springs Intend, and using economics analysis formula (Donley M G, Lim T C, Steyer G C.Dynamic analysis ofautomotive gearing systems.SAE International Congress and Exposition, Detroit,Michigan:1992,920762.) the equivalent meshing parameter of gear is calculated, the influence of gear tooth friction is not accounted for, and Influence of the nonlinear gear system contact performance to transmission system can not accurately be embodied.

Still lacking at present a kind of can accurately consider the transmission system method for numerical simulation of gear and bearing Non-linear coupling.

The content of the invention

Regarding to the issue above, gear and bearing Non-linear coupling can accurately be considered it is an object of the invention to provide a kind of Transmission system method for numerical simulation.

To achieve the above object, the present invention takes following technical scheme：A kind of consideration gear and bearing Non-linear coupling Transmission system method for numerical simulation, comprises the following steps：

1) the body unit FEM model of drivetrain components is set up：

The body unit FEM model of input shaft, output shaft and casing is set up, wherein, pass through between input shaft and output shaft Gearing relationships are connected, and input shaft and output shaft are connected by bearing with casing respectively, and input torque is applied on input shaft, Output torque is applied on output shaft, and above-mentioned connection is included in the body unit FEM model of input shaft, output shaft and casing Boundary node is set up at the position of relation and loading, it is specific as follows：

1. input shaft：Bearing centre node, input torque loading node, the equivalent working pitch point of driving gear；

2. output shaft：Bearing centre node, output torque loading node, the equivalent working pitch point of driven gear；

3. casing：Bearing centre node.

Each boundary node is comprising 6 frees degree, respectively with rigid coupling unit by boundary node and the body of each part Unit FEM model is connected, and the body unit FEM model according to actual conditions in casing applies external constraint.

2) the transmission system reduced model comprising non-linear bearing unit and the equivalent engaging element of gear is set up：

Bearing is simulated using the non-linear bearing unit of analytical form, if bearing is axially its local coordinate system z The stiffness matrix K of direction of principal axis, then non-linear bearing unit_bIt is expressed as formula (1)：

In above formula, F_xAnd F_yFor the radial load suffered by bearing unit along local coordinate system x-axis and y-axis direction；F_zFor bearing list Axial force suffered by member；M_xAnd M_yBearing unit is around the radial direction moment of flexure suffered by local coordinate system x-axis and y-axis direction；δ_xAnd δ_yFor axle Bearing unit is along local coordinate system x-axis and the radial deformation in y-axis direction；δ_zFor the axial deformation of bearing unit；θ_xAnd θ_yFor bearing list Corner of the member around local coordinate system x-axis and y-axis direction is deformed.

In order to realize to comprising non-linear bearing unit actuation system models carry out nonlinear iteration numerical solution, it is necessary to Polycondensation conversion is carried out to the body unit model comprising the great deal of nodes free degree, realized using Guyan condensation methods, only retains border section The point free degree, calculates the polycondensation stiffness matrix of each part：

Wherein, the polycondensation stiffness matrix K of input shaft_iIt is expressed as formula (2)：

In above formula, k_iaaFor the corresponding stiffness matrix of the input shaft border free degree；k_ibbFor input shaft internal degree of freedom correspondence Stiffness matrix；k_iabAnd k_ibaFor the stiffness coupling of the input shaft border free degree and internal degree of freedom.

The polycondensation stiffness matrix K of output shaft_oIt is expressed as formula (3)：

In above formula, k_oaaFor the corresponding stiffness matrix of the output shaft border free degree；k_obbFor output shaft internal degree of freedom correspondence Stiffness matrix；k_oabAnd k_obaFor the stiffness coupling of the output shaft border free degree and internal degree of freedom.

The polycondensation stiffness matrix K of casing_hIt is expressed as formula (4)：

In above formula, k_haaFor the corresponding stiffness matrix of the casing border free degree；k_hbbIt is corresponding just for the box house free degree Spend matrix；k_habAnd k_hbaFor the stiffness coupling of the casing border free degree and internal degree of freedom.

One end connection input shaft or the bearing centre node of output shaft reduced model of bearing unit, the other end are connected therewith The bearing centre node of corresponding casing reduced model.

The engagement of gear pair is simulated using equivalent engaging element, the equivalent mesh stiffness matrix K of gear_mRepresent For formula (5)：

In above formula, N is equivalent engagement forced direction unit vector；k_mFor mesh stiffness coefficient.

One end of equivalent engaging element connects the equivalent working pitch point of driving gear of input shaft reduced model, other end connection The equivalent working pitch point of driven gear of output shaft reduced model.

By the polycondensation stiffness matrix K of above-mentioned input shaft_i, output shaft polycondensation stiffness matrix K_o, casing polycondensation stiffness matrix K_h, non-linear bearing unit stiffness matrix K_b, the equivalent mesh stiffness matrix K of gear_mAccording to annexation group collection, that is, obtain The stiffness matrix K of whole transmission system reduced model is expressed as formula (6)：

3) the nonlinear static mechanics of transmission system reduced model is solved and linear equivalence bearing rigidity is calculated：

Apply input torque T in input torque loading node₁, and it is free to constrain the axial rotation of output torque loading node Degree, to eliminate the rigid body free degree of transmission system, using Newton-Raphson method to including the biography of non-linear bearing unit Dynamic system reduced model carries out nonlinear iteration solution, static balance shape that can be in the hope of transmission system under correspondence load working condition State, while obtaining the stiffness matrix K of the non-linear bearing unit shown in poised state following formula (1)_b, because the stiffness matrix is base In the iterative calculation tangent stiffness tried to achieve of the non-linear bearing unit under correspondence load working condition of analytical form, so can not be straight Connect is used for marine hydrostatic calculation as linear rigidity.

In order to realize the simulation to bearing rigidity, this hair in the transmission system finite element contact analysis model subsequently set up It is bright to propose a kind of linear equivalence Rigidity Matrix of Bearings K_be, for the stiffness matrix of the non-linear bearing unit shown in equivalent expression (1) K_bState in transmission system static balance, linear equivalence Rigidity Matrix of Bearings K_beIt is expressed as formula (7)：

From formula (7), linear equivalence Rigidity Matrix of Bearings K_beIn only include diagonal term, each diagonal term represents bearing The direction institute it is loaded with deformation the ratio between, i.e., bearing is equivalent to the Hookean spring of all directions independence, and meet formula (8)：

From formula (8), linear equivalence Rigidity Matrix of Bearings K_beCorrespondence bearing deformation δ_bThe bearing load F of generation_bWith it is non- The load that linear bearing unit iteration is produced when balancing is identical, i.e., under corresponding load working condition, linear equivalence bearing is firm Spend matrix K_beThe stiffness matrix K of non-linear bearing unit when effect played in actuation system models and balance_bIt is identical.

4) the transmission system finite element contact analysis model for including linear equivalence bearing rigidity is set up：

Transmission system finite element contact analysis model is set up in ABAQUS softwares, including：Input shaft, output shaft, casing, Driving gear and driven gear.Wherein, the finite element contact analysis model of driving gear and driven gear is joined according to design of gears Number is set up, the finite element contact analysis model and step 1 of input shaft, output shaft and casing) the middle body unit finite element mould set up Type is identical, and the same position comprising annexation and loading in input shaft, output shaft, the finite element contact analysis model of casing Vertical boundary node is set up, it is specific as follows：

1. input shaft：Bearing centre node, input torque loading node；

2. output shaft：Bearing centre node, output torque loading node；

3. casing：Bearing centre node.

With step 1) it is identical, respectively with rigid coupling unit by the finite element contact analysis model of boundary node and each part Connect, and the finite element contact analysis model according to actual conditions in casing applies external constraint.

Because including the corresponding wheel tooth model of actual design parameter, the meshing relation of gear in finite element contact analysis model By face define embody, so no longer establishment step 1) in the equivalent working pitch point of gear.By the driving gear of foundation The a certain position of engagement is adjusted to the finite element contact analysis model of driven gear, and is defined between the flank of tooth for having meshing relation Surface-to-surface contact relation and gear tooth friction coefficient.The finite element contact analysis model of driving gear is contacted with the finite element of input shaft Analysis model is connected, and the finite element contact analysis model of driven gear and the finite element contact analysis model of output shaft are connected.

With step 3) in the linear equivalence Rigidity Matrix of Bearings K that tries to achieve_beThe bearing centre node and bridge of coupled drive shafting The bearing centre node of shell, i.e., in ABAQUS softwares, set up Bushing mono- between the bearing centre node of pairing successively Member, and stiffness parameters are defined as K_be。

Loading and constraint and step 2) identical, input torque T₁Input torque loading node is applied to, and constrains output and is turned Square loads the axial rotation free degree of node.Can be tried to achieve using ABAQUS softwares between each pair of gear teeth being meshed equivalent is connect Touch vector f_iWith equivalent action point coordinates r_i, it is expressed as formula (9) and formula (10)：

f_i=[f_xi,f_yi,f_zi]^T (9)

r_i=[x_i, y_i,z_i]^T (10)

In above formula, f_xi,f_yi,f_ziFor force component of the equivalent contact force in global coordinate system between the i-th pair gear teeth；x_i, y_i,z_iFor coordinate components of the equivalent operating point between the i-th pair gear teeth in global coordinate system.

Total gear mesh force vector F_mIt is expressed as formula (11)：

F_m=[F_mx,F_my,F_mz]^T (11)

In above formula,N is the gear teeth that engage between gear pair Logarithm.

The equivalent engagement forced direction unit vector N of gear is expressed as formula (12)：

N=[n_mx,n_my,n_mz]^T (12)

In above formula, n_mx=F_mx/|F_m|；n_my=F_my/|F_m|；n_mz=F_mz/|F_m|。

Coordinate position of the equivalent working pitch point in global coordinate system is expressed as formula (13)：

R=[x_m, y_m,z_m]^T (13)

In above formula, x_m, y_m,z_mFor coordinate components.

Can be in the hope of as the equivalent and equalising torque relation of the power shown in formula (14)~formula (16)：

z_m=(M_my+F_mzx_m)/F_mx (15)

y_m=(M_mx+F_myz_m)/F_mz (16)

In above formula, M_myAnd M_mxRespectively gear mesh force is to global coordinate system X-axis and the torque of Y-axis.

5) Gear Contact state and the equilibrium iteration of non-linear bearing rigidity are calculated：

Because there is Non-linear coupling between Gear Contact state and non-linear bearing rigidity, it is necessary to be balanced iteration meter Calculate, process is as follows：

1. the equivalent engagement point coordinates R of gear tried to achieve with the economics analysis formula using propositions such as Donley⁰Nibbled with equivalent Force action direction unit vector N⁰As initial value, according to step 2) in method set up include non-linear bearing unit and gear The transmission system reduced model of equivalent engaging element, as the initial model of iterative calculation, and based on step 3) in formula (7) meter Calculate the initial model corresponding linear equivalence bearing rigidity under correspondence load working condition

2. according to step 4) in method set up transmission system finite element contact analysis model, it is linear etc. by what is tried to achieve in 1. Imitate bearing rigidityFor in contact analysis model, calculating contact condition of the gear under the bearing rigidity, and based on step 4) In formula (9)~formula (16), calculate the equivalent engagement point coordinates R of new gear¹With equivalent engagement forced direction unit vector N¹；

3. by the equivalent engagement point coordinates R of the gear tried to achieve in 2.¹With equivalent engagement forced direction unit vector N¹For 1. In the equivalent engaging element of gear, be again based on step 3) in formula (7) calculate new linear equivalence bearing rigidity

4. repeat 1.~3., when the equivalent meshing parameter that adjacent iteration twice is tried to achieve meets tolerance, calculate convergence, convergence Condition is expressed as formula (17)：

||N^k-N^k-1||+||R^k-R^k-1||<ε (17)

In above formula, N^k-1And N^kThe equivalent engagement forced direction list of gear that respectively kth -1 time and kth time iteration are tried to achieve Bit vector；R^kAnd R^k-1The equivalent working pitch point coordinate of gear that respectively kth -1 time and kth time iteration are tried to achieve；ε is convergence tolorence, For a less positive number, convergence tolorence is taken as 0.01 during actual calculating.

Said process can accurately try to achieve the poised state that any gear engages orientation lower transmission system, and obtain whole pass The stress and deformation result of dynamic various parts, by changing the orientation of gear, can also simulate complete Meshing Process of Spur Gear, Analyzed so as to transmission system more comprehensively calculate.

The present invention is due to taking above technical scheme, and it has advantages below：1st, the invention provides a kind of non-linear axle The linear equivalence method of rigidity is held, bearing tangent stiffness during transmission system reduced model static balance is equivalent to linear rigidity, And for transmission system finite element contact analysis model, can accurately consider bearing rigidity when gear finite element is contacted and calculated Influence.2nd, the present invention combines finite element contact analysis model and comprising non-linear bearing unit and the equivalent engaging element of gear Advantage of the reduced model in transmission system numerical simulation, two are set up using the equivalent meshing parameter of gear and equivalent bearing rigidity The coupled relation of model is planted, by equilibrium iteration, the transmission system number for considering gear and bearing Non-linear coupling is accurately realized Value simulation.3rd, the present invention can accurately calculate complete power train system and stress of wherein any part under different loads operating mode and Deformation characteristicses, design setting model and analysis for transmission system, which are checked, provides effective means.

Brief description of the drawings

Fig. 1 is the schematic flow sheet of the inventive method；

Fig. 2 is main reducing gear of drive axle actuation system models floor map；

Fig. 3 is driving gear shaft FEM model schematic diagram；

Fig. 4 is differential carrier FEM model schematic diagram；

Fig. 5 is axle housing FEM model schematic diagram；

Fig. 6 is drive axle finite element contact analysis model；

Fig. 7 is the linear equivalence rigidity in the translational degree of freedom direction of bearing 2 with gear corner change curve.

Embodiment

The present invention is described in detail with reference to the accompanying drawings and examples.It should be appreciated, however, that the offer of accompanying drawing is only For a better understanding of the present invention, they should not be interpreted as limitation of the present invention.

By taking main reducing gear of drive axle Hypoid Gear Drives system as shown in Figure 2 as an example, driving gear shaft is input End, differential carrier is output end, and global coordinate system origin O is defined as differential mechanism center, driving gear axis and global coordinate system X Axle is parallel, and driven gear axis is parallel with global coordinate system Y-axis, and driving gear is biased along global coordinate system Z axis, hypoid gear Take turns parameter as shown in table 1.

The hypoid gear parameter of table 1

As shown in figure 1, the transmission system method for numerical simulation for considering gear and bearing Non-linear coupling that the present invention is provided, It comprises the following steps：

1) the body unit FEM model of transmission system is set up：

The body unit FEM model of driving gear shaft, differential carrier and axle housing is set up, it is single respectively as shown in Fig. 3~Fig. 5 Element type is four node tetrahedron elements, and boundary node is set up in the position for having annexation in each part：

1. driving gear shaft model has 5 boundary nodes：3 bearing centre node (b_i1,b_i2,b_i3), 1 input torque Load node (l₁), the equivalent working pitch point (m of 1 driving gear₁)；

2. differential carrier model has 4 boundary nodes：2 bearing centre node (b_i4,b_i5), 1 output torque loading section Point (l₂), the equivalent working pitch point (m of 1 driven gear₂)；

3. axle housing model has 5 boundary nodes：Bearing centre node (b_o1,b_o2,b_o3,b_o4,b_o5)。

The initial equivalent meshing parameter of hypoid gear is first calculated using economics analysis formula, the equivalent engagement of gear is tried to achieve The coordinate position R of node⁰For (- 181.025, -18.145,24.931), equivalent engagement forced direction unit vector N⁰For (- 0.7062, -0.2466, -0.6637), the equivalent working pitch point of driving gear and the equivalent working pitch point of driven gear of initial model Coordinate position be disposed as R⁰.Respectively with rigid coupling unit by the node coupling on boundary node and each part body unit model Close, and constrain the translational degree of freedom of the spring block upper surface node of axle housing, simulate constraint of the actual leaf spring to axle housing.

2) the drive axle reduced model comprising non-linear bearing unit and the equivalent engaging element of gear is set up：

Polycondensation conversion is carried out to the body unit FEM model of each part using Guyan condensation methods, only retains boundary node The free degree, calculates the polycondensation stiffness matrix of each part, wherein, driving gear shaft polycondensation stiffness matrix K_iComprising 30 frees degree, Differential carrier stiffness matrix K_oInclude 24 frees degree, axle housing stiffness matrix K_hInclude 30 frees degree.

Rolling bearing, bearing unit non-linear tangent stiffness matrix point are simulated using the non-linear bearing unit of analytical form Wei not K_b1,K_b2,K_b3,K_b4,K_b5, the bearing centre node of bearing unit one end connection drive shaft system, other end connection axle housing Bearing centre node, i.e., as shown in Figure 2：K is used successively_{b j}Coupling boundary node b_{i j}And b_{o j}, j is bearing number, j=1~5.

Simulate the meshing relation of gear using equivalent engaging element, driving gear is equivalent nibbles for the connection of equivalent engaging element one end Close node, the other end connection equivalent working pitch point of driven gear, i.e., in fig. 2, with the equivalent mesh stiffness matrix K of gear_mCoupling Boundary node m₁And m₂。

The method programmed using MATLAB is by above-mentioned polycondensation stiffness matrix K_i,K_o,K_h, non-linear Rigidity Matrix of Bearings K_b1, K_b2,K_b3,K_b4,K_b5, the equivalent mesh stiffness matrix K of gear_mAccording to annexation group collection, driving bridge system reduced model is obtained Stiffness matrix K, drive axle reduced model has 84 frees degree.

3) the nonlinear static mechanics of transmission system reduced model is solved and linear equivalence bearing rigidity is calculated：

In the input torque loading node l of driving gear shaft₁Apply input torque T₁=1616Nm, and constrain differential mechanism Shell output torque loading node l₂The axial rotation free degree, using Newton-Raphson method carry out nonlinear iteration ask Solution, 7 calculating convergences of iteration, takes 10 seconds, tries to achieve static balance state of the driving bridge system under correspondence load working condition, and try to achieve Linear equivalence Rigidity Matrix of Bearings K_be1,K_be2,K_be3,K_be4,K_be5, as shown in 2~table of table 6：

The linear equivalence Rigidity Matrix of Bearings of the bearing 1 of table 2

Kbe1	δX/mm	δY/mm	δZ/mm	θX/rad	θY/rad	θZ/rad
							FX/N	0	0	0	0	0	0
FY/N	0	4.12E+05	0	0	0	0
							FZ/N	0	0	4.05E+05	0	0	0
MX/N·mm	0	0	0	0	0	0
							MY/N·mm	0	0	0	0	1.31E+07	0
MZ/N·mm	0	0	0	0	0	5.80E+06

The linear equivalence Rigidity Matrix of Bearings of the bearing 2 of table 3

Kbe2	δX/mm	δY/mm	δZ/mm	θX/rad	θY/rad	θZ/rad
							FX/N	6.50E+05	0	0	0	0	0
FY/N	0	1.03E+06	0	0	0	0
							FZ/N	0	0	4.77E+06	0	0	0
MX/N·mm	0	0	0	0	0	0
							MY/N·mm	0	0	0	0	3.20E+09	0
MZ/N·mm	0	0	0	0	0	-1.34E+09

The linear equivalence Rigidity Matrix of Bearings of the bearing 3 of table 4

Kbe3	δX/mm	δY/mm	δZ/mm	θX/rad	θY/rad	θZ/rad
							FX/N	-2.08E+04	0	0	0	0	0
FY/N	0	5.58E+04	0	0	0	0
							FZ/N	0	0	7.79E+04	0	0	0
MX/N·mm	0	0	0	0	0	0
							MY/N·mm	0	0	0	0	8.35E+07	0
MZ/N·mm	0	0	0	0	0	1.81E+08

The linear equivalence Rigidity Matrix of Bearings of the bearing 4 of table 5

The linear equivalence Rigidity Matrix of Bearings of the bearing 5 of table 6

Kbe5	δX/mm	δY/mm	δZ/mm	θX/rad	θY/rad	θZ/rad
							FX/N	4.93E+05	0	0	0	0	0
FY/N	0	-3.76E+04	0	0	0	0
							FZ/N	0	0	4.79E+05	0	0	0
MX/N·mm	0	0	0	-1.87E+09	0	0
							MY/N·mm	0	0	0	0	0	0
MZ/N·mm	0	0	0	0	0	-6.49E+09

4) the gear train assembly finite element contact analysis model for including equivalent bearing rigidity is set up：

Gear train assembly finite element contact analysis model is set up in ABAQUS softwares, including：Driving gear shaft, differential Device shell, axle housing, driving gear, driven gear.The finite element contact analysis model of driving gear shaft, differential carrier and axle housing and step The rapid body unit FEM model 1) set up is identical, and boundary node is set up in the position for equally having annexation in each part, wherein：

1. driving gear shaft model has 4 boundary nodes：3 bearing centre node (b_i1,b_i2,b_i3), 1 input torque Load node (l₁)；

2. differential carrier model has 3 boundary nodes：2 bearing centre node (b_i4,b_i5), 1 output torque loading section Point (l₂)；

3. axle housing model has 5 boundary nodes：Bearing centre node (b_o1,b_o2,b_o3,b_o4,b_o5)。

With step 1) it is identical, respectively with rigid coupling unit by boundary node and each part finite element contact analysis model Node coupling, and the translational degree of freedom of the spring block upper surface node of axle housing is constrained, to simulate actual leaf spring to axle housing Effect of contraction.

Because including the corresponding wheel tooth model of actual design parameter in finite element contact analysis model, step is not resettled It is rapid 1) in the equivalent working pitch point of gear.

The finite element contact analysis model of driving gear and driven gear is set up, cell type is hexahedral element, such as Fig. 6 It is shown, gear pair is adjusted to a certain position of engagement, surface-to-surface contact relation is defined between the flank of tooth for having meshing relation, and it is fixed Artificial tooth face coefficient of friction is 0.15.During actual Modeling Calculation, consider analysis and require, computational accuracy and calculate cost, can be with Only the gear teeth for participating in contacting calculating are modeled, and by the control of flank of tooth size of mesh opening within 1mm, to ensure that contact calculates essence Degree and convergence.Adopted between the finite element contact analysis model of driving gear and the finite element contact analysis model of driving gear shaft It is connected with TIE units, between the finite element contact analysis model of driven gear and the finite element contact analysis model of differential carrier Also it is connected using TIE units.

With step 3) in the linear equivalence bearing rigidity K that tries to achieve_beThe bearing centre node of coupled drive shafting and axle housing Bearing centre node, i.e., in ABAQUS softwares, successively in node b_{i j}And bo_jBetween set up Bushing units, and by rigidity Parameter is defined as K_{be j}, j is bearing number, j=1~5.

Loading and constraint and step 2) identical, input torque T₁Input torque loading node is applied to, and constrains output and is turned Square loads the axial rotation free degree of node.

5) Gear Contact state and the equilibrium iteration of non-linear bearing rigidity are calculated：

Because there is Non-linear coupling between Gear Contact state and non-linear bearing rigidity, it is necessary to be balanced iteration meter Calculate, process is as follows：

1. the equivalent meshing parameter R of gear tried to achieve with use economics analysis formula⁰、N⁰As initial value, set up comprising non-linear The transmission system reduced model of bearing unit and the equivalent engaging element of gear, as the initial model of iterative calculation, and calculating should Model corresponding linear equivalence bearing rigidity under correspondence load working condition

2. transmission system finite element contact analysis model is set up, and the equivalent bearing rigidity that will be tried to achieve in 1.For having In the first contact analysis model of limit, contact condition and new gear equivalent meshing parameter R of the gear under the bearing rigidity are calculated¹、 N¹；

3. the equivalent meshing parameter R of gear that will be tried to achieve in 2.¹、N¹For the equivalent engaging element of gear in 1., calculate new Linear equivalence bearing rigidity

4. repeat 1.~3., when the equivalent meshing parameter that adjacent iteration twice is tried to achieve meets tolerance, calculate convergence.

Convergence tolorence is taken as 0.01 during actual calculating, altogether 4 convergences of iteration, the equivalent engagement ginseng of the gear that each iteration is tried to achieve Number is as shown in table 7 and table 8, and it is the initial results tried to achieve using economics analysis formula that wherein iterations 0 is corresponding, due to not examining Consider gear tooth friction and the influence of true Gear Contact state can not be embodied, economics analysis formula result of calculation contacts meter with finite element Calculate result and there is significant difference.

The equivalent working pitch point coordinate of the gear of table 7

Iterations	X-coordinate/mm	Y coordinate/mm	Z coordinate/mm
				0	-181.025	-18.145	24.931
1	-179.246	-14.407	31.272
				2	-179.479	-14.409	31.22
3	-179.491	-14.407	31.222
				4	-179.492	-14.406	31.225

The equivalent engagement forced direction unit vector of the gear of table 8

Iterations	X-component	Y-component	Z component
				0	-0.7062	-0.2466	-0.6637
1	-0.6451	-0.1685	-0.7452
				2	-0.6477	-0.1624	-0.7443
3	-0.6481	-0.1617	-0.7442
				4	-0.6481	-0.1616	-0.7442

The further driving gear shaft and differential side bearing center boundary node of contrast iteration convergence latter two model Displacement result of calculation as shown in Table 9 and Table 10, coincide, and illustrates that the stress of two kinds of models is consistent with deformation state, tests by result of calculation The correctness that driving bridge system static balance state is solved is demonstrate,proved.

The reduced model boundary node displacement of table 9

Boundary node	X displacements/μm	Y coordinate/μm	Z coordinate/μm
				bi1	-14	4	-44
bi2	-13	-7	-57
				bi3	-14	-12	-66
bi4	26	30	13
				bi5	17	31	16

The finite element contact analysis model boundary modal displacement of table 10

Boundary node	X displacements/μm	Y coordinate/μm	Z coordinate/μm
				bi1	-15	4	-43
bi2	-13	-7	-58
				bi3	-14	-11	-65
bi4	26	30	13
				bi5	17	31	16

Said process can engage the poised state of the system under orientation in the hope of any gear, by the side for changing gear Position, can also try to achieve the static balance state of complete gear wheel engagement process, and obtain the change of the bearing rigidity in Meshing Process of Spur Gear Curve, by taking the bearing 2 in Fig. 2 as an example, try to achieve its translational degree of freedom equivalent stiffness with gear corner change curve as shown in fig. 7, Bearing rigidity is in cyclically-varying trend, and wherein Z-direction rigidity fluctuation is more obvious.

In summary, method proposed by the present invention can accurately realize the power train for considering gear and bearing Non-linear coupling System numerical simulation, can be widely applied to the design analysis of all kinds of machine driven systems such as drive axle, gearbox.

The various embodiments described above are merely to illustrate the present invention, wherein the structure of each part, connected mode and manufacture craft etc. are all It can be varied from, every equivalents carried out on the basis of technical solution of the present invention and improvement should not be excluded Outside protection scope of the present invention.

Claims (6)

1. a kind of transmission system method for numerical simulation for considering gear and bearing Non-linear coupling, it is characterised in that this method bag Include following steps：

1) the body unit FEM model of transmission system is set up；

2) the transmission system reduced model comprising non-linear bearing unit and the equivalent engaging element of gear is set up；

3) the nonlinear static mechanics of transmission system reduced model is solved and linear equivalence bearing rigidity is calculated；

4) the transmission system finite element contact analysis model for including linear equivalence bearing rigidity is set up；

5) Gear Contact state and the equilibrium iteration of non-linear bearing rigidity are calculated.

2. a kind of transmission system method for numerical simulation for considering gear and bearing Non-linear coupling as claimed in claim 1, its It is characterised by, above-mentioned steps 1) be specially：

The body unit FEM model of input shaft, output shaft and casing is set up, wherein, pass through gear between input shaft and output shaft Meshing relation is connected, and input shaft and output shaft are connected by bearing with casing respectively, and input torque is applied on input shaft, output Torque is applied on output shaft, and above-mentioned annexation is included in the body unit FEM model of input shaft, output shaft and casing And set up boundary node at the position of loading：

1. input shaft：Bearing centre node, input torque loading node, the equivalent working pitch point of driving gear；

2. output shaft：Bearing centre node, output torque loading node, the equivalent working pitch point of driven gear；

3. casing：Bearing centre node；

Each boundary node is comprising 6 frees degree, respectively with rigid coupling unit by boundary node and the body unit of each part FEM model is connected, and the body unit FEM model according to actual conditions in casing applies external constraint.

3. a kind of transmission system method for numerical simulation for considering gear and bearing Non-linear coupling as claimed in claim 2, its It is characterised by, above-mentioned steps 2) be specially：

Bearing is simulated using the non-linear bearing unit of analytical form, if bearing is axially its local coordinate system z-axis side To the then stiffness matrix K of non-linear bearing unit_bIt is expressed as formula (1)：

$K_b={\left[\begin{array}{cccccc}\frac{\&part;F_x}{\&part;{\mathrm{\&delta;}}_x}& \frac{\&part;F_x}{\&part;{\mathrm{\&delta;}}_y}& \frac{\&part;F_x}{\&part;{\mathrm{\&delta;}}_z}& \frac{\&part;F_x}{\&part;{\mathrm{\&theta;}}_x}& \frac{\&part;F_x}{\&part;{\mathrm{\&theta;}}_y}& 0\\ \frac{\&part;F_y}{\&part;{\mathrm{\&delta;}}_x}& \frac{\&part;F_y}{\&part;{\mathrm{\&delta;}}_y}& \frac{\&part;F_y}{\&part;{\mathrm{\&delta;}}_z}& \frac{\&part;F_y}{\&part;{\mathrm{\&theta;}}_x}& \frac{\&part;F_y}{\&part;{\mathrm{\&theta;}}_y}& 0\\ \frac{\&part;F_z}{\&part;{\mathrm{\&delta;}}_x}& \frac{\&part;F_z}{\&part;{\mathrm{\&delta;}}_y}& \frac{\&part;F_z}{\&part;{\mathrm{\&delta;}}_z}& \frac{\&part;F_z}{\&part;{\mathrm{\&theta;}}_x}& \frac{\&part;F_z}{\&part;{\mathrm{\&theta;}}_y}& 0\\ \frac{\&part;M_x}{\&part;{\mathrm{\&delta;}}_x}& \frac{\&part;M_x}{\&part;{\mathrm{\&delta;}}_y}& \frac{\&part;M_x}{\&part;{\mathrm{\&delta;}}_z}& \frac{\&part;M_x}{\&part;{\mathrm{\&theta;}}_x}& \frac{\&part;M_x}{\&part;{\mathrm{\&theta;}}_y}& 0\\ \frac{\&part;M_y}{\&part;{\mathrm{\&delta;}}_x}& \frac{\&part;M_y}{\&part;{\mathrm{\&delta;}}_y}& \frac{\&part;M_y}{\&part;{\mathrm{\&delta;}}_z}& \frac{\&part;M_y}{\&part;{\mathrm{\&theta;}}_x}& \frac{\&part;M_y}{\&part;{\mathrm{\&theta;}}_y}& 0\\ 0& 0& 0& 0& 0& 0\end{array}\right]}_{6\&times;6}---\left(1\right)$

In above formula, F_xAnd F_yFor the radial load suffered by bearing unit along local coordinate system x-axis and y-axis direction；F_zFor bearing unit institute The axial force received；M_xAnd M_yBearing unit is around the radial direction moment of flexure suffered by local coordinate system x-axis and y-axis direction；δ_xAnd δ_yFor bearing list Member is along local coordinate system x-axis and the radial deformation in y-axis direction；δ_zFor the axial deformation of bearing unit；θ_xAnd θ_yFor bearing unit around The corner deformation in local coordinate system x-axis and y-axis direction；

Polycondensation conversion is carried out to the body unit model comprising the great deal of nodes free degree using Guyan condensation methods, only retains border section The point free degree, and calculate the polycondensation stiffness matrix of each part：

Wherein, the polycondensation stiffness matrix K of input shaft_iIt is expressed as formula (2)：

$K_i=k_{iaa}-k_{iab}{k}_{ibb}^{-1}{k}_{iba}---\left(2\right)$

In above formula, k_iaaFor the corresponding stiffness matrix of the input shaft border free degree；k_ibbIt is corresponding just for input shaft internal degree of freedom Spend matrix；k_iabAnd k_ibaFor the stiffness coupling of the input shaft border free degree and internal degree of freedom；

The polycondensation stiffness matrix K of output shaft_oIt is expressed as formula (3)：

$K_o=k_{oaa}-k_{oab}{k}_{obb}^{-1}{k}_{oba}---\left(3\right)$

In above formula, k_oaaFor the corresponding stiffness matrix of the output shaft border free degree；k_obbIt is corresponding just for output shaft internal degree of freedom Spend matrix；k_oabAnd k_obaFor the stiffness coupling of the output shaft border free degree and internal degree of freedom；

The polycondensation stiffness matrix K of casing_hIt is expressed as formula (4)：

$K_h=k_{haa}-k_{hab}{h}_{hbb}^{-1}{k}_{hba}---\left(4\right)$

In above formula, k_haaFor the corresponding stiffness matrix of the casing border free degree；k_hbbFor the corresponding rigidity square of the box house free degree Battle array；k_habAnd k_hbaFor the stiffness coupling of the casing border free degree and internal degree of freedom；

One end connection input shaft or the bearing centre node of output shaft reduced model of bearing unit, other end connection are corresponded to therewith Casing reduced model bearing centre node；

The engagement of gear pair is simulated using equivalent engaging element, the equivalent mesh stiffness matrix K of gear_mIt is expressed as formula (5)：

$K_m=k_m{\left[\begin{array}{cc}{N}^{T}N& 0\\ 0& 0\end{array}\right]}_{6\&times;6}---\left(5\right)$

In above formula, N is equivalent engagement forced direction unit vector；k_mFor mesh stiffness coefficient；

One end of equivalent engaging element connects the equivalent working pitch point of driving gear of input shaft reduced model, other end connection output The equivalent working pitch point of driven gear of axle reduced model；

By the polycondensation stiffness matrix K of above-mentioned input shaft_i, output shaft polycondensation stiffness matrix K_o, casing polycondensation stiffness matrix K_h、 The stiffness matrix K of non-linear bearing unit_b, the equivalent mesh stiffness matrix K of gear_mAccording to annexation group collection, that is, obtain complete pass The stiffness matrix K of dynamic system reduced model is expressed as formula (6)：

$K=\left[\begin{array}{ccc}{K}_i+K_m+K_b& -K_m& -K_b\\ -K_m& K_o+K_m+K_b& -K_b\\ -K_b& -K_b& K_h+K_b\end{array}\right].---\left(6\right)$

4. a kind of transmission system method for numerical simulation for considering gear and bearing Non-linear coupling as claimed in claim 3, its It is characterised by, above-mentioned steps 3) be specially：

Apply input torque T in input torque loading node₁, and the axial rotation free degree that output torque loads node is constrained, with The rigid body free degree of transmission system is eliminated, using Newton-Raphson method to including the transmission system of non-linear bearing unit Reduced model carries out nonlinear iteration solution, tries to achieve static balance state of the transmission system under correspondence load working condition, obtains simultaneously The stiffness matrix K of non-linear bearing unit under poised state_b；

Using linear equivalence Rigidity Matrix of Bearings K_beCarry out the stiffness matrix K of equivalent non-linear bearing unit_bIt is quiet flat in transmission system State during weighing apparatus, linear equivalence Rigidity Matrix of Bearings K_beIt is expressed as formula (7)：

$K_{be}={\left[\begin{array}{cccccc}\frac{F_x}{{\mathrm{\&delta;}}_x}& 0& 0& 0& 0& 0\\ 0& \frac{F_y}{{\mathrm{\&delta;}}_y}& 0& 0& 0& 0\\ 0& 0& \frac{F_z}{{\mathrm{\&delta;}}_z}& 0& 0& 0\\ 0& 0& 0& \frac{M_y}{{\mathrm{\&theta;}}_x}& 0& 0\\ 0& 0& 0& 0& \frac{M_y}{{\mathrm{\&theta;}}_y}& 0\\ 0& 0& 0& 0& 0& 0\end{array}\right]}_{6\&times;6}---\left(7\right)$

From formula (7), linear equivalence Rigidity Matrix of Bearings K_beIn only include diagonal term, each diagonal term represents bearing at this Direction institute it is loaded with deformation the ratio between, i.e., bearing is equivalent to the Hookean spring of all directions independence, and meet formula (8)：

$K_{be}{\mathrm{\&delta;}}_b=\left[\begin{array}{cccccc}\frac{F_x}{{\mathrm{\&delta;}}_x}& 0& 0& 0& 0& 0\\ 0& \frac{F_y}{{\mathrm{\&delta;}}_y}& 0& 0& 0& 0\\ 0& 0& \frac{F_z}{{\mathrm{\&delta;}}_z}& 0& 0& 0\\ 0& 0& 0& \frac{M_x}{{\mathrm{\&theta;}}_x}& 0& 0\\ 0& 0& 0& 0& \frac{M_y}{{\mathrm{\&theta;}}_y}& 0\\ 0& 0& 0& 0& 0& 0\end{array}\right]\left[\begin{array}{c}{\mathrm{\&delta;}}_x\\ {\mathrm{\&delta;}}_y\\ {\mathrm{\&delta;}}_z\\ {\mathrm{\&theta;}}_x\\ {\mathrm{\&theta;}}_y\\ 0\end{array}\right]=\left[\begin{array}{c}{F}_x\\ F_y\\ F_z\\ M_x\\ M_y\\ 0\end{array}\right]=F_b---\left(8\right)$

From formula (8), linear equivalence Rigidity Matrix of Bearings K_beCorrespondence bearing deformation δ_bThe bearing load F of generation_bWith non-linear axle The load that bearing unit iteration is produced when balancing is identical, i.e., under corresponding load working condition, linear equivalence Rigidity Matrix of Bearings K_beThe stiffness matrix K of non-linear bearing unit when effect played in actuation system models and balance_bIt is identical.

5. a kind of transmission system method for numerical simulation for considering gear and bearing Non-linear coupling as claimed in claim 4, its It is characterised by, above-mentioned steps 4) be specially：

Transmission system finite element contact analysis model is set up in ABAQUS softwares, including：Input shaft, output shaft, casing, active Gear and driven gear；Wherein, the finite element contact analysis model of driving gear and driven gear is built according to design of gears parameter It is vertical, finite element contact analysis model and the step 1 of input shaft, output shaft and casing) the middle body unit FEM model phase set up Together, built at the position comprising above-mentioned annexation and loading in input shaft, output shaft, the finite element contact analysis model of casing Vertical boundary node：

1. input shaft：Bearing centre node, input torque loading node；

2. output shaft：Bearing centre node, output torque loading node；

4. casing：Bearing centre node；

Boundary node is connected with the finite element contact analysis model of each part with rigid coupling unit respectively, and according to actual feelings Condition applies external constraint in the finite element contact analysis model of casing；

The finite element contact analysis model of the driving gear of foundation and driven gear is adjusted to a certain position of engagement, and nibbled Surface-to-surface contact relation and gear tooth friction coefficient are defined between the flank of tooth of conjunction relation；The finite element contact analysis model of driving gear It is connected with the finite element contact analysis model of input shaft, the finite element contact analysis model and the finite element of output shaft of driven gear Contact analysis model is connected；

With the linear equivalence Rigidity Matrix of Bearings K tried to achieve_beThe bearing centre node and the bearing centre of axle housing of coupled drive shafting Node, i.e., in ABAQUS softwares, set up Bushing units, and rigidity is joined between the bearing centre node of pairing successively Number is defined as K_be；

Input torque T₁Input torque loading node is applied to, and constrains the axial rotation free degree that output torque loads node；Profit The equivalent contact force vector f between each pair of gear teeth being meshed is tried to achieve with ABAQUS softwares_iWith equivalent action point coordinates r_i, respectively It is expressed as formula (9) and formula (10)：

f_i=[f_xi,f_yi,f_zi]^T (9)

r_i=[x_i,y_i,z_i]^T (10)

In above formula, f_xi,f_yi,f_ziFor force component of the equivalent contact force in global coordinate system between the i-th pair gear teeth；x_i,y_i,z_i For coordinate components of the equivalent operating point between the i-th pair gear teeth in global coordinate system；

Total gear mesh force vector F_mIt is expressed as formula (11)：

F_m=[F_mx,F_my,F_mz]^T (11)

In above formula,N is the gear teeth logarithm that engages between gear pair；

The equivalent engagement forced direction unit vector N of gear is expressed as formula (12)：

N=[n_mx,n_my,n_mz]^T (12)

In above formula, n_mx=F_mx/|F_m|；n_my=F_my/|F_m|；n_mz=F_mz/|F_m|；

Coordinate position of the equivalent working pitch point in global coordinate system is expressed as formula (13)：

R=[x_m,y_m,z_m]^T (13)

In above formula, x_m,y_m,z_mFor coordinate components；

Tried to achieve as the equivalent and equalising torque relation of the power shown in formula (14)~formula (16)：

$x_m=\frac{\underset{i=1}{\overset{n}{\&Sigma;}}{x}_i\left|f_i\right|}{\underset{i=1}{\overset{n}{\&Sigma;}}\left|f_i\right|}---\left(14\right)$

z_m=(M_my+F_mzx_m)/F_mx (15)

y_m=(M_mx+F_myz_m)/F_mz (16)

In above formula, M_myAnd M_mxRespectively gear mesh force is to global coordinate system X-axis and the torque of Y-axis.

6. a kind of transmission system method for numerical simulation for considering gear and bearing Non-linear coupling as claimed in claim 5, its It is characterised by, above-mentioned steps 5) be specially：

1. the equivalent engagement point coordinates R of gear tried to achieve with the economics analysis formula using propositions such as Donley⁰With equivalent engagement masterpiece With direction unit vector N⁰As initial value, by step 2) in set up non-linear bearing unit and the equivalent engaging element of gear biography Initial model of the system reduced model as iterative calculation is moved, and the initial model is calculated in correspondence load working condition based on formula (7) Under corresponding linear equivalence bearing rigidity

2. the linear equivalence bearing rigidity that will be tried to achieve in 1.For step 4) the middle transmission system finite element contact analysis set up In model, contact condition of the gear under the bearing rigidity is calculated, and based on formula (9)~formula (16), calculate new gear equivalent Engage point coordinates R¹With equivalent engagement forced direction unit vector N¹；

3. by the equivalent engagement point coordinates R of the gear tried to achieve in 2.¹With equivalent engagement forced direction unit vector N¹Used in 1. The equivalent engaging element of gear, is again based on formula (7) and calculates new linear equivalence bearing rigidity

4. repeat 1.~3., when the equivalent meshing parameter that adjacent iteration twice is tried to achieve meets tolerance, calculate convergence, the condition of convergence It is expressed as formula (17)：

||N^k-N^k-1||+||R^k-R^k-1||<ε (17)

In above formula, N^k-1And N^kThe equivalent engagement forced direction Unit Vector of gear that respectively kth -1 time and kth time iteration are tried to achieve Amount；R^kAnd R^k-1The equivalent working pitch point coordinate of gear that respectively kth -1 time and kth time iteration are tried to achieve；ε is convergence tolorence.