Using a Point-Line-Plane Representation for Unified Simulation of Planar and Spherical Mechanisms

Kinematic simulation of a mechanism requires calculation of the position and orientation of all of its constituent links during its entire range of motion. Simulation methodologies can be broadly classified into three categories: graphical, analytical, and numerical [1]. The graphical analysis method is based on dyadic decomposition, i.e., identification of four-bar loops in mechanisms [2]. Although this approach is prominently used in simulation packages like Linkages [3] and PMKS [4], its limitations are well known [5], e.g., they are unable to handle complex mechanisms like a double butterfly mechanism. Analytical methods involve solving a loop closure constraint-based system of non-linear equations [6]. Most analytical methods use the polynomial continuation method [7,8], elimination method, or Grobner bases [9] to solve the simulation problem. These methods are able to find all the possible assembly configurations of a given mechanism. However, they are not general in nature since the motion equations need to be derived for each type of mechanism separately. As a result, these approaches cannot be used to simulate n-bar planar or spherical mechanisms without manually deriving equations on a case-by-case basis.

Numerical simulation methods are iterative in nature and can handle extremely complex mechanism [10,11]. They use numerical methods such as the Newton–Raphson method to solve the system of non-linear equations for one solution only instead of all possible ones [12]. These methods accept the mechanism joints and link information as inputs. Subsequently, the algorithm repeatedly solves the finite displacement problem, i.e., the input link is iteratively moved with finite displacement, and consequently, positions of remaining links are calculated. As a result, the entire range of motion for the specified mechanism is calculated.

Hernández and Petuya have proposed a geometrical-iterative method that performs better than Newton–Rhapson method [13]. However, the approach is limited to n-bar planar mechanism with revolute joints only. Radhakrishnan and Campbell [14] have created a computational tool for planar mechanism, which carries out position analysis of planar mechanisms using a geometrically iterative algorithm. However, due to the use of dyadic decomposition, it shares the limitations of graphical methods and is limited to the planar mechanisms.

Commercial CAD software such as autodesk inventor, solidworks, adams, etc. solve differential-algebraic equations numerically to provide multibody simulation capability [15–17].² However, their use is more prominent during detail design stage rather than the conceptual design stage. Creation of feature-based assembly of planar and spherical mechanisms and initializing constraints on these systems is a nontrivial task. Changing the type of joints or the number of links for a mechanism is also more involved than carrying out the same operation on purely kinematic simulators like PMKS [4]. Additionally, their solvers model the motion problem as a set of coupled differential and algebraic equations. This type of model is more suited for dynamic simulations rather than kinematic simulations, which involves purely algebraic constraints. Also, the algebraic equations for commercial softwares are modeled using reference point representation that leads to more number of constraints when compared with other representations. Thus, use of these softwares for concept design is not ideal. SAM and GIM are two more packages that support n-bar simulation for planar linkages with both revolute and prismatic joints [18].³ Furlong et al. [19] have demonstrated a virtual reality environment for simulating spherical four-bar mechanisms, and in the academic domain, sphinx, isis, and osiris are softwares that enable the analysis and synthesis of spherical mechanisms [20–22].

However, currently, there are no approaches that unify n-bar planar and spherical mechanism analysis and can be demonstrated to more complex linkage systems. The proposed approach hopes to bridge this gap and augment the capability of pure kinematic design systems like MotionGen [23]. In this paper, planar and spherical mechanisms are represented as a collection of geometric constraints spanned by points, lines, and planes. The geometric mechanism representation enables the design of unified constraint equations that are easily programmed. As a result, a simple generalized real-time framework for mechanism simulation is achieved.

Our previous work on mechanism synthesis problems includes the creation of geometric constraint equations for four-bar mechanisms with revolute or prismatic joints [24–28]. In Ref. [29], we have demonstrated an algebraic fitting-based approach in the space of planar quaternions to simulate planar four-bar linkages. However, that approach does not scale for more complex planar or spherical linkages. In this paper, we show that by using homogeneous coordinates, we can derive unified geometric constraint equations for both planar and spherical linkages, which simplifies the simulation without resorting to calculations for individual types of mechanisms. The rigidity constraints imposed by the links are modeled as simple geometric constraints using points, lines and planes. Once the mechanism is specified, the solver proceeds with iteratively perturbing the input and solving the constraints for other links. To the best of authors’ knowledge, this work is the first attempt at using a point-line-plane mechanism representation and presenting unified geometric constraints for simulation.

The major intellectual contributions of this paper can be summarized as (1) presenting generalized constraint equations for planar and spherical mechanisms using point-line-plane representation and (2) enabling real-time simulation of n-bar planar or spherical mechanisms.

Rest of the paper is organized as follows. Section 2 discusses the representation and constraints required to describe the motion of a general planar and spherical mechanism. Section 3 demonstrates the algorithm required to simulate a mechanism using the iterative numerical approach. Finally, in Sec. 4, we present a few examples to demonstrate the use of proposed algorithm.

Selection of an apt mechanism representation and constraints is important as it has a profound effect on algorithm’s simplicity and efficiency. Conventionally, a multibody system has been specified using multiple representations, namely: relative coordinates, reference point coordinates, and natural coordinates [10]. Relative coordinates are based on parameters specifying one link relative to another; reference point coordinates are based on specifying absolute position of each link independently; while the natural coordinates are based on each link being specified by two points. Using relative coordinates enables a scalable representation while reference point coordinates tend to be more computationally efficient. Natural coordinates provide a compromise between the two approaches in terms of simplicity and efficiency. Most commercial softwares use reference point coordinate representation that usually leads to maximum number of constraint equations and subsequently high computation time.

Figure 1 shows an RRPR (R: Revolute, P: Prismatic) four-bar mechanism and its specification using different representations. For the relative coordinate representation, there are three unknown coordinates, i.e., ψ₁, ψ₂, L₃. For the reference point coordinate representation, the mechanism has nine unknown variables, i.e., location and orientation of each link x_i, y_i, ψ_i. Similarly, for the natural coordinate representation, there are six unknown variables namely x₁, y₁, x₂, y₂, x₃, y₃. Since the four-bar mechanisms are a single degree-of-freedom mechanisms, each of the representation requires two, eight, and five constraint equations, respectively, to fully define the motion. In this paper, we derive unified constraint equations for all types of planar and spherical linkages consisting of both revolute and prismatic joints. We use homogeneous coordinates to write geometric constraints on points, lines, and planes. For example, our representation for the shown RRPR mechanism will require using four unknown point and line coordinates, i.e., x₁, y₁, a₂, b₂ since we can set homogenizing factors z₁ = 1 and c₂ = 1 without any loss in generality. Such a representation would keep the number of unknown variables smaller while also enabling construction of simpler geometric constraint equations. Computational efficiency aside, this representation also naturally maps to the geometric constraints, which for the shown RRPR mechanism is a circle-constraint on the moving pivot and a line-constraint on the fixed pivot of the RP dyad. It is well known that the time complexity of multidimensional Newton–Raphson method, which is used in this paper, is at-least O(n²) for a single iteration where n is the number of unknowns in state variable [11]. This is a direct consequence of a dense Jacobian matrix having n² elements that need to be calculated after every iteration. Thus, more unknowns result in needing to calculate larger Jacobian matrices, which is computationally expensive.

Planar mechanisms can be uniquely specified using their joint and link data. A joint can be prismatic or revolute, which naturally associates with points and lines, respectively. We use homogeneous coordinates to represent both points and lines. Thus, a point P is given by homogeneous coordinates (x, y, z) whose Affine coordinates are given as (x/z, y/z), while a line L is also represented using homogeneous coordinates (a, b, c), where equation of the line passing through the point P in the projective plane is given by ax + by + cz = 0. Depending on the constraint being expressed, this line can be fixed or floating in the plane. A link can be represented by a subset of joints. The link can be binary, ternary, or n-ary depending on the number of joints it contains. An example six-bar planar mechanism is displayed in Fig. 2. Its joints are represented as points and lines while its links are defined as a group of joints as shown in Table 1.

Similarly, spherical mechanisms can also consist of revolute and prismatic joints. A spherical prismatic joint constrains the link movement along a circular arc instead of a line. We represent a spherical revolute joint as a point P in terms of its homogeneous coordinates (x, y, z, w) with respect to the center of the unit sphere such that its Affine coordinates are (x/w, y/w, z/w). A spherical prismatic joint is defined as a plane Pl : (a, b, c, 0) passing through the center of the sphere and is given by the equation ax + by + cz = 0. The intersection of the plane and the unit sphere defines the great circle along which the constituent links are constrained to move for a spherical prismatic joint. Similar to the lines for planar mechanisms, this plane can be fixed or moving depending on the geometric constraint being expressed. An example RRPR spherical mechanism is displayed in Fig. 3. Its joints are represented as points and planes while its links are defined as a group of joints as shown in Table 2.

During the motion, a mechanism is subjected to a set of constraints imposed by the rigidity of its links. Thus, to simulate a mechanism, these constraint equations need to be formulated. For planar and spherical mechanisms, modeling three constraint equations are sufficient for simulation. We propose a unique constraint equation for each of the binary links RR, RP, PR, and PP. But, we will see that all of these constraints can be expressed in a single equation. Any link with more than two joints can easily be reduced to an equivalent collection of binary links. For example, a ternary link can be treated as three binary links. Thus, these constraints can successfully be used to enforce the rigidity of any link in a general mechanism. Figures 4 and 5 show different planar and spherical binary links, which are building blocks for any planar and spherical mechanism and are discussed below.

It can be seen that Eqs. (2), (4)–(7) are all degenerate case of the Eq. (1). In the projective plane for the planar geometric constraints, the lines and points are dual to each other; thus, their meanings can be interchanged without changing the underlying structure of the equations. In the projective three-space for the spherical constraints, the points and planes are dual to each other and thus their meanings can be interchanged. Thus, Eq. (1) is the single equation that unifies all the geometric constraints associated with all types of links for both planar and spherical mechanisms. This facilitates creation of the following metrics for computation: (1) distance between two points in space, (2) perpendicular distance between a point and a plane, and (3) angle between two planes.

In this section, we discuss the algorithmic steps required to solve the kinematic simulation problem. The general approach is to iteratively perturb the input links by a finite displacement and find the new position of the mechanism.

The simulation process involves iteratively perturbing the input link by a finite displacement. Depending on the actuating joint being revolute or prismatic, the displacement could be translation or rotational in nature. In this paper, we restrict ourselves to consider actuation at the fixed joints. The relations governing the motion of input link are derived in this subsection.

With these expressions, we can successfully calculate the location of input link after imparting it a discrete perturbation. The next step is to find the coordinates of all the other unknown joint coordinates that are compatible with the rigidity constraints imposed on the mechanism during simulation.

To automate the calculation of residual vector Φ(q_i) and the Jacobian matrix [J(q_i)], the constraints are handled in a sequential manner. While creating the residual vector in our implementation, first the rigidity constraints for each link are calculated and then the constraints for joints are calculated. Similarly, the Jacobian matrix is created in a column-first manner, i.e., all the partial differential equations with respect to an unknown state variable are calculated before progressing to the next variable. The outlined method is just one way of calculating Φ(q_i) and [J(q_i)] since their values are independent of the sequence adopted to calculate each element.

Thus, using the constraint equations and their first-order partial derivatives, it is possible to solve iteratively for the solution using Newton–Raphson method. The iterations are continued until a solution within desired accuracy is calculated.

For some input link perturbations, the Newton–Raphson method might fail to converge even after many iterations. In these instances, there does not exist a mechanism state which fulfills all the constraint equations. As a result, this input perturbation is outside the possible limits of motion of mechanism.

Thus, by iteratively perturbing the input link and solving the constraints for other joint coordinates, we are able to simulate any general planar or spherical mechanism. Numerous techniques exist that can improve the convergence and efficiency of the Newton–Raphson method. However, the basic method suffices to achieve real-time simulation. The complete algorithm has been described in

Algorithm 1.

This section presents sample examples to demonstrate the use of proposed algorithm for mechanism simulation. The simulation has been carried out in matlab on a PC running Core i5-7300 at 2.6 GHz with 8 GB RAM. The simulation is carried out within seconds for residual value of 1.0e − 8. Each closed-loop output curve is made up of 180 points while open-loop curves have less than 180. Each point corresponds to 2π/180 radian or units input perturbations.

To compare the speed of commercially available CAD systems and proposed algorithm, a planar four-bar crank rocker mechanism with revolute joints is modeled and simulated. Autodesk Inventor 2020 with educational license is used as reference commercial CAD system. Simulation is performed for 180 time-steps with one-degree input link perturbation for each time-step. The simulation takes 5 s to complete Inventor while it finishes in 1.4 s when the proposed algorithm is used. Thus, even for a simple mechanism, there is a significant speed difference between the commercial solvers and proposed methodology.

A planar Stephenson-II six-bar linkage is simulated in this example. This linkage does not have a four-bar linkage and proves to be challenging to simulate using dyadic decomposition based approaches. However, our approach handles these non-dyadic mechanisms without issues.

The six-bar mechanism is displayed in Fig. 2, and its joint and link data are given in Table 1. The mechanism has J₁, J₇ as the fixed joints, J₂ as the perturbed joint and J₃, J₄, J₅, J₆, J₈ as the unknown joints defining the 11-dimensional state vector. The mechanism consists of ten rigidity constraint equations and one homogeneous coordinate equation for prismatic joint. The simulation algorithm successfully solves these constraints and plots the trajectory of the coupler point J₈ as shown in Fig. 2. The run-time of this simulation was 3.79 s.

In this example, a planar-modified Theo Jansen linkage with one of its revolute joints replaced by a floating prismatic joints is simulated. The eight-bar mechanism is displayed in Fig. 6 and its joint and link data is given in Table 3. J₁, J₅ are the fixed joints, J₂ is the perturbed joint and the state vector consists coordinates of J₃, J₄, J₆, J₇, J₈. This results in a 11-dimensional state vector. Ten rigidity constraint equations for links and one homogeneous coordinate equation for prismatic joint are available for this mechanism. The simulation algorithm plots the trajectory of the coupler point J₈ as shown in Fig. 6. Note, the length of stride for this modified mechanism is larger than that of the conventional Theo Jansen mechanism which has revolute joints only. As a result, this mechanism is a prospective candidate for walking robots. The run-time of this simulation was 3.81 s.

This example presents the simulation of a spherical RRPR mechanism which is the spherical analog of the Whitworth quick-return mechanism. The four-bar mechanism is displayed in Fig. 3, and its joint and link data are given in Table 2. The mechanism has J₁, J₄ as the fixed joints, J₂ as the perturbed joint, and J₃, J₅ as the unknown joints defining the six-dimensional state vector. The mechanism consists of four rigidity constraint equations for output and coupler links which can be described using Eqs. (1) and (4). Also, one unit sphere equation for revolute joint using Eq. (9) and one homogeneous coordinate equation for prismatic joint using Eq. (8) can be written. Once the simulation is completed, the trajectory of coupler point J₅ can be plotted as shown in Fig. 3. The run-time of this simulation was 1.49 s.

In this example, a spherical Watt-I six-bar linkage with prismatic input joint is simulated. Spherical Watt I type linkages have been used to design door hinges for spatial movement. The six-bar mechanism is shown in Fig. 7, and its link and joint data are given in Table 4. From the data, it is known that J₁, J₆ are the fixed joints, J₂, J₃ are the perturbed joints, and J₄, J₅, J₇, J₈ are the unknown joints representing the 12-dimensional state vector. The mechanism can be described using eight rigidity constraints for links and four unit circle constraints. Perturbing the input link along the input prismatic joint results in the motion of coupler point J₈ as shown in Fig. 7. The run-time of this simulation was 3.33 s.

In this section, we demonstrate the scalability of proposed algorithm to spatial mechanisms by simulating a 5-SS platform linkage. A 5-SS mechanism consists of five binary spherical-spherical (SS) links connected to a floating coupler link on one end and the ground link on the other [30]. Although using the Grüebler criterion [2], the mobility of this mechanism is six, five of the rotational degrees for each binary link are redundant, and therefore, 5-SS platform linkage is a one degree-of-freedom mechanism.

An example 5-SS mechanism is displayed in Fig. 8 and its joint and link data is given in Table. 5. Since the input link perturbation approaches outlined in Sec. 3.1 only cover revolute and prismatic joints, a new approach is needed for spherical input links. We attach a linear actuator between J₁ and J₇ to actuate the mechanism [30]. This results in the mechanism having J₁, J₂, J₃, J₄, J₅ as fixed joints and J₆, J₇, J₈, J₉, J₁₀, J₁₁ as the unknown joints defining the 18-dimensional state vector. Five rigidity constraints for binary links, 12 independent rigidity constraints for coupler and one additional constraint for input actuator length is available for this mechanism. All these constraints are modeled using Eq. (1) and the simulation is carried out to generate a spatial trajectory of coupler point J₁₁. The input prismatic link is perturbed by 0.01 units and run-time of this simulation was 0.28 s.

In this paper, we have presented unified equations for motion simulation of planar and spherical n-bar mechanisms and an efficient algorithm for computation to enable real-time, interactive simulation. The approach is general and uses simple geometric primitives, such as point, line, and planes to represent the constraints inherent in mechanisms. A 5-SS mechanism is simulated to demonstrate the scalability of the proposed approach to spatial mechanisms. Once the mechanism is simulated and the path of coupler point determined, velocity and acceleration curves can easily be determined using numerical differentiation. Future research would involve finding appropriate representation and rigidity constraints for cylindrical and helical joints to further unify spatial synthesis.

This work has been financially supported by The National Science Foundation under a research grant # CMMI-1563413 to Stony Brook University. All findings and results presented in this paper are those of the authors and do not represent those of the funding agencies.