Design General Cam Proﬁles Based on Finite Element Method

: This paper presents a general framework to design a cam proﬁle using the ﬁnite element method from given displacements of the follower. The arbitrarily complex cam proﬁle is described by Lagrangian ﬁnite elements, which are formed by the connectivity of nodes. In order to obtain the desired proﬁle, a penalty-type functional that enforces the prescribed displacement of the follower is proposed. Additionally, in order to ensure convexity of the functional, a numerical stabilization scheme is used. The nodal positions are then obtained by solving a nonlinear system of equations resulting from minimizing the total functional. The geometrical accuracy of the cam proﬁle can be controlled by the number of ﬁnite elements. A case study is considered to illustrate the ﬂexibility, accuracy, and robustness of the proposed approach.


Introduction
Cam mechanisms are widely used in various mechanical machines, such as machine tools, sewing machines, and engines due to their unlimited motions, operation speed, motion accuracy, and structural rigidity. In cam design, the follower motion affects the kinematic and dynamic characteristics of cam systems. Various traditional functions can be used to describe the follower motion such as harmonic, cycloidal, trapezoidal, polynomial, and piecewise polynomial [1,2]. So far, spline functions have also been used to develop the cam curves, for example, the Bezier curve [3], the B-spline curve [4][5][6], and the NURBS curve [7][8][9][10]. Besides, kinematics and dynamics have been considered to optimize the peak values of the displacement, velocity, acceleration, and jerk with referring to several functions, e.g., harmonic, cycloidal, and polynomial [11]. Qin and Chen [12] also proposed the kinematic optimization of cam motions for engine valve trains regarding the characteristics of polynomial cam curves. In order to achieve motion curves with good dynamic characteristics under arbitrary design conditions, the parameters of the motion model were obtained by optimizing the dynamic performance [13,14]. Related to the cam design process, the cam size optimization has been also investigated [15][16][17]. This plays an important role due to affect the cost of cam systems.
The design of the cam profile is one of the main tasks in designing a cam mechanism. The profile of the cam must be designed correctly to ensure the characteristics of kinematics and dynamics of cam mechanisms. However, when a priori displacement of the follower is given, obtaining the corresponding curve equations for cam profiles is nontrivial in general. Many approaches to design the profile of a cam have been proposed. The graphical method has been established by several studies for disk cams [1,18,19]. In this method, kinematic inversion is used. To increase the accuracy, the analytical method using a Cartesian coordinate system is presented by [20][21][22]. Since then, the analytical method has been continuously advanced, among others by Shigley and Jr. [1] and Norton [23]. In these studies, the equations of cam profile have been formulated by the loop-closure equation and using complex polar notation. Wu et al. [24] presented a cam profile of a translating cam driven. In this approach, the follower motion can be arbitrarily chosen for only the forward or the return stroke. The desired cam profile is obtained by using the theory of envelopes [25]. To this end, the envelopes are given by fitting a tangent curve to the successive follower positions, and the cam profile is achieved by the boundary of the family. Biswas et al. [26] presented two approximate analytical methods for determining disk cam profiles of roller followers. In the first method, the average location is defined by two lines tangent of three roller positions. The cam profile is obtained by the interpolation of the roller positions. In the second method of this research, the cam profile is achieved by the intersection of the curvature center and the roller center. A systematic method for the design and analysis of cams with three circular-arc profiles is presented by Jung-Fa Hsieh [27]. This study is based on constructing a generic geometric model of cam profile and deriving the corresponding equations of the constituent circular arcs using a coordinate transformation. Moreover, using coordinate transformation, the cam profile is obtained by [28]. This study proposes a simple yet comprehensive method for the design and analysis of an indexing cam mechanism with parallel axes.
In this paper, we propose a finite element method to design a cam profile from a given displacement of the follower. To the best of our knowledge, this is the first time a finite element description is used directly for obtaining the cam profile. Compared to earlier works, the presented method contains the following merits and novelties:

•
A penalty type functional to enforce the prescribed displacement of the follower is proposed; • A stabilization scheme is added to guarantee convexity of the functional; • For the Newton-Raphson solution method, the consistent linearization of the weak form is presented; • The arbitrary desired accuracy of the cam profile can be achieved by increasing the number of finite elements; • The resulting profile represented by finite elements can be used for cam manufacturing.
The remainder of the paper is structured as follows. In Section 2, the geometrical description of cam profiles is presented. Section 3 introduces the synthetical kinematics. The formulation of the general framework for the synthesis of cam mechanisms is presented in Section 4. The finite element (FE) formulation is discussed in Section 5. Section 6 then presents a solution of the finite element. The case study for this method is presented in Section 7. Finally, conclusions are drawn in Section 8.

Geometrical Description of Cam Profiles
In this section, we present the geometrical description of the cam profile. The differential geometry (see [29]) is used for an arbitrarily complex cam surface. The cam surface, denoted by S (see Figure 1), is represented in curvilinear coordinates as where ξ ξ 1 , ξ 2 are the curvilinear coordinates. Accordingly, the covariant and contravariant tangent vectors, a α and a α , at point x depending on S are calculated by and respectively, where a αβ := [a αβ ] −1 , a αβ := a α · a β are the contravariant and covariant components of the metric tensor. Here and henceforth, the summation convention is implied for repeated Greek indices, e.g., α = 1, 2. The unit normal vector n at x is defined by [30] n = a 1 ×a 2 Appl. Sci. 2021, 11, x FOR PEER REVIEW 3 of 16 where   ,  are the curvilinear coordinates. Accordingly, the covariant and contravariant tangent vectors, and , at point depending on S are calculated by and respectively, where a : a , a : are the contravariant and covariant components of the metric tensor. Here and henceforth, the summation convention is implied for repeated Greek indices, e.g., = 1, 2. The unit normal vector at is defined by [30] ‖ ‖ .
It can be shown that The curvature tensor can be computed as (see, e.g., [29] where b : , .

Synthetical Kinematics
In this section, the established approach to describe the deformation of membranes (see, e.g., [29]) is employed in the context of synthetical kinematics for designing cam It can be shown that a 1 × a 2 = deta αβ .

Synthetical Kinematics
In this section, the established approach to describe the deformation of membranes (see, e.g., [29]) is employed in the context of synthetical kinematics for designing cam profiles. To obtain the desired cam profile, a given initial configuration is selected. The initial curve is then "deformed" to the desired configuration. The differences between the initial (reference) configuration S 0 and the current configuration S are shown in Figure 1. Both configurations are described by the geometrical quantities discussed in Section 2. For the reference configuration S 0 we use the upper-case symbols X, A α , A αβ , A α , and N. The corresponding lower-case symbols x, a α , a αβ , a α , and n are used in the current configuration S. In order to characterize the deformation between the reference configuration S 0 and the current configuration S, a line element in S is considered as and its corresponding line dX = A α dξ α in S 0 . So that we can write dξ α = A α dX. With this, dx in the current configuration can be mapped to dX in S 0 by Here the tensor denotes the deformation gradient. Vice versa maps the dX in S 0 to dx in S. Through F we thus have the following transformation Finally, the stretch between configurations S 0 and S is defined by Here, "det" denotes the surface determinant.

Governing Equation for the Synthesis of Cam Profile
In this section, the governing equation for designing cam profile is proposed. The strong form, weak form, and stabilization scheme are presented.

Strong form Statement
The point-wised governing equation for the synthesis of cam profile is given by subjected to the boundary condition where s 0 and s denote the reference and current displacement of the follower, respectively. s 0 and s will be explicitly defined for a particular type of cam. Γ u is a subset of S. From Equation (15), the variation δg can be expressed as where The strong form equation as shown in Equations (15) and (16) is not a well-posed statement since it is still missing an in-plane constraint equation. In another word, the positions of two arbitrary material points x 1 ∈ S\Γ u and x 2 ∈ S\Γ u can be exchanged without affecting the solution of Equation (15). The system is thus unstable and needs stabilization, which will be discussed in Section 4.4.

Penalty-Type Functional
The point-wise condition in Equation (15) can be enforced by a penalty-type functional. The penalty-type functional can be expressed as where is the penalty parameter. The variation of the penalty-type functional is computed as where λ := g (21) represents the pressure of the enforcement. It should be noted that the potential in Equation (19) exactly enforces the condition at the limit, i.e., → ∞ . At → 0 , Π p does not affect the system. It turns out that this property is advantageous for the Newton solution method since it allows for the application of s in an increment manner until a desired accuracy is reached.

Augmented Lagrange Functional
In order to enforce the condition in Equation (15) stronger, the Augmented Lagrange method can be used instead of the Penalty method. Accordingly, the pressure λ (see Equation (21)) is further corrected by using extra loops called augmentation steps. For this, we consider the functional where λ is a non-negative parameter. Thus, its variation reads: From Equation (23), one finds the updated new parameter λ new for Equation (20) as

Stabilization Governing Functional
As noted above, the system is unstable and thus cannot be solved uniquely with FEM. In this paper, a stabilization method is considered as found in [31]. The stabilization functional can be expressed by where µ is a given parameter. The variation of Π stab can be found as where, τ αβ = A αβ − a αβ .

Stabilized weak Form Statement
The cam surface is assumed to be governed by the weak form where Π p and Π stab are defined in Equations (19) or (22), and (25), respectively. The cam profile is obtained by minimizing potential (27). I.e., Here, δΠ p and δΠ stab are defined in Equations (20) or (23), and (26), respectively.

Linearized Weak Form
The linearization is necessary for the calculation of tangent matrices. The linearization of the potential energy as shown in Equation (20) can be computed by Here ∆δg can be expressed in the form

Linearized Stabilization Term
From Equation (26), the linearization of the stabilization equation can be computed by

Finite Element Discretization
Equation (28) is solved by the finite element method. The surface S is therefore discretized into a set of finite elements Ω e that are defined by nodal points x i which corresponds to the Lagrangian finite element description (see Figure 2).
The variation of the unit normal vector is given by [30] δ δ , with  .
The linearization of vector is followed from Equation (37) as where ∆a α is can be found in [31] ∆ ∆ with a   , and ∆ is the same manner as δ from Equation (36). The variation of J is given respectively in [29] δJ J . δ .
Further, the variation and the linearization of a and a are:

Discretized Synthetical Kinematic Quantities
Within elements Ω e , the geometry is approximated by the nodal interpolations where N i = N i (ξ α ) and n denote the shape function and number of nodes, respectively. For shorthand notation, Equations (32) and (2) can be expressed as follows: The variation of the unit normal vector is given by [30] with The linearization of vector n is followed from Equation (37) as where ∆a α is can be found in [31] ∆a α = B αβ ∆a β with B αβ = a αβ (n ⊗ n) − a β ⊗ a α , and ∆a α is the same manner as δa α from Equation (36). The variation of J is given respectively in [29] δJ = J a α · δa α .
Further, the variation and the linearization of a αβ and a αβ are:

Discretized Weak Form of the Governing Equation
With the expressions in Section 5.1, the weak form is now discretized. Accordingly, Equation (20) can be written as follows: where f e p : = In Equation (46), g p is to be defined for a given type of cam. The linearization can be expressed as where kp is the tangent matrix. It can be written as

Discretized Stabilization Terms
The weak form of the stabilization equation is discretized. The formula as shown in Equation (26) can be expressed as The linearization of the stabilization functional can be expressed as where the tangent matrix can be written by

Solution Procedure
The discretized weak form as shown in Equations (45) and (49) of the governing equation and the stabilization equation can be written as The nonlinear in Equation (53) can be solved by numerical methods. The Newton-Raphson method is used to solve this equation. Equation (53) can be described as At iteration i, an approximation solution of cam profile can be sought and denoted by x i . The method can be described by Gouri Dhatt, et al. [32] The algorithm is obtained by developing this residual into a Taylor series in the vicinity of x i−1 and By ignoring terms greater than one, we get where is the tangent matrix.

Numerical Quadrature
For finite element computations, the integral for the weak form and the tangent matrices must be evaluated. These integrations can be performed on the element level. Since the integration is carried out in the finite element analysis of the reference configuration, the parameter domain coordinate ξ α ∈ [−1, +1] must be transformed to this configuration (see Figure 3).

Numerical Quadrature
For finite element computations, the integral for the weak form and the tangent matrices must be evaluated. These integrations can be performed on the element level. Since the integration is carried out in the finite element analysis of the reference configuration, the parameter domain coordinate  ∈ [−1, +1] must be transformed to this configuration (see Figure 3). The Jacobian can be calculated as follows [30]: The integration form is then carried out with the standard Gaussian quadrature on the parameter domain. Thus, from Equations (46)   The Jacobian can be calculated as follows [30]: The integration form is then carried out with the standard Gaussian quadrature on the parameter domain. Thus, from Equations (46)  (63) The integration will be done numerically. n gp is representative for the number of Gauss points. Hence, the integral in Equation (63) will be approximated by Here, W p is weighting factors and ξ gp denotes the coordinates of the evaluation points. The integration for the whole system can be computed as follows: (65)

Computational Algorithm
This section presents a computational to solve a cam profile. With the above-discussed framework, the computational algorithm of a cam profile that satisfies a prescribed displacement of the follower is given in Table 1. Step 1: Load step for penalty method • Step 2: Reduce parameter µ • Step 3: Update new parameter λ of Augmented Lagrangian method following as equation

Case Study
This section presents a case study to demonstrate the procedure of this work. We consider the cam mechanism with a flat-face follower for illustration. The finite element is therefore reduced to line elements. With the general designing framework discussed above, an explicit expression for the follower displacement is required. The expression is, however, generally dependent on the type of cam. The cam curves resulted from the proposed framework are compared with the analytical solution.

Computation of Cam Flat-Face Follower
With the flat-face follower, as shown in Figure 4, the function of the follower displacement s 0 can be written as follows: where ϕ is the angle of the camshaft, and the displacement function f(ϕ) is shown in Appendix A.
Appl. Sci. 2021, 11, x FOR PEER REVIEW 12 of 16 From the calculation above, can be simply obtained by where, ≔ . The tangent matrix is .
Here ϵ   dA For cam plate, we can describe the tangent vector and others with α = 1 and β = 1. ξ 1 is denoted by ξ. From the variation of the stabilized equation (see Equation (26)), f e stab can be computed as follows: Likewise, the tangent matrix is calculated by where k e stab1 = Ω e 0 2µ a 11 N T , a 1 a 1 N T ,ξ dA e , As shown in Figure 4, the displacement of the follower can be computed as Taking the variation of the current displacement in Equation (70) is presented by The angle of the camshaft can be written as (see Figure 4) ϕ = arccos(n ·e 1 ).
The variation of s 0 can be also computed by In Equation (73), d c is calculated in Appendix A. From the calculation above, f e p can be simply obtained by where, R A := R 1 A . The tangent matrix is k e p = k e p1 + k e p2 . Here and Tensor b in Equation (76) and tensors G 1 , G 2 in Equation (77) are defined in Appendix A.
In the example, the cycloid function for the displacement of the follower as shown in Equations (A1) and (A2) in Appendix A is used for the design of the cam profile. The types of motion programs are rise, fall, and dwell. These motions are defined by the angle of segments, i.e., β 1 = β 2 = β 3 = 120 • and the lift L = 5 mm. By using the computational algorithm in Section 6.3, the cam profile is shown in Figure 5a and the kinematic characteristics of the follower including displacement, velocity, and acceleration are presented in Figure 5b-d. The finite element solution of the cam profile is compared with the analytical solution as shown in Figure 5. From this figure, it can be seen that there are very small differences between the FE solution and the analysis solution.
In order to determine the minimum base circle, the parameters including ρ , s @ , and a need to be computed. The minimum radius of curvature ρ is calculated from Equation (7). Values of s @ and a are determined from the displacement and acceleration functions. Substituting values of the parameters into Equation (78), the minimum base circle is 3.9878 mm. In addition, with the cam flat-face follower, the radius of curvature must be positive to maintain the contact between the cam and the follower.

Influence of Parameters and μ
Next, the influence of penalty parameters ϵ and μ is discussed in this section. The accuracy of the cam profile is estimated by parameter as shown in Figure 6a. It is clearly seen that an accurate solution is obtained for a sufficiently large ϵ. However, it should be noted that for the augmented method, a much smaller penalty parameter is already sufficient for the comparative accuracy due to the extra correction loop. Figure 6b further shows the convergence of the finite element calculation with mesh refinement for three minimum values of the parameter μ = 50, 25, and 10. It is observed In this example, 128 elements are used to compute the cam profile. In the penalty step, the maximum difference between the FE displacement and the input displacement is small, i.e., ε = 1.9 × 10 −3 . For the second step, the stabilization parameter µ is reduced from 100 to 25, which increases the accuracy in the displacement to ε = 7.0282 × 10 −5 . For the augmented step, the accuracy is further improved with ε = 2.6133 × 10 −5 .
Next, the size of the cam is considered in the cam design process since it decides the cost of the cam follower. The cam size is influenced by the base circle r O . Usually, a minimum base circle, denoted by r O min , can be defined by (see Robert L. Norton [23]) where ρ min , a min , s @ amin are the minimum radius of curvature, minimum acceleration, and the displacement value at a min , respectively. In order to determine the minimum base circle, the parameters including ρ min , s @ amin , and a min need to be computed. The minimum radius of curvature ρ min is calculated from Equation (7). Values of s @ amin and a min are determined from the displacement and acceleration functions. Substituting values of the parameters into Equation (78), the minimum base circle is 3.9878 mm. In addition, with the cam flat-face follower, the radius of curvature must be positive to maintain the contact between the cam and the follower.

Influence of Parameters and µ
Next, the influence of penalty parameters and µ is discussed in this section. The accuracy of the cam profile is estimated by parameter as shown in Figure 6a. It is clearly seen that an accurate solution is obtained for a sufficiently large . However, it should be noted that for the augmented method, a much smaller penalty parameter is already sufficient for the comparative accuracy due to the extra correction loop.
Appl. Sci. 2021, 11, x FOR PEER REVIEW 14 of 16 from the figure that the rate of convergence increases with decreasing the parameter μ.
On the other hand, for the fine meshes, parameter μ affects the convergent rate significantly, i.e., the smaller stabilization parameter, the more accuracy is obtained.

Conclusions
This paper presents a new and general computational formulation for the synthesis of cam mechanisms using the finite element method. The formulation is based on a general theoretical description of the cam curve in curvilinear coordinates. The governing equation for the synthesis of cam profile is proposed. The penalty-type functional is proposed for the weak form. The stabilization functional is further added for convexity of the total potential. The discretization of the weak form is derived in the framework of the finite element. The Gaussian quadrature is used to integrate the finite element computation. The cam profile is obtained by solving a nonlinear system of equations resulting from minimizing the total functional.
The case study of the flat-face cam mechanism demonstrates that the formulation is able to handle all aspects accurately, robustly, and efficiently. The proposed framework can be applied for more general cam mechanisms, which is a subject of future work.  Acknowledgments: This work is supported by Thai Nguyen University of Technology (TNUT).

Conflicts of Interest:
The authors declare no conflicts of interest. They do not know of any competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Appendix A
This appendix presents the displacement function of the follower f φ , vector and tensors , in Section 7. Figure 6b further shows the convergence of the finite element calculation with mesh refinement for three minimum values of the parameter µ = 50, 25, and 10. It is observed from the figure that the rate of convergence increases with decreasing the parameter µ. On the other hand, for the fine meshes, parameter µ affects the convergent rate significantly, i.e., the smaller stabilization parameter, the more accuracy is obtained.

Conclusions
This paper presents a new and general computational formulation for the synthesis of cam mechanisms using the finite element method. The formulation is based on a general theoretical description of the cam curve in curvilinear coordinates. The governing equation for the synthesis of cam profile is proposed. The penalty-type functional is proposed for the weak form. The stabilization functional is further added for convexity of the total potential. The discretization of the weak form is derived in the framework of the finite element. The Gaussian quadrature is used to integrate the finite element computation. The cam profile is obtained by solving a nonlinear system of equations resulting from minimizing the total functional.
The case study of the flat-face cam mechanism demonstrates that the formulation is able to handle all aspects accurately, robustly, and efficiently. The proposed framework can be applied for more general cam mechanisms, which is a subject of future work.