Contact Stiffness Study: Modelling and Identification

Computerized fixture design &analysis has become means of providing solutions in production operation improvement. Although fixtures can be designed by using CAD functions, a lack of scientific tool and systematic approach for evaluating the design performance makes them rely on trial-and-errors, which leads to several problems, for instance, over design in functions, which is very common and sometimes depredates the performance (e.g., unnecessary heavy design); the quality of design that cannot be ensured before testing; the long cycle time of fixture design, fabrication, and testing, which may take weeks if not months; a lack of technical evaluation of fixture design in the production planning stage.


Introduction
In machining processes, fixtures are used to accurately position and constrain workpiece relative to the cutting tool.As an important aspect of tooling, fixturing significantly contributes to the quality, cost, and the cycle time of production.Fixturing accuracy and reliability is crucial to the success of machining operations.
Computerized fixture design &analysis has become means of providing solutions in production operation improvement.Although fixtures can be designed by using CAD functions, a lack of scientific tool and systematic approach for evaluating the design performance makes them rely on trial-and-errors, which leads to several problems, for instance, over design in functions, which is very common and sometimes depredates the performance (e.g., unnecessary heavy design); the quality of design that cannot be ensured before testing; the long cycle time of fixture design, fabrication, and testing, which may take weeks if not months; a lack of technical evaluation of fixture design in the production planning stage.
Over past two decades, Computerized Aided Fixture Design (CAFD) has been recognized as an important area and studied from fixture planning, fixture design to fixturing analysis/verification.The fixture planning is to determine the locating datum surfaces and locating/clamping positions on the workpiece surfaces for a totally constrained locating and reliable clamping.The fixture design is to generate a design of fixture structure as an assembly, according to different production requirements such as production volume and machining conditions.And the design verification is to evaluate fixture design performances for satisfying the production requirements, such as completeness of locating, tolerance stack-up, accessibility, fixturing stability, and the easiness of operation.
For many years, fixture planning has been the focus of fixture related academic research with significant progress in both theoretical and practical studies.Most analyses are based on strong assumptions, e.g., frictionless smooth surfaces in contact, rigid fixture body, and single objective function for optimization.Fixture design is a complex problem with considerations of many operational requirements.Four generations of CAFD techniques and systems have been developed: group technology (GT)-based part classification for fixture design and on-screen editing, automated modular fixture design, permanent fixture design with predefined fixture components types, and variation fixture design for part families.The study on a new generation of CAFD just started to consider operational requirements.Both geometric reasoning, knowledge-based as well as case-based reasoning (CBR) techniques have been intensively studied for CAFD.How to make use of the best practice knowledge in fixture design and verify the fixture design quality under different conditions has become a challenge in the fixture design &analysis study.
In fixture design verification, it was proved that when the fixture stiffness and machining force are known as input information, the fixturing stability problem could be completely solved.However most of the studies were focused on the fixtured workpiece model, i.e., how to configure positions of locators and clamps for an accurate and secured fixturing.FEA method has been extensively used to develop fixtured workpiece model (e.g., Fang, 2002;Lee, 1987;Trappey, 1995) with an assumption of rigid or linear elastic fixture stiffness.The models and computational results cannot represent the nonlinear deformation in fixture connections identified in previous experiments.As Beards (1983) pointed out, up to 60% of the deformation and 90% of the damping in a fabricated structure can arise from various connections.The determination of fixture contact stiffness is the key barrier in the analysis of fixture stiffness.The existing work is very preliminary, by either simply applying the Hertzian contact model or considering the effective contact area.
The development of fixture design &analysis tools would enhance both the flexibility and the performance of the workholding systems by providing a more systematic and analytic approach to fixture design.Fixture stationary elements, such as locating pads, buttons, and pins, immediately contact with the workpiece when loading the workpiece.Subsequent clamping (by moveable elements) creates pre-loaded joints between the workpiece and each fixture component.Besides, there may be supporting components and a fixture base in a fixture.In fixture design, a thoughtful, economic fixture-workpiece system maintains uniform maximum joint stiffness throughout machining while also providing the fewest fixture components, open workpiece cutting access, and shortest setup and unloading cycles.Both static and dynamic stiffness in this fixture-workpiece system rely upon the component number, layout and static stiffness of the fixture structure.These affect fixture performance and must be addressed through appropriate design solutions integrating the fixture with other process elements to produce a highly rigid system.This requires a fundamental understanding of the fixture stiffness in order to develop an accurate model of the fixture -workpiece system.

Computer-aided fixture design with predictable fixture stiffness
The research on fixture-workpeice stiffness is a crucial topic in fixture design field.Currently, based on the elastic body assumption, using FEA method to predict the fixture stiffness has been widely accepted.With the consideration on the contact and friction conditions, the validity and accuracy of the methodology was been illustrated by two cases simulation and experimental comparison (Zhu, 1993).
The following is an introduction on the general methodology.
First the stiffness of typical fixture units is studied with considerations of contact friction conditions.The results of the fixture unit stiffness analysis are integrated in fixture design as a database with variation capability driven by parametric representations of fixture units.When a fixture is designed using fixture design &analysis tool, the fixture stiffness at the contact locations (locating and clamping positions) to the workpiece can be estimated and/or designed based on the machining operation constraints (e.g., fixture deformation and dynamic constraints).Fig. 1 shows a diagram of the integrated fixture design system.In order to study the fixture stiffness in a general manner, fixture structure is decomposed into functional units with fixture components and functional surfaces (Rong, 1999).In a fixture unit, all components are connected one to another where only one is in contact directly with fixture base and one or more in contact with the workpiece serving as the locator, clamp, or support.Fig. 2 shows a sketch of the fixture units in a fixture design.When a workpiece was located and clamped in the fixture, the fixture units are subjected to the external loads that pass through the workpiece.If the external load is known and acting on a fixture unit, and the displacement of the fixture unit at the contact position is measured or calculated based on a finite element (FE) model, the fixture unit stiffness can be determined.
The fixture unit stiffness is defined as the force required for a unit deformation of the fixture unit in normal and tangential directions at the contact position with workpiece.The stiffness can be static if the external load is static (such as clamping force), and dynamic if the external load is dynamic (such as machining force).It is the key parameter to analyze the relative performance of different fixture designs and optimize the fixture configuration.
Analysis of fixture unit stiffness may be divided into three categories: analytical, experimental and finite element analysis (FEA).Conventional structural analysis methods may not work well in estimating the fixture unit stiffness.Preliminary experimental study has shown the nature of fixture deformation in T-slot based modular fixtures (Zhu, 1993).An integrated model of a fixture-workpiece system was established for surface quality prediction (Liao, 2001) based on the experiment results in (Zhu, 1993), but combining zhu's experimental work and finite element analysis (FEA).Hurtado used one torsional spring, and two linear springs, one in the normal direction and the other in the tangential direction, to model the stiffness of the workpiece, contact and fixture element.(Hurtado, 2002) FEA method has not been studied for fixture unit stiffness due to the complexity of the contact conditions and the large computation effort for many fixture components involved.

Finite element formulation
Consider a general fixture unit with two components I and J, as shown in Fig. 3 (Zheng, 2005(Zheng, , 2008)).For multi-component fixture units, the model can be expanded.The fixture unit is discretized into finite element models using a standard procedure, except for the contact surfaces, where each nodes on the finite element mesh for the contact surface is modelled by a pair of nodes at the same location belonging to components I and J, respectively, which are connected by a set of contact elements.The basic assumptions include that material is homogenous and linearly elastic, displacements and strains are small in both components I and J, and the frictional force acting on the contact surface follows the Coulomb law of friction.
The total potential energy p of a structural element is expressed as the sum of the internal strain energy U and the potential energy Ω of the nodal force; that is, It is well known that the element strain energy can be expressed as,

 
where   K is the element stiffness matrix; and   q is the element nodal displacement vector.
The potential energy of the nodal force is, where   R is the vector of the nodal force.It includes internal force and external force.
When the two components I and J are in contact, a number of three-dimensional contact elements are in effect on the contact surfaces.It should note that the problem is strongly nonlinear, partially due to the fact that the number of contact elements may vary with the change of contact condition.The original contacting nodes might separate or recontact after separation, based on the deformation condition on the contact surface; also contact stiffness may not constant either.The contact elements are capable of supporting a compressive load in the normal direction and tangential forces in the tangential directions.When the two components are in contact, and the displacements in the tangential directions and normal direction are assumed as independent, the element itself can be treated as three independent contact springs: two having stiffness k t and k in the tangential directions of the contact surface at the contact point and one having stiffness k n in the normal direction.
Usually, there are two methods used to include the contact condition in the energy equation: the Lagrange multiplier and the penalty function methods.In order to understand these methods, a physical model of the contact conditions is presented, shown in Fig. 4. When two contact surfaces of fixture components, i.e., body J and I, are loaded together, they will contact at a few asperities.

Fig. 3. Contact Model of Two Fixture Components
The contact criteria can be written as: Where,  is distance from a contact point i in body I to a contact point j on the body J in the normal direction of contact; f ni is the contact force acting on point i of body I in the normal direction.

325
It shows the kinematic condition of no penetration and the static condition of compressive normal force.To prevent interpenetration, the separation distance  for each contact pair must be greater or equal to zero.If >0, the contact force f ni =0.When =0, the points are in contact and f ni <0.If <0, penetration occurs.In real physics, the actual contact area increases, and contact stiffness is enhanced when the load increases.Therefore, the contact deformation is nonlinear as a function of the preload as shown Figure 4(e).In the Lagrange multiplier method, the function w(, f ni ) represents the constraint, which prevents the penetration between contact pairs.In the penalty function method, an artificial penalty parameter is used to prevent the penetration between contact pairs.
Body I State after loading d) Normal contact deformation curve e)

Fig. 4. Physical Model of the Contact Conditions
In the penalty function method, the contact condition is represented by the constraint equation, Where {t} is the constraint equation,   C K is the contact element stiffness matrix,   Q is the contact force vector of the active contact node pairs.When     0 t  , it means that the constraints are satisfied.So the constraint equation Eq. 4 becomes The total potential energy Π p in Eq. 1 can be augmented by a penalty function    where    is a diagonal matrix of penalty value i  .The total potential energy in the penalty function method becomes The minimization of Π pP with respect to   where is the penalty matrix.
On the other hand, in the Lagrange multiplier method, the contact constraint equation can be written as: where the components of the row vector i  (i=1, 2, …, N), are often defined as Lagrange multipliers i  .
Adding Eq. 8 to the potential energy in Eq. 1, we have the total energy in the Lagrange multiplier method, The minimization of Π pL with respect to   q and    requires that In a matrix form, Eqs. 10 and 11 can be expressed as, While the constraints in Eq. 8 can be satisfied, the Lagrange multiplier method has disadvantages.Because the stiffness matrix in Eq. 12 may contain a zero component in its diagonal, there is no guarantee of the absence of the saddle point.In this situation, the computational stability problem may occur.In order to overcome that difficulty, a perturbed Lagrange multiplier method was introduced (Aliabadi, 1993).
where  is an arbitrary positive number.At the limit  goes to  , the perturbed solutions converge to the original solutions.The introduction of  will maintain a small force across and along the interface.This will not only maintain stability but also avoid the stiffness matrix being singular, due to rigid body motion.Similarly, the minimization of p pL  with respect to   q and    results in the following matrix, Eq. 14 can be expressed as: Substitute Eq. 16 into Eq.15, For simplicity, let all i in [] of penalty function equal to  , i.e.  i =  .Thus, the perturbed Lagrange multiplier is equivalent to the penalty function method.
In the Lagrange multiplier method, both displacement and contact force are regarded as independent variables; thus, the constraint (contact) conditions can be satisfied and the contact force can be calculated.It has disadvantages.The stiffness matrix contains zero www.intechopen.com Finite Element Analysis -From Biomedical Applications to Industrial Developments 328 components in its diagonal, and the Lagrange multiplier terms must be treated as additional variables.This leads to the construction of an augmented stiffness matrix, the order of which may significantly exceed the size of the original problem in the absence of constraint equations (Aliabadi, 1993).In comparison with the Lagrange multipliers method, the implementation of the penalty function method is relatively simple and does not require additional independent variables.It is often adopted in the practical analysis because of its simple implementation.

Contact conditions
Based on an iterative scheme (Mazurkiewicz, 1983), the contact conditions in FEA model are classified into the following three cases: 1. Open condition: gap remains open; 2. Stick condition: gap remains closed, and no sliding motion occurs in the tangential directions; and 3. Sliding condition: gap remains closed, and the sliding occurs in the tangential directions.
Let f ji and u ji be the contact nodal load vector and the nodal displacement, respectively, which are defined in the local coordinate system, where the subscript j indicates the component number ( j = I or J), and i indicates the coordinate (i = n, t, τ), as shown in Fig. 5.By equilibrium of the contact element, 0 where x is the node number of body I or J.
The displacement and force must satisfy the equilibrium equations in the three contact conditions (note that {n, t, τ} is the local coordinate system).

Open condition
When the normal nodal force F n is positive (tension), the contact is broken, and no force is transmitted.The displacement change in the normal and tangential directions, denoted respectively by

 
,, where u Jn and u In are the current displacements of node J and node I in a normal direction, respectively.For each structural contact element, stiffness and forces are updated, based upon current displacement values, in order to predict new displacements and contact forces.

Stick condition
The force in the tangential direction   S F , which is the composition of the nodal force in t and  directions (F t and F  ), is defined only when , where  is the Coulomb friction coefficient, there is no slide-motion in the interface, and the contact element responds like a spring.The stick condition exists if where k t and k  are the tangential contact stiffness in t and  directions, respectively.In the analysis of fixture unite stiffness, set t kk   .

Sliding condition
Slide-motion will occur when the absolute value of S F is more than || n F


. The slidemotion may occur in both the element t and  directions.That is, if where   mean the maximum friction force in t and  directions.

Solution procedure
The model presented in the previous section can be implemented to determine the fixture unit stiffness in clamping and machining.Because the model involves high nonlinearity, the Newton-Raphson (N-R) approach is used to solve the problem.Considering the full Newton-Raphson iteration it is recognized that in general the major computational cost per iteration lies in the calculation and factorization of the stiffness matrix.Since these calculations can be quite expensive when large-order systems are considered, the modified Newton-Raphson algorithm is used in this research (Bathe, 1996).
Given the applied load R and the corresponding displacement u, the applied load is divided into a series of load increments.At each load step, the contact stiffness and contact conditions remain constant.And several iterations may be necessary to find a solution with acceptable accuracy.The modified Newton-Raphson method is used first to evaluate the initial out-of-balance load vector at the beginning of the iteration at each load step.The outof-balance load vector is defined as the difference between the applied load vector R and the vector of restoring loads r i R .When the out-of-balance load is non-zero, the program performs a linear solution, using the initial out-of-balance loads, and then checks for convergence.If the convergence criteria are not satisfied, the out-of-balance load vector is reevaluated, the new contact conditions and the stiffness matrix are updated, and a new solution is obtained.This iterative procedure continues until the solution converges.The modified Newton-Raphson method and its flowchart are outlined by Fig. 6.

Contact stiffness identification using a dynamic approach
First, the dynamic method is studied for use in the estimation of normal contact stiffness.
The results of the dynamic methods are compared with the results based on the static test of normal contact stiffness; then the validated dynamic test method is used in estimation of tangential contact stiffness.

Theoretical formulation of 1-D normal contact stiffness
The idea behind the identification of normal contact stiffness is that the contact interface is modeled by a discrete linear spring.When the preload is changed, contact stiffness will change.When body I is in contact with the ground, the dynamic model of the entire structure can be shown as in Fig. 7.According to this theoretical model the relationship between natural frequencies and normal contact stiffness can be established.When natural frequencies are obtained from impact test, along with a theoretical model, normal contact stiffness can be estimated.
In the one-dimensional model of body I, m is the mass of body I, k is the contact stiffness, p is the preload, f(t) is impulse excitation, u(x,t) the longitudinal displacement of the bar at distance x from a fixed reference.With use of a bar in Fig. 7, the governing equation of the longitudinal vibration of the bar can be expressed as The boundary conditions of the bar are: and at x=l: Initially, the system starts from rest, from the static equilibrium position of the bar, such that the initial displacement condition is: The response of a system to an impulsive force can also be obtained by considering that the impulse produces an instantaneous change in the momentum of the system before any appreciable displacement occurs.The second initial condition is - 2 is called the separation constant and is designated to be negative (De Silva, 1999).
Therefore, the mode shapes X (x) satisfies   whose general solution is According to the general solution and the modal boundary conditions, one can get Set the structure stiffness as * EA k l  and the ratio of the stiffness as . Since the structure stiffness k* is constant and known, the ratio of the stiffness  is proportional to the contact stiffness k n .Therefore Eq. 11 can be expressed as This transcendental equation has an infinite number of solutions  i (i=1,2,…)that correspond to the modes of vibration.When  is changed, the solution of i will change.When  is changed from 0.1 to 10, one can get the corresponding  i l as shown in Fig. 8.The natural frequencies can also be obtained using In an experimental study, the natural frequencies can be obtained by an impact test. i can calculated from Eq.31 since the natural frequencies are related to the system characteristics.
Then  can be determined from Eq. 30.Finally, the contact stiffness, k n , can be estimated based on the definition of .According to the assumption that contact stiffness is a function of the preload, the natural frequencies can be determined in experiments under different preloads.The change of contact stiffness can then be identified based on the change of the through measurement of the natural frequency variation.It should be noted that although any mode of the natural frequency can used to estimate the contact stiffness, some modes might be more sensitive than others to the change of the preloads.

Experimental procedure and results
The experiments were conducted in order to verify the method of identifying contact stiffness in the normal direction (Zheng, 2005(Zheng, , 2008)).The measurement instrumentation includes the proximity, the impact hammer with a load cell, power supply, and a Fast Fourier transformation (FFT) analyzer, as shown in Fig. 9.The experimental procedure can be expressed as follows: 1. Frequency response function (FRF) of the bar is measured by using the hammer to excite the system.Thus, the natural frequencies of the bar can be obtained.

According to the natural frequency equation
3. Based on the relationship between  i l and  of the first three modes in Fig. 8, the  can be inferred from the comparison of experimental results and theoretical results.Then the normal contact stiffness can be obtained from the equation When the natural frequencies are obtained from the experiment, along with the curves of the relationships between  i l and , contact stiffness can be determined from each mode of vibration.However, when the preload changes, the natural frequencies may not necessarily change significantly with the change of normal load for certain modes.Contact stiffness should be identified from the mode most sensitive to changes of a preload.Fig. 10 shows the FRF of the test system under different preload.Fig. 11 shows the relationships between the natural frequencies and preload.The natural frequency of the third mode f3 is the most sensitive to changes in a preload.It can be seen that the results from the dynamic tests are consistent with the numerical calculation results based on the static test results.When the results of dynamic test are consistent with the static test results, the dynamic test method can be used in identifying tangential contact stiffness, for which the static tests are too difficult to conduct.

Theoretical formulation of tangential contact stiffness
Two fixture components are in contact at a certain number of asperities due to the inherent roughness of the surface.When they are subjected to tangential forces, the components are mutually constrained through frictional contacts.A friction model based on the Coulomb friction theory is shown in Fig. 13.The tangential contact stiffness results from the elasticity of asperities of the contact surfaces, and the total resulting stiffness of these contact surfaces depends on their statistical topographical parameters.

Slide area
Stick area Slide area

Fig. 13. A Friction Model
Consider that body I is brought into contact with the flat surface of the support under a uniform preload, P, and is subjected to an small excitation, F, as shown in Figure 14.It is assumed that the tangential contact stiffness will change as the preload increases.The friction at each contact point is governed by Coulomb's law.When force is applied in the tangential direction, the asperities in body I will also deform until the shear stress between the asperities exceeds the limit, then the contact surface will slide each other.The friction model of body I in contact is shown in the Fig. 15.The friction force is given by Eq. 32.
The idea of the identification of the tangential contact stiffness is to compare the two sets of system natural frequencies: one set is identified from the measured impulse response in tangential direction under different preloads, and the other set is calculated from the FEA model of the system.Based on the numerical simulation, a relationship between tangential contact stiffness and the natural frequencies can be established.If the natural frequencies are measured in the experiments under different preloads, the contact stiffness can be calculated from the relationship obtained by the numerical simulation.In order to do the numerical simulation, the effect of the contact force needs to be included into FEA model of the system.The additional contact stiffness matrix will be introduced in the general FEA model.The derivation of contact stiffness matrix is briefly given as follows.
Consider an elastic body I in Fig. 15, the kinetic, strain, and potential energies of the system can be expressed respectively as: where K is the kinetic energy; {u} is the displacement vector; V is the volume of the elastic body I;  is the mass density of the material; U is the strain energy;    and {}; are the strain and stress components, respectively.W is the potential energy of external forces;   F is the external surface force vector specified on the boundary S 1 ;   ˆc R is the contact force vector on the contact surface S c .Note that ; The body force is ignored.Using the above energy expressions the total potential energy of the system is Based on the well-known Hamilton's principle, a discretized FEA formulation for a typical element can be expressed as To obtain the matrix form, the displacement field of a typical element {u}, which is a function of both space and time, can be written as: where [N(x)] is a vector of the space function; and {d(t)} is the nodal response vector.Using the interpolation relationship the element, where [D c ] is the contact property matrix.In the section, the displacements of contact element in the normal direction are assumed to keep stick.Therefore, the normal contact stiffness becomes infinity.The tangential contact stiffness is considered.
The derived contact stiffness matrix should be added to the general FEA model for the fixture stiffness analysis to take into account the effects of the contact force.Followed the standard procedure of the eigenvalue problem, the system natural frequencies can be obtained using the FEA method to establish the relationship between the tangential contact stiffness and natural frequencies.For example, a specimen that has the dimensions 530.75 in was used to measure dynamic characteristics.Fig. 16 shows the FEA model of the specimen.Contact elements were modeled as separate springs on the top and bottom surfaces of the specimen.There are two nodes for each contact element.One node is on the contact surface of the specimen.The other node is constrained at all degrees of freedom.The impulse force was applied at the side of the specimen.The response was obtained at point M, at the other side of the specimen.

Conclusions
Forces in a workpiece-fixture system have a crucial impact on the deformation and accuracy of the system.In this chapter, an FEA model of fixture unit stiffness is proposed.A contact model between fixture components are utilized for solving the contact problem encountered in the study of fixture unit stiffness.By several simple experiments and comparison with the corresponding analytical solution and experimental results in the literature, this methodology is validated.This analytic approach also can be extended in the research of complex fixture system with multiple units and components, which will lead to a new progress in the design and verification of fixture-workpiece system study.

Fig. 5 .
Fig. 5. Sketch of Contact Force on the Contact Surface

n
is the gap between a pair of the potential contact points.In each increment of load, the gap status and the stiffness values are iteratively changed until convergence.As the load is increased, n  will change and hence should be adjusted initial gap before any deformation and T n  is the gap change caused by the total combined normal movement at the pair of points.

Fig. 8 .
Fig. 8. Relationships between the Nondimensional Natural Frequencies and the Stiffness Ratio  in the First Four Modes

Fig. 14 .
Fig. 14.Body I on the Support

Finite
] is the geometry matrix.Comparing to the standard FEA formulation an additional term of C K  , referred as the contact stiffness matrix in included in Eq. (37).The term stems from the work done by the contact force on the contact surface.A brief derivation is presented as follows.The work done by the contact force on the contact surface can be written as Eqs.(38) and (42) into (41) yields www.intechopen.com

Fig. 16 .
Fig. 16.Finite Element Model of Specimen Fig.17 shows the relationships between tangential contact stiffness and natural frequencies of the first two vibration modes.The results are obtained through numerical simulation.From experiments, the frequency response is measured under the different preloads.The contact stiffness can be determined based on the relationships shown in Fig.17.
Fig.17 shows the relationships between tangential contact stiffness and natural frequencies of the first two vibration modes.The results are obtained through numerical simulation.From experiments, the frequency response is measured under the different preloads.The contact stiffness can be determined based on the relationships shown in Fig.17.

Fig. 17
Fig. 17.The Relationship of Tangential Stiffness vs. the First Two Natural Frequencies Fig. 17.The Relationship of Tangential Stiffness vs. the First Two Natural Frequencies