An Evolutionary Analytic Method of Multi-DOF Nonlinear Coupling Dynamic Model for Controllable Close-Chain Linkage Mechanism System

The 2-DOF controllable close-chain linkage mechanism is investigated in this paper. Based on the characteristics of the multi-DOF nonlinear coupling dynamic equation of the system established by the finite element method, an analytic method of multiple-scales Newmark is presented after thinking about the method of perturbation and the method of numerical analysis. Firstly, the firstorder approximate solution of the dynamic responses of the system at the time of t is calculated by the multiple scales method. Then, taken the first-order approximate solution as the initialization of the generalized coordinate of the system, the stable dynamic response of the system is obtained by the implicit Newmark method. The simulation and experimental results are given in the end. The studies indicate that the method of multiple-scales Newmark is correct and practicable to study the dynamic characteristics of such kind of multi-DOF nonlinear coupling system.


Introduction
Multi-DOF controllable linkage mechanism, which can accurately actualize the given trajectory, velocity, and acceleration, will has a wide outlook of application in robots, automatic production lines, and so on 1-5 .The dynamic equation of such kind of linkage mechanism is nonlinear coupling time-variant from the dynamics modeling process of the linkage system.The numerical Newmark method can calculate the dynamic responses, but that method could not analyze internal relation between the dynamic characteristics and scale and electromagnetism parameters.Though the multiple-scales method can analyze the vibratory mechanism, it usually adopts first approximations and quadratic or higher order approximation is very complex 6-10 .An evolutionary analytic method of multiple-scales-Newmark method is firstly presented to study the dynamic characteristics of such kind of multi-DOF coupling system synthetically using the property that the implicit Newmark method can calculate the unconditioned stable solution when the Newmark parameters are specific constants the features of the multiple scales method.

Nonlinear Dynamic Equation of System
The 2-DOF controllable close-chain linkage mechanism is investigated in this paper.The analysis diagram of system is shown in Figure 1.Based on the air-gap field of nonuniform airspace of controllable motors of the linkage mechanism caused by the eccentricity of rotor, the controllable motor element 11 as shown in Figure 2, which defined the transverse vibration and torsional vibration of the controllable motors as its nodal displacement, was established.In the diagram, numbers 1, 2, 3, and 4 denote four nodes of the element.So the transverse vibration and the torsional vibration can be expressed by the generalized coordinate vector u 1 U 1 U 2 U 3 U 4 T , and XY Z is the coordinate system of controllable motor element.
According to the mechatronics analysis dynamics, the air-gap eccentric vibration is shown in Figure 3, where point O is the inner circle geometric center of the motor stator, point O 1 is the outer circle geometric center of the rotor journal, point O 2 is the outer circle center of the journal under the deformation of shaft or bearing, the coordinate of point O 3 x, y is the outer circle geometric center of the rotor, δ is the length of air gap, and e 1 is the air-gap eccentricity.e 01 e 2 11 e 2 12 is the static eccentricity, which is caused by the rotor gravity and mismachining tolerance of motor.ε 01 is the rotational eccentricity, which is caused by centering error between the outer circle center of the journal and the outer circle geometric center of the rotor.The static eccentricity and the rotational eccentricity are considered at the same time in this paper.The static and rotational eccentricities also can be ignored, and only the eccentricity caused by rotor vibration can be considered in engineering application latterly.Then, e 1 x 2 y 2 , x U 1 e 11 ε 01 cos β 1 , y U 2 e 12 ε 01 sin β 1 ,

2.1
where u 1 and u 2 are the components of vibration eccentricity in the x and y direction, and is the rotational angle of rotor with respect to the stator of controllable motor.ω 01 is the synchronous speed of rotation of controllable motor.s 1 is the slide ratio.
According to the theory of electromechanics and electromechanical analysis dynamics, as far as controllable motor is concerned, the voltage and current between the two windings are asymmetric, that is to say, the controllable motor works in elliptic magnetic field.According to this real running state of controllable motor, the kinetic energy and potential energy of controllable motor can be obtained as follows: where W 1 x, t and V 1 x, t the transverse displacements and the angle of elastic torsion of any points in the controllable motor element, respectively see the appendix , and l l 1 l 2 l 3 is the length of controllable motor shaft see the appendix ; m 1 x , including the rotor mass m 0 which is at x l 1 , is the mass distribution function of the controllable motor shaft; J 01 x , including the moment of inertia of rotor J 0 which is at x l 1 , is the moment of inertia distribution function of the controllable motor shaft.R 01 is the inner radius of the controllable motor stator, L 01 is the effective length of the rotor, Λ 01 μ 0 /σ is the even air gap permeance of the controllable motor, μ 0 is the magnetic permeability coefficient of air, σ k μ δ 0 , k μ is saturation, δ 0 is the uniform air-gap size, k μ 1 δ F e /k 1 δ 0 , k 1 is the calculation air-gap coefficient of the even air gap, δ F e is the equivalent air-gap of ferromagnetic materials, ϕ 10 is the phase angle of the positive-sequence current of rotor lagging behind the positive-sequence current of stator, and ϕ 20 is the phase angle of the negative-sequence current of rotor lagging behind the negative-sequence current of stator.F s , F −s , F r , and F −r are the positive-sequence and negative-sequence components of the magnetomotive amplitude of stator and rotor respectively see the appendix .M 1 is the mass 4 × 4 matrix of controllable motor element, K 11 is the stiffness 4 × 4 matrix of controllable motor element in connection with the structural parameters of the rotor, K 12 is the stiffness 4 × 4 matrix of controllable motor element in connection with the electromagnetic parameters of the rotor, and k 01 is the 4-order vector in connection with the eccentric motor see the appendix .e e 11 ε 01 cos β 1 e 12 ε 01 sin β 1 0 0 T is a matrix in connection with the static and rotary eccentricity of the controllable motor.
In general, the links of the 2-DOF controllable linkage mechanism are slim bars, so they are adapted to be simulated using beam element as shown in Figure 4.In dynamic analysis of the beam element, the coupling terms of the elastic motion and the rigid body motion in the Coriolis acceleration and transport acceleration are neglected in studying the absolute acceleration of any point in the beam element.In calculation of strain energy, the shearing deformation energy and yield deformation energy are also omitted.The material of the link is adopted as metal.Therefore, the kinetic energy and potential energy, respectively, are as follows: where ρ is the average mass density of beam element, A x is the cross-section area function, and V 3 x, t and W 3 x, t are, respectively, the longitudinal displacement and the transversal displacement of any point in the beam element see the appendix .E is the modulus of elasticity, J x is the moment of inertia distribution function of the element.M 3 and K 3 are the mass matrix and stiffness matrix of the beam element, respectively, and they are 8 × 8 matrix see the appendix .We regard the links BC and CD as two beam elements, respectively, and regard the crank AB and DE as one element respectively as shown in Figure 1.i 3 , i 4 , i 5 , i 6 , i 7 and i 8 are the serial number of the beam element as shown in Figure 5. U 11 , U 12 , U 15 , U 16 , U 19 , U 20 , U 23 , U 24 , U 27 , and u 28 are elastic displacement, U 9 , U 13 , U 14 , U 17 , U 21 , U 22 , U 25 , U 29 , U 30 , and u 31 are elastic rotational angle, u 10 , u 18 , u 26 , and u 32 are the curvature.The serial numbers of the controllable motor elements are i 1 and i 2 .The global coordinate system is established, and the system elastic displacement vector in the global coordinate system is defined as u.Assume that R i is the transformation matrix of the ith element between the element coordinates and the global coordinates, and B i is the coordinate matrix of the ith element between the local numbering and the system numbering see the appendix .
The linear viscous damper damping is adopted in this paper, according to the second Lagrange equation The electromechanical coupling nonlinear dynamic equation of the system can be obtained as where u is the generalized coordinates array of the system, M and C are, respectively, the n × n mass matrix and damping matrix of the system, K and K e 1 are the n × n stiffness matrixes in connection with the structural parameters of system, K 0 and K e 2 are the n × n stiffness matrixes in connection with the electromagnetic parameters of system, k 0 is the n orders array in connection with the electromagnetic parameters of system, P is the external force array of the system, ür is the rigid acceleration array of system in the global coordinates, e is an array in connection with the eccentricity of rotor, ε is the nonlinear term and also is a small parameter, and G n , K n , G kl and K kl are the n × n matrixes in connection with the structural parameters of the linkage mechanism, n is the number of the elements, and here n 8. γ i are the shape functions see the appendix .In 2.7 , −M ür − K e 1 K e 2 e − k 0 are periodic terms.They can be expanded to the form of Fourier series as follows 12 : where v 1 is the working frequency of rotor of controllable motor 1, v 2 is the working frequency of rotor of controllable motor 2, and v 3 is the smallest common multiple of v 1 and v 2 .v 4 v 1 , v 5 v 2 , v 6 is the rotation frequency of magnetic field of stator of the controllable motor 1, v 7 is the rotation frequency of magnetic field of stator of the controllable motor 2, v 5 is the smallest common multiple of v 6 and v 7 , and v 8 is the smallest common multiple of v 1 , v 2 , v 6 , and v 7 .F fki express the amplitude values, τ fki are the corresponding phase angles, and m is the number of terms of the Fourier expansion formula.

First Approximate Solutions
Substituting 2.10 into 2.7 and assuming that the linear transfer functions are u φη, 3.
where φ is the modal transfer matrix, η is the modal coordinate vector corresponding to φ. Substituting 3.1 into 2.7 , and premultiply the equation by φ T , then 2.7 can be transferred as where α s , δ st , and γ stu are the coefficients of η s , η s η t , and η s η t η u , respectively.ξ is the nth-order vector, and the rth element of ξ is 1 and the other elements are zero.Therefore, the system is affected by multifrequency excitations.The method of multiple time scales is employed to study the nonlinear 3.4 , and the frequent factors of system can be obtained.and some resonance phenomena will take place in the system under certain conditions.
Considering the first two orders modal of vibration, the first approximation solutions of the multi-DOF nonlinear coupling system in the generalized coordinate are also obtained by the method of multiple scales as follows: where a r and θ r can be resolved by the method of Newton-Raphson, ω r is the r-order instantaneous natural frequency mean of the system, and Λ rfk Γ rfk exp iτ fki and . Therefore, the resonance characteristics of the multi-DOF nonlinear coupling system can be analyzed using the method of multiple scales, and the first approximate solutions of the dynamic model of the system also can be obtained by that method.

Newmark Stable Numerical Solution
Because the implicit Newmark integration method can calculate the unconditioned stable solution when the Newmark parameters ς and ϑ are specific constants 13, 14 , the method is adopted to obtain the stable solution of the multi-DOF nonlinear coupling system.Therefore, the stable solutions can be obtained using the implicit Newmark integration method based on the first approximate solutions calculated by the method of multiple scales.
In order to calculate the dynamic response of the system as shown in 2.7 , the period t of motion of the system is divided into S time steps Δt under given trajectory, that is to say, T S • Δt.

4.1
And the first approximate solution of dynamic response of the system is solved by the multiple-scales method firstly as shown in 3.7 .The first approximate solution when t is at a moment is looked as the initial value,and then the dynamic response of system is solved by immediate integration implicit Newmark method.When we calculate the response of the system by implicit Newmark method, the hypothesis about the response of some moment must be taken as where ς and ϑ are the Newmark parameters, namely, the constant of integration and ς ≥ 0.5 and ϑ ≥ 0.25 0.5 ς 2 .
Because the multi-DOF nonlinear coupling dynamic 2.7 is time variant, the equation at the moment of t Δt can be expressed as

4.4
Substituting üt Δt calculated according to 4.3 into 4.4 , the u t Δt expressed by u t , ut and üt can be obtained.The calculation accuracy is related to time stepΔt.

Simulation and Experimental Results
Every link of the linkage mechanism is homogeneous.The width and thickness are 30 mm and 2 mm, respectively.The lengths of each links are, crank L 1 200 mm and L 4 150 mm, coupler L 2 L 3 400 mm, and frame L 5 400 mm as shown in Figure 1.The material of link is aluminium.The density of aluminium ρ 2700 kg/m 3 and the Young's modulus E 70 GPa.The lumped mass of the intersection between the crank and the coupler is m 01 0.142 kg.The lumped mass of the intersection between the two couplers is m 02 0.092 kg.

Mathematical Problems in Engineering
The controllable motor 1 is 90ZYT motor and motor 2 is YS8024 motor.The motors are custom made by motor manufacturer.The parameters of motors are offered by the manufacturer.The parameters of the motors are as follows: i 1 The Parameters of the Control Motor 1.The rated power of the 90ZYT motor, P N 0.75 kW, the rated voltage, U kn 220 V, the rated rotational speed of the motor is 1500 r/min, and the stall torque is 2.0 N • m.The static-geometric eccentricity of the motor, e 01 0.75 mm and the rotational eccentricity, ε 01 0.5 mm.The magnetic permeability coefficient of air is μ 0 4π × 10 −7 H/m, the length of the even air is gap δ 0 0.25 mm, and the saturation is The rated power of the YS8024 motor is P N 0.75 Kw, the rated current is I N 3.48/2.01A, the rated voltage is U kn 220 V, and the rated rotational speed is n N 1440 r/ min.The static-geometric eccentricity is e 02 0.73 mm and the rotational eccentricity is ε 02 0.5 mm.The magnetic permeability coefficient of air is μ 0 4π × 10 −7 H/m the length of the even air gap is δ 0 0.25 mm, and the saturation is k μ 1.2.The number of excitation windings of the motor is W 824 and the coefficient is K w 1.The peak value of field current is I 3.58 A. The number of the magnetic pole-pair of compounded magnetic field is p 2. m 1 3 and m 2 0.5 are the number of phases of the stator and rotor respectively.The reactance of field windings is x m 600 Ω.The reduction value of resistance and equivalent self-induction reactance of rotor, respectively, are r 30 Ω and x 0.02 Ω.The slide ratio is s 0.15.The control voltage is U k 95 V.The mass of the motor rotor is m 0 2.93 kg.The moment of inertia of the motor rotor J 0 0.021 kg • m 2 .The length of the motor shaft l 208 mm l 1 100 mm, l 2 33 mm, and l 3 175 mm , the effective length of the rotor L 01 100 mm and the inner radius of the motor stator is R 01 38 mm.
As shown in Figure 1, The initial angles of the two cranks are 0 • .The calculation initial value of the response is the one when the two cranks move to the initial position of them after the system comes to the stabilized state.Given ς 0.25 and ϑ 0.5, the period of motion of the linkage mechanism is divided into S 100 time steps.
The responses of transversal and longitudinal displacement of the midpoints of the links L 2 , and L 3 can be calculated by the method of multiple-scales Newmark mentioned above, then the dynamic responses on the direction perpendicular to the links' axis can be simulated by the Multiple-Scales-Newmark method mentioned above as shown in Figure 6, and the corresponding experimental curve of dynamic response of the links' midpoint on the same parameter conditions can be obtained through experiment as shown in Figure 6.The dynamic responses of midpoints of the links are measured by the dynamic test system in experiment.Comparing the simulation with experimental figures, one can find that the multiple-scales Newmark method studied in this paper is correct and practicable.

Conclusions
The numerical method can calculate the dynamic responses, but that method could not analyze internal relationship between the dynamic characteristics and the scale and electromagnetism parameters.Though the multiple-scales method can analyze the dynamics mechanism, it usually adopt first approximations and second, or higher-order approximation is very complex.An evolutionary analytic method of multiple-scales Newmark is firstly presented in this paper to study the dynamic characteristics of such kind of multi-DOF nonlinear coupling system synthetically using the property that the implicit Newmark method can calculate the unconditioned stable solution when the parameters of Newmark are specific constants with the features of the multiple scales method.The numerical simulation and experimental results indicate that the multiple-scales Newmark method mentioned in this paper is correct and practicable.The study provides the basis of further study on the dynamic control of such kind of mechanism.

Appendix
According to the finite element method, W 1 x, t and V 1 x, t can be expressed as where ψ i x are the shape functions, Mathematical Problems in Engineering where e 1 x/l 1 , e 2 l 1 l 2 − x /l 2 , e 3 l 1 l 2 l 3 − x /l 3 , l 1 , l 2 , and l 3 are, respectively, the length between point 1 and point 2, point 2 and point 3, point 3 and point 4. U i t is the nodal displacement as shown in Figure 2.
The mass matrix M 1 of the electromotor element is 4 × 4 matrix, and the other components are equal to zeros.The stiffness matrixes of the electromotor element are A.4 According to the finite element method, V 3 x, t and W 3 x, t can be expressed as where x is the coordinate of beam element in local coordinate system, and the shape functions are as follows: where I is the moment of inertia of cross-section.The positive-sequence and negative-sequence components of the magnetomotive amplitude of stator and rotor, respectively, are where F m 0.9 WK w /pI, W is the number of turns of field winding, K w is the coefficient of field winding, p is the number of pole pairs, I is the peak value of field current.a e U k /U kn is the effective signal coefficient, U kn is the rated voltage, and U k is the actual control voltage.m 1 and m 2 are the number of phases of the stator and rotor, respectively, x m is the reactance of field windings, r and x are respectively the reduction value of resistance and equivalent self-induction reactance of rotor windings, and the slip ratio s is the difference between the rotor speed and forward revolving field. 2 − s is the slip ratio between the rotor speed and backward revolving field.The transformation matrix of the ith element R i is as follows

Figure 3 :
Figure 3: Diagram of air-gap eccentric vibration of controllable motor.

Figure 4 :
Figure 4: Beam element model of linkage mechanism.

Figure 5 :
Figure 5: Diagram of finite element analysis of mechanism.

Figure 6 :
Figure 6: Simulation and experiment curve of dynamic response for midpoint of link L 2 and L 3 .

γ 2 x 1 − 14 Mathematical Problems in Engineering γ 3 x L × e − 6e 3 8e 4 − 3e 5 ,
10e 3 15e 4 − 6e5 , e x/L, L is the length of the beam element.The mass matrix and stiffness matrix of the beam element are as follows:
g, l 2, 3, 4, 6, 7, 8; r, s, t, u 1, 2, . . ., n , The number of excitation windings of the motor is W 924, and the coefficient is K w 0.92.The peak value of field current, I 3.58 A. The number of magnetic polepairs of the motor is p 1. m 1 2 and m 2 2 are the number of phases of the stator and rotor, respectively.The reactance of field windings is x m 594.35 Ω.The reduction value of resistance and equivalent self-induction reactance of rotor, respectively, are r 27.24 Ω and x 0.0196 Ω.The slide ratio is s 0.0713.The control voltage is U k 26 V.The mass of the motor rotor is m 0 2.2 kg.The moment of inertia of the motor rotor is J 0 0.018 kg • m 2 .The length of the motor shaft is l 363 mm l 1 130 mm, l 1 53 mm, and l 1 180 mm , the effective length of the rotor L 01 140 mm, and the inner radius of the motor stator is R 01 23 mm.
i 2 The Parameters of the Control Motor 2.