A Methodology for Quantifying the Kinematic Position Errors due to Manufacturing and Assembly Tolerances

A systematic and general methodology for kinematic position errors analysis of multibody systems is investigated throughout this work, taking into account the influence of the manufacturing and assemble tolerances on the performance of planar mechanisms. The generalized Cartesian coordinates are used to mathematically formulate kinematic constraints and equations of motion of the multibody systems. Thus, the systems are defined by a set of generalized coordinates, which represents the instantaneous positions of all bodies, together with a set of generalized dimensional parameters that defines the functional dimensions of the system. These generalized dimensional parameters take into account the tolerances associated with the lengths, fixed angles, diameters and distance between centers, among others. This paper highlights the relation among kinematic constraints, dimensional parameters and output kinematic parameters. Based on the theory of dimensional tolerances, the variation of the geometrical dimensions is regarded as a tolerance grade with an interval associated with each dimension and, consequently, a kinematic amplitude variation for the bodies’ position. The methodology presented is implemented in a computational code developed for kinematic analysis of multibody systems, capable of automatically generating and solving the equations of motion for general multibody systems. Finally, a slider-crank mechanism is used as a numerical example to demonstrate the accuracy of the presented methodology, as well as to discuss the main assumptions and procedures adopted in this work. ©2011 Journal of Mechanical Engineering. All rights reserved.


INTRODUCTION
It is well known that the tolerances play a key role in the modern design process by introducing quality improvements and limiting manufacturing costs.According to standard ANSI Y14.5M-1994, tolerances are used to define the allowable limits of geometric variation that are inherent in manufacturing and assemble processes [1].Thus, the assignment of geometric tolerances is always a trade-off between two distinct situations: (i) a part with tight tolerances is good for assembly, but the cost to manufacture the part is increased; (ii) alternatively, loose tolerance in one part may make the whole assemble infeasible.The tolerance analysis deals with the study of the aggregate behaviors of given individual tolerances [2] and [3].
Over the last few decades a number of research papers have been published on the influence of the tolerance and clearance on the kinematic performance of multibody systems [4] to [14].However, most of these research works lack generality, that is, they are developed for specific mechanisms and situations.Garrett and Hall [4] defined the concept of mobility bands to study the effect of tolerance and clearance in the design of linkages.Dhande and Chakraborty [5] presented a stochastic model for the analysis and synthesis of the four-bar mechanism considering tolerances and clearances.Hummel and Chassapis [6] described an approach to the design and optimization of Cardan joints with manufacturing tolerances.Based on the reliability concept, Shi [7] presented and developed a probabilistic model of mechanical error in spatial mechanisms.Choi et al. [8] presented an analytical approach to tolerance optimization for planar mechanisms with lubricated joints based on mechanical error analysis.Wittwer et al. [9] established the direct linearization method applied to position error in kinematic linkages due to the link-length and angle variation.Fogarasy and Smith [10] presented a complete investigation on the study of the influence of manufacturing tolerances on the kinematic response of mechanisms.However, this approach required that constraint equations to be known and independent.Flores et al. [11] and [12] studied the influence of clearance in joints on the kinematic and dynamic performance of multibody systems.These works are valid for both planar and spatial systems and for dry and lubricated joints; however, they do not account for tolerance effects.Dong and Ye [13] modeled and studied a reheat-stop mechanism considering the effects of tolerance, misalignment and thermal action.A good survey on the research work developed in the field of tolerance analysis of kinematic mechanism is provided by Chase and Parkinson [14].
By and large, there are two main approaches to study the effect of the manufacturing tolerances on the kinematic position errors, namely, the deterministic and probabilistic methods.The deterministic method involves fixed values or constraints that are used to find an exact solution.These methods are used mostly when tolerances are known and the worst position error is to be determined.In contrast, the probabilistic or statistical methods deal with random variables that result in a probabilistic response.The statistical approaches are utilized when dimensions have some type of random distribution and the probability of being within a given tolerance band is to be evaluated.
Chase and Greenwood [15] introduced a statistical model to sum tolerances of characteristics affected by variability with mean shift.This model weights the values of tolerances sum between two extremes cases: (i) the worst case and (ii) the root sum square case.The worst case model assumes that all the component dimensions occur, in each assembly, at their extreme and worst limit simultaneously.The root sum square method considers that the component dimensions occur statistically having a Gaussian distribution.From statistical point of view the worst model is the most pessimistic evaluation of the sum variability, meanwhile the root sum square is the more optimistic evaluation of the sum variability.
In general, the manufacturing cost increases geometrically for uniform incremental tightening of tolerances [16].The cost is also related to the characteristics of manufacturing processes used, and the degree of maturity of the workers.The allocated tolerances should be as large as possible for the sake of economy and ease of manufacture.However, large tolerances usually increase mechanical errors.Thus, designers should allocate tolerances to minimize the manufacturing cost while keeping mechanical errors below a certain specific limit [17].
The main purpose of this research work is to present a general and systematic approach to quantify the kinematic position errors due to manufacturing and assemble tolerances.Based on the worst case the deterministic method is utilized.The kinematic constraints and equations of motion of the multibody systems are formulated under the framework of multibody systems methodologies.The system is defined by a set of generalized coordinates, which represents the instantaneous positions of all bodies, together with a set of generalized dimensional parameters that defines the functional dimensions of the system.The generalized dimensional parameters take into account the tolerances associated with the lengths.The relation between the kinematic constraints, dimensional parameters and output kinematic parameters is demonstrated.Finally, the proposed methodology is applied to an elementary planar multibody system in order to demonstrate its features.

KINEMATIC ANALYSIS
The kinematic analysis is the study of the motion of a system, independently of the forces that produce it.Since in the kinematic analysis the forces are not considered, the motion of the system is specified by driving elements that govern the system motion during the analysis, while the position, velocity and acceleration of the remaining elements are defined by kinematic constraint equations that describe the system topology.It is clear that in the kinematic analysis, the number of driver constraints must be equal to the number of degrees of freedom of the multibody mechanical system.In short, the kinematic analysis is performed by solving a set of equations that result from the kinematic and driver constraints.
When the configuration of a multibody system is described by n Cartesian coordinates, then a set of algebraic kinematic independent holonomic constraints Φ can be written in a compact form as [18] to [20], where q is the vector of generalized coordinates and t is the time variable, usually associated with the driving elements.
The velocities and accelerations of the system elements are evaluated through the velocity and acceleration constraint equations.Thus, the first time derivative with respect to time of Eq. ( 1) provides the velocity constraint equations: where Φ q is the Jacobian matrix of the constraint equations, that is, the matrix of the partial derivates, ∂Φ/∂q,  q is the vector of generalized velocities and υ is the right hand side of velocity equations, which contains the partial derivates of Φ with respect to time, ∂Φ/∂t.It should be noticed that only rheonomic constraints, associated with driver equations, contribute with non-zero entries to the vector υ.Furthermore, it is assumed that this vector does not present any dependency on the vector of coordinates.
A second differentiation of Eq. ( 1) with respect to time leads to the acceleration constraint equations, obtained as: where  q is the acceleration vector and γ is the right hand side of acceleration equations, i.e., the vector of quadratic velocity terms, which contains the terms that are only function of velocity, position and time.In the case of scleronomic constraints, that is, when Φ is not explicitly dependent on the time, the term Φ t in Eq. ( 2) and the Φ qt and Φ tt terms in Eq. ( 3) vanish.
The constraint equations represented by Eq. ( 1) are, in general, non-linear in terms of q and are usually solved by the Newton-Raphson method.Eqs. ( 2) and ( 3) are linear in terms of  q and  q , respectively, and can be solved by any usual method for linear equations' systems.Thus, the kinematic analysis of a multibody system can be carried out by solving the set of Eqs. ( 1) to ( 3).The necessary steps to perform this analysis are sketched in Fig. 1, and described as, (i) Specify the initial conditions for positions q 0 and initialize the time counter t 0 .(ii) Evaluate the position constraint Eq. ( 1) and solve it for positions, q. (ii) Evaluate the velocity constraint Eq. ( 2) and solve it for velocities,  q .(iv) Evaluate the acceleration constraint Eq. (3) and solve it for accelerations,  q .
Fig. 1.Flowchart of computational procedure for kinematic analysis of a multibody system (v) Increment the time.If the time is smaller than final time, go to step ii), otherwise stop the analysis.

MATHEMATICAL FORMULATION OF THE KINEMATIC POSITION ERRORS DUE TO TOLERANCES
In order to evaluate, in a systemic and general way, the influence of the manufacturing and assemble tolerances on the kinematic position errors, special attention needs to be given to the mathematical formulation of the description of the systems' configuration.Thus, according to the previous section, the equations of constraints can be written as: where q 1 , q 2 , ..., q n represent the generalized vector of coordinates that define the kinematic system's configuration at any instant, and d 1 , d 2 , ..., d m are the generalized vectors of the dimensional parameters defining the functional dimensions of the system.It should be noted that Eq. ( 4) represents the kinematic system's constraints, which can easily be written using, for instance, Cartesian coordinates.Furthermore, the number of generalized coordinates, n, and the number of generalized dimensional parameters, m, must be adequately selected bearing in mind the correct system's description and system's degrees of freedom.The kinematic analysis of any multibody system implies the resolution of Eq. ( 4) for q 1 , q 2 , ..., q n and their derivatives, according to what was presented in the previous section.In this process, it is assumed that vectors d 1 , d 2 , ..., d m are constants, meaning that there is no variation of the dimensional parameters and, consequently, that they do not affect the global system's performance.However, it is well known that this is not the case in practical engineering design and manufacturing processes.
Since, in general, multibody systems are conducted by driving elements, excluding these elements and considering that the kinematic constraints are independent, the Jacobian matrix can be written as follows: where indices k and l represent, respectively, the n-dr kinematic constraints and n Cartesian coordinates.The number of driving elements is represented by variable dr.Considering all coordinates and dimensional parameters as global system's variables, the variation of the constraint equation is expressed as: In a compact form, Eq. ( 6) is be written as: , (k = 1, …, n-dr) (7) in which, δq i is the variation of the generalized system's coordinates and δd j is the variation of the dimensional parameters.This last term represents the manufacturing tolerances of the corresponding functional dimensions, such as the lengths of the multibody system parts.In a matrix form, Eq. ( 7) is expressed as: where Φ Φ q i is the Jacobian matrix and Φ Φ d i represents the derivative of the constraint equations in relation to the dimensional parameters.

Fig. 2. Double pendulum
Since in the kinematic analysis of multibody systems the Jacobian matrix is known, as illustrated in the previous section, specifying the manufacturing tolerances, δd j , only the matrix Φ Φ d i needs to be evaluated in order to obtain the kinematic position errors of all system's bodies, by solving Eq. ( 8).It should be noted, that with this approach very low computational effort is added to the standard kinematic analysis procedure.Moreover, it should be highlighted that Eq. ( 8) represents a linear system of equations that can easily be solved by employing any numerical method, such as the LU factorization procedure, available in the thematic literature [21] and [22].
With the intent to illustrate the application of this approach, a double pendulum system demonstrative example follows.Fig. 2 shows the double pendulum, where the body numbers, local and global coordinate systems are illustrated.Since this simple multibody system has two ideal revolute joints, the corresponding constraint equations expressed in Cartesian coordinates are written as: where x 2 , y 2 , θ 2 , x 3 , y 3 and θ 3 are the global system coordinates and r 2 and r 3 are the selected dimensional tolerance parameters.The input parameters corresponding to the driving elements are the angles θ 2 and θ 3 , being x 2 , y 2 , x 3 , and y 3 the output parameters.The differentiation of the constraints' Eqs. ( 9) to ( 12) yields the variation of the constraints as follows: Rearranging Eqs. ( 13) to ( 16) results: ---- which represents a linear system of equations that can be solved for the unknowns δx 2 , δy 2 , δx 3 and δy 3 , since the values of the variables θ 2 and θ 3 are specified as inputs and δr 2 and δr 3 represent the tolerance amplitude for the double pendulum arm lengths r 2 and r 3 .
In a compact form Eq. ( 17) is be written as: Again, from Eq. ( 18) great simplicity and generality of the proposed methodology for the study of kinematic position errors due to manufacturing tolerances is evident.In Eq. ( 18) the term Φ q , which corresponds to the partial derivatives of the constraint equations with respect to the dimensional parameters, represents the quantitative influence of the individual variation of the selected dimensional tolerance parameters on the kinematic accuracy of the output parameters.Fig. 3 shows the flowchart of the computational procedure to perform kinematic analysis of a multibody system, in which the kinematic position errors owing to tolerances are evaluated.
The basic and fundamental idea of this study is to present a simple and general approach to kinematic position error analysis of mechanisms.This analysis is simple and general because it is valid for any planar or spatial mechanism.Furthermore, the proposed approach is systematic as it is developed under the framework of multibody system methodologies being included in a standard procedure for kinematic analysis.The main reason for that is based on the fact that the additional computational effort is not significant because almost all the necessary parameters to solve Eq. ( 18) are already known from the standard kinematic analysis.However, this is not necessary to embed this kinematic position error analysis in the general kinematic analysis procedure.Obviously, it is possible to complete the standard kinematic analysis, and then go back at each step of time simulated and evaluate the position error analysis be solving Eq. ( 18).However, this way requires to store the history response, being easier and computationally more efficient to perform the position error analysis at each step of the standard kinematic procedure.The problem of sensitivity analysis of mechanical systems has grasped the attention of several authors such as Arora and Haug [23], Lee et al. [24], Schulz and Brauchli [25], only to mention a few.Furthermore, the problem of optimization and parameter identification of multibody system models using gradient based optimization methods has been developed by Bock and his co-workers [26] and [27].
In order to quantify the kinematic accurate position of a multibody system link, it is first necessary to define the amount of tolerance allowed for each of the dimensional parameters considered.Thus, according to standard ISO 286-1, for common mechanical operating conditions the IT grades are usually in the range IT8 to IT11 [17].For the manufacturing tolerances of the dimensional parameters of a multibody system, the bilateral tolerances specified in ISO 286-1 are commonly used.Hence, in Eq. ( 18), the variation of the functional parameters δd can be regarded as such tolerance fields.Therefore, it is possible to write the following relation: where T represents the total manufacturing tolerance range corresponding to the dimensional parameters.

MANUFACTURING AND ASSEMBLY TOLERANCE ANALYSIS
It is known that tolerances are used to define the allowable limits of geometric variations that are inherent in the manufacturing and assembly processes [1].Broadly, there are two approaches to solve this problem; worst case assemblability and statistical assemblability.In the first case, it is assumed that all the component dimensions Fig. 3. Flowchart of the computational procedure to perform kinematic analysis of a multibody system including the evaluation of the kinematic position errors due to manufacturing tolerances occur, in each assembly, at their extreme and worst limit simultaneously.When this approach is employed, the designer desires to ensure that the components can always be assembled.This means that the probability of having a kinematic error exceeding the specific limits in a particular system is null.On the other hand, the statistical assemblability can be used to take advantage of statistical averaging over of components, allowing for the use of less restrictive tolerances in exchange for admitting the small probability of non-assembly.In general, the standard process is defined at the confidence level corresponding to ±3σ interval [28] and [29].From the statistical point of view the worst case model is the most pessimistic, while the statistical assemblability is the most optimistic case.In practical situations it is expected that they fall between the worst and statistical models.These two approaches are discussed in detail in the following paragraphs.Using Eq. ( 18) as reference, the variation of the generalized system's coordinates can be rewritten as: where s represents the sensitive coefficients given by: s Considering a system with m generalized dimensional parameters and that the maximum tolerance or error is specified, based on the worst scenario assemblability, it is possible to evaluate the maximum error of a general output coordinate as: The average tolerance can be determined by the following expression: In a similar way, the tolerance associated with any system component can be given by: It is worth noting that T maximum is a given quantity, allowing the evaluation of average tolerance values, which ultimately can be used as a guiding reference to specify the manufacturing precision requirements of the system components by selecting the IT grades.
Tolerance is statistical in nature since the output of a random variable is, in general, normally or Gaussian distributed, with the level of confidence three-sigma considered.This means that only two or three cases in a thousand have the probability to be outside the ±3σ range.This confidence level of tolerance becomes an important design parameter to be evaluated and optimized.Thus, using the statistical approach, the root mean square considers that the component dimensions occur statistically having a Gaussian distribution and can be expressed by: In the statistical approach, the average tolerance can be determined by the following expression: In a similar way, the tolerance associated with any system component can be given by: In order to convert the tolerance value to the standard process-tolerance, the γ factor is introduced [3]: Solving Eq. ( 28) for σ yields: Consequently, the admissible tolerance of any component can be expressed by: In Table 1, the γ factor is listed as a function of confidence p value of the tolerance.This table was constructed by integrating the standard normal distribution function f(γ), which can be defined by the error integral in the form:  An elementary slider-crank mechanism is used here to show the influence of manufacturing tolerances on kinematic performance.Fig. 4 shows the configuration of the mechanism, comprising four rigid bodies that represent the crank, connecting rod, slider and ground, three revolute joints and one ideal translational joint.The inertia properties and length characteristics of each body, as well as the associated tolerance according to ISO 286-1 standard are given in Table 2.In the present example the tolerance grade IT 10 was selected establishing the number of generalized dimensional parameters equal to two (m = 2), which are related to the lengths of the crank (r 2 ) and connecting rod (r 3 ).In the kinematic simulation, the crank is the driving element and rotates at a constant angular velocity of 500 rpm clockwise.The initial system configuration corresponds to the top dead point.In the numerical simulation, r 2 and r 3 were selected as dimensional tolerance parameters, being δx 4 and δθ 3 the output parameters.Thus, applying the methodologies presented in the previous sections, Figs. 5 and 6 show, respectively, the maximum absolute errors on the linear slider position and the angular position of the connecting rod, when the worst case approach is considered.These maximum position errors were evaluated at 25 crank angular positions.By observing the obtained global results it can be concluded that the maximum position errors vary during the computational simulation of the slider crank mechanism.

Fig. 6. Maximum angular position error of the connecting rod evaluated over a full crank cycle
In order to clearly show the differences of the worst case and the statistical model, let us consider the mathematical equation that allows the evaluation of the position error of the slider, which can be written as [30]: or in another form: where the manufacturing tolerances on dimensions r 2 and r 3 are represented by T 1 and T 2 , and the sensitive coefficients s 1 and s 2 can be written as follows: From the analysis of Fig. 5, it can be observed that the most critical situation occurs when the crank and the connecting rod are collinear, that is, when θ 2 = 0 or 180º and θ 3 = 180º.In these circumstances, the sensitive coefficients s 1 and s 2 are equal to 1.When the worst case is considered, the maximum tolerance can be calculated using Eq. ( 22) yielding: 1 50 1 70 120 µm.
It is obvious that this value is relatively high from a practical engineers view point.On the other hand, when the statistical model is considered, the maximum tolerance is determined using Eq.(25), that is: 1 50 1 70 86 µm.
This value is still high for practical purposes.Therefore, considering, for example that the maximum admissible tolerance is equal to ±55 μm, then the γ factor can be evaluated using Eq. ( 30): This value corresponds to a level of confidence less than 95%, which is far too high.With the intent to increase the confidence level, the sensitive coefficients associated with the tolerance dimensions should have larger values.This desideratum can be achieved by reducing the tolerance grade IT, for instance from 10 to 9 [17].Consequently, according to ISO 286-1 standard the corresponding tolerance ranges of grade IT 9 for crank and connecting rod lengths are ±31 μm and ±44 μm.Thus, the maximum tolerance is now: µm .
Hence, the γ factor is evaluated as: The corresponding confidence level is over 99.7%, which can be considered to be clearly satisfactory.

CONCLUDING REMARKS
In this paper, a general and systematic methodology for kinematic positional error analysis of multibody systems was investigated, taking into account the influence of the manufacturing and assembly tolerances on the performance of planar mechanisms.In the process, the main aspects for kinematic analysis of multibody systems were revised.Based on the theory of dimensional tolerances, the variation of the geometrical dimensions is regarded as a tolerance grade with an interval associated with each dimension and, consequently, a kinematic amplitude variation for the positions.The presented deterministic method evaluates the relation between variations in the dimensional parameters and variation in the generalized coordinates.The statistical approach based on the confidence level three-sigma was also studied.The methodologies proposed have been exemplified through the application of kinematics to a slidercrank mechanism.The simplicity and generality of the proposed methodology for the study of kinematic position errors due to manufacturing tolerances was thus demonstrated.It should be highlighted that in this paper, only dimensional parameters such as length of links have been addressed.However, other parameters can be integrated in the general methodology presented throughout this work, namely those related to roundness of circular surfaces and clearances in joints.

ACKNOWLEDGMENTS
This work is supported by the Portuguese Foundation for the Science and Technology under the project (PTDC/EME-PME/099764/2008).

Table 2 .
Geometric properties of the slider-crank mechanism