## Abstract

This paper presents a geometric constraints driven approach to unified kinematic simulation of *n*-bar planar and spherical linkage mechanisms consisting of both revolute and prismatic joints. Generalized constraint equations using point, line, and plane coordinates have been proposed which unify simulation of planar and spherical linkages and are demonstrably scalable to spatial mechanisms. As opposed to some of the existing approaches, which seek to derive loop-closure equations for each type of mechanism separately, we have shown that the simulation can be made simpler and more efficient by using unified version of the geometric constraints on joints and links. This is facilitated using homogeneous coordinates and constraints on geometric primitives, such as point, line, and plane. Furthermore, the approach enables simpler programming, real-time computation, and ability to handle any type of planar and spherical mechanism. This work facilitates creation of practical and intuitive design tools for mechanism designers.

## 1 Introduction

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].^{2} 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].^{3} 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.

## 2 Mechanism Representation and Constraints

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., *ψ*_{1}, *ψ*_{2}, *L*_{3}. 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*_{1}, *y*_{1}, *x*_{2}, *y*_{2}, *x*_{3}, *y*_{3}. 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*_{1}, *y*_{1}, *a*_{2}, *b*_{2} since we can set homogenizing factors *z*_{1} = 1 and *c*_{2} = 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*^{2}) 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*^{2} 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.

Joint | Type | Coordinates |
---|---|---|

J_{1,input} | Revolute | 0, −1 |

J_{2} | Revolute | 1, 0.5 |

J_{3} | Prismatic | −0.17, 0.98, −4.28 |

J_{4} | Revolute | 3.25, 1.4 |

J_{5} | Revolute | 7.72, 1.44 |

J_{6} | Revolute | 11.66, 4.17 |

J_{7} | Prismatic | 0, 1, 1.24 |

J_{8} | Coupler point | 6, −2 |

Joint | Type | Coordinates |
---|---|---|

J_{1,input} | Revolute | 0, −1 |

J_{2} | Revolute | 1, 0.5 |

J_{3} | Prismatic | −0.17, 0.98, −4.28 |

J_{4} | Revolute | 3.25, 1.4 |

J_{5} | Revolute | 7.72, 1.44 |

J_{6} | Revolute | 11.66, 4.17 |

J_{7} | Prismatic | 0, 1, 1.24 |

J_{8} | Coupler point | 6, −2 |

Link | Constituent joints |
---|---|

L_{1} | J_{1}, J_{2} |

L_{2} | J_{2}, J_{3}, J_{4} |

L_{3} | J_{3}, J_{6} |

L_{4} | J_{4}, J_{5}, J_{8} |

L_{5} | J_{5}, J_{6}, J_{7} |

L_{6,ground} | J_{1}, J_{7} |

Link | Constituent joints |
---|---|

L_{1} | J_{1}, J_{2} |

L_{2} | J_{2}, J_{3}, J_{4} |

L_{3} | J_{3}, J_{6} |

L_{4} | J_{4}, J_{5}, J_{8} |

L_{5} | J_{5}, J_{6}, J_{7} |

L_{6,ground} | J_{1}, J_{7} |

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.

Joint | Type | Coordinates |
---|---|---|

J_{1,input} | Revolute | 0.94, 0.24, 0.24 |

J_{2} | Revolute | 0.80, 0.27, 0.53 |

J_{3} | Prismatic | 0.68, −0.68, 0.26 |

J_{4} | Revolute | −0.38, 0.76, 0.53 |

J_{5} | Coupler point | 0.50, −0.21, 0.84 |

Joint | Type | Coordinates |
---|---|---|

J_{1,input} | Revolute | 0.94, 0.24, 0.24 |

J_{2} | Revolute | 0.80, 0.27, 0.53 |

J_{3} | Prismatic | 0.68, −0.68, 0.26 |

J_{4} | Revolute | −0.38, 0.76, 0.53 |

J_{5} | Coupler point | 0.50, −0.21, 0.84 |

Link | Constituent joints |
---|---|

L_{1} | J_{1}, J_{2} |

L_{2} | J_{2}, J_{3}, J_{5} |

L_{3} | J_{3}, J_{4} |

L_{4,ground} | J_{1}, J_{4} |

Link | Constituent joints |
---|---|

L_{1} | J_{1}, J_{2} |

L_{2} | J_{2}, J_{3}, J_{5} |

L_{3} | J_{3}, J_{4} |

L_{4,ground} | J_{1}, J_{4} |

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.

*a*

_{1},

*a*

_{2},

*a*

_{3},

*a*

_{4}) and floating point (

*b*

_{1},

*b*

_{2},

*b*

_{3},

*b*

_{4}), where

*a*

_{4}and

*b*

_{4}are homogenizing factors. The RR link imposes the constraint that the distance between two points remains constant, i.e.,

*dist*(

*R*

_{1},

*R*

_{2}) =

*r*in Figs. 4 and 5. The constraint equation is given as

*a*

_{0}is given as

*r*is the radius of the sphere formed by the RR link with the center given by (

*a*

_{1},

*a*

_{2},

*a*

_{3},

*a*

_{4}). When the

*z*-coordinate is set to zero, the constraint equation degenerates into the one for a planar RR link. The constraint equation for a planar RR link represented by points (

*a*

_{1},

*a*

_{2},

*a*

_{4}) and (

*b*

_{1},

*b*

_{2},

*b*

_{4}) is thus given as

*a*

_{0}is given as

*a*

_{1},

*a*

_{2},

*a*

_{3},

*a*

_{4}) and (

*L*

_{1},

*L*

_{2},

*L*

_{3},

*L*

_{4}), respectively. An RP or PR link imposes the constraint that the distance between a point and a line (planar case) or a point and a plane (spherical case) is constant, i.e.,

*dist*(

*R*,

*P*) =

*d*in Figs. 4 and 5. RP and PR links are inversions of each other and are expressed by the same constraint. The general constraint equation for a spherical RP link is given as

*d*is the signed perpendicular distance between the revolute joint and prismatic joint. For spherical linkages, the prismatic plane always passes through the origin, i.e.,

*L*

_{4}= 0. Thus, the spherical RP link constraint equation reduces to

*z*-coordinates is set to zero, the general constraint equation degenerates into a planar case. Thus, constraint equation for a planar RP or PR link represented by a point (

*a*

_{1},

*a*

_{2},

*a*

_{4}) and a line (

*L*

_{1},

*L*

_{2},

*L*

_{4}) is given as

*d*becomes zero, the equation describes a line passing through a point, i.e., constraint equation of PR or RP links. This is usually the case with the PR links where the moving joint point is constrained to be on the fixed line of the prismatic joint and in case of the RP link where the moving line is constrained to pass through the point of the fixed joint.

*L*

_{1},

*L*

_{2},

*L*

_{3}, 0) and (

*M*

_{1},

*M*

_{2},

*M*

_{3}, 0). For the PP link, the angle between two lines (planar case) or two planes (spherical case) remains constant, i.e.,

*angle*(

*P*

_{1},

*P*

_{2}) = cos

^{−1}(

*k*) in Figs. 4 and 5. The constraint equation is given as

*k*represents the cosine of angle between two prismatic joints. Similarly, for planar PP binary link represented by two lines (

*L*

_{1},

*L*

_{2},

*L*

_{4}) and (

*M*

_{1},

*M*

_{2},

*M*

_{4}), the constraint equation degenerates to

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.

*λ*to prismatic coordinates (

*L*

_{1},

*L*

_{2},

*L*

_{3}) does not change the coordinates. Thus, the magnitude of this vector can be fixed to unity without losing generality and another constraint can be written as

*a*

_{1},

*a*

_{2},

*a*

_{3}) are the coordinates of any revolute joint on a spherical mechanism. Thus, the rigidity constraints described in Eqs. (1), (2), (4)–(9) are sufficient to uniquely determine the unknown coordinates of a

*n*-bar planar or spherical mechanism. This concludes our discussion on representation and constraints for a generalized planar or spherical mechanism.

## 3 Solving Constraint Equations

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.

### 3.1 Input Link Perturbation.

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.

*x*

_{1},

*y*

_{1}) and moving joint (

*x*

_{2},

*y*

_{2}), the new coordinates of moving revolute joint can be given as

*X*

_{2},

*Y*

_{2}) represent the moving joint after perturbation and

*θ*is the angle through which the input link is perturbed.

*x*

_{1},

*y*

_{1}) and moving joint (

*a*,

*b*,

*c*), the new coordinates of moving line representing the prismatic joint can be given as

*A*,

*B*,

*C*) are the moving line coordinates after perturbation, and

**T**and

**R**are the translation and rotation matrices as described in Eq. (11).

_{x}*a*,

*b*,

*c*) and moving joint (

*x*,

*y*), the new coordinates of translating revolute joint can be given as

*X*,

*Y*) are the moving joint coordinates after perturbation and

*d*is the distance through which the prismatic joint is moved along the fixed line. It can be seen that in Eq. (13), the actuating line coordinate

*c*does not affect the new position of moving joint coordinates as new position only depends on direction cosines.

*a*

_{1},

*b*

_{1},

*c*

_{1}) and moving joint (

*a*

_{2},

*b*

_{2},

*c*

_{2}), the new coordinates of translating prismatic joint can be given as

*A*

_{2},

*B*

_{2},

*C*

_{2}) are the moving prismatic joint coordinates after perturbation and

*d*is the distance through which the input link has been perturbed.

*l*,

*m*,

*n*) is given by

**R**

_{x},

**R**

_{y}, and

**R**

_{z}are the rotation matrix around

*x*-,

*y*-, and

*z*-axis and

*θ*is the angle by which the link is rotated around the axis.

*x*

_{1},

*y*

_{1},

*z*

_{1}) are the fixed joint coordinates, (

*x*

_{2},

*y*

_{2},

*z*

_{2}) are the moving joint coordinates before perturbation, and (

*X*

_{2},

*Y*

_{2},

*Z*

_{2}) are the moving joint coordinates after perturbation.

*x*,

*y*,

*z*) are the fixed joint coordinates, (

*a*,

*b*,

*c*) are the moving joint coordinates before perturbation, (

*A*,

*B*,

*C*) are the moving joint coordinates after perturbation, and as described above.

*a*,

*b*,

*c*), its pole coordinates are also given as (

*a*,

*b*,

*c*). Thus, for a spherical RP link with fixed prismatic joint, the coordinates of moving revolute joint can be given as rotation, i.e.,

*a*,

*b*,

*c*) are the fixed prismatic joint coordinates, (

*x*,

*y*,

*z*) are the moving revolute joint coordinates before perturbation and (

*X*,

*Y*,

*Z*) are the moving revolute joint coordinates after perturbation.

*a*

_{1},

*b*

_{1},

*c*

_{1}) are the coordinates of fixed prismatic joint, (

*a*

_{2},

*b*

_{2},

*c*

_{2}) are moving prismatic joint coordinates before perturbation, and (

*A*

_{2},

*B*

_{2},

*C*

_{2}) are moving prismatic joint coordinates after perturbation.

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.

### 3.2 Numerical Non-linear System of Equation Solving.

**q**is the state vector that consists of all the unknown coordinates. The well-known Newton-Raphson method can be used to solve this non-linear system of equation. It is featured by quadratic convergence in the neighborhood of the solution. Since the input link is perturbed by a small finite displacement, the previous state of mechanism serves as a good initial approximation. The number of constraint equations should be equal to or greater than the number of unknowns for this approach to work. For planar and spherical mechanisms, it is always possible to satisfy this criterion using the constraints outlined in section.

**q**

_{i}is the state vector at

*i*th iteration,

**Φ**(

**q**

_{i}) is the vector of residuals at

**q**=

**q**

_{i}, and [

*J*^{−1}(

**q**

_{i})] is the inverse of Jacobian matrix evaluated at

**q**=

**q**

_{i}. The Jacobian matrix is of the following form

*m*is the number of constraints and

*n*is the number of unknown coordinates. Thus to calculate the Jacobian matrix, relations describing the first-order partial derivatives of constraint equations are required.

*a*

_{0}and

*b*

_{4}has been assumed as unity without loss in generality. For planar RR link, only ∂

*C*

_{RR}/∂

*a*

_{1}and ∂

*C*

_{RR}/∂

*a*

_{2}exists.

*a*

_{0}has been assumed as unity without loss in generality. For planar RR constraint given in Eq. (5), the first-order differentials are

*L*

_{3}= 0 and

*M*

_{3}= 0 and ∂

*C*

_{P,PP}/∂

*L*

_{4}= 0.

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

## 4 Examples

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.

### 4.1 Speed Comparison With a Commercial Software.

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.

### 4.2 Planar Stephenson-II Linkage.

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*_{1}, *J*_{7} as the fixed joints, *J*_{2} as the perturbed joint and *J*_{3}, *J*_{4}, *J*_{5}, *J*_{6}, *J*_{8} 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*_{8} as shown in Fig. 2. The run-time of this simulation was 3.79 s.

### 4.3 Planar-Modified Theo Jansen Linkage.

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*_{1}, *J*_{5} are the fixed joints, *J*_{2} is the perturbed joint and the state vector consists coordinates of *J*_{3}, *J*_{4}, *J*_{6}, *J*_{7}, *J*_{8}. 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*_{8} 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.

Joint | Coordinates |
---|---|

J_{1,input} | 2.77, 2.31 |

J_{2} | 2.17, 3.33 |

J_{3} | −0.50, 0.87, −4.80 |

J_{4} | 0.66, −1.3 |

J_{5} | −0.22, 1.72 |

J_{6} | −3.17, 0.66 |

J_{7} | −2.08, −2.24 |

J_{8} | 2.54, −4.64 |

Joint | Coordinates |
---|---|

J_{1,input} | 2.77, 2.31 |

J_{2} | 2.17, 3.33 |

J_{3} | −0.50, 0.87, −4.80 |

J_{4} | 0.66, −1.3 |

J_{5} | −0.22, 1.72 |

J_{6} | −3.17, 0.66 |

J_{7} | −2.08, −2.24 |

J_{8} | 2.54, −4.64 |

Link | Constituent joints |
---|---|

L_{1} | J_{1}, J_{2} |

L_{2} | J_{2}, J_{3} |

L_{3} | J_{2}, J_{4} |

L_{4} | J_{3}, J_{5}, J_{6} |

L_{5} | J_{5}, J_{4} |

L_{6} | J_{6}, J_{7} |

L_{7} | J_{4}, J_{7}, J_{8} |

L_{8,ground} | J_{1}, J_{5} |

Link | Constituent joints |
---|---|

L_{1} | J_{1}, J_{2} |

L_{2} | J_{2}, J_{3} |

L_{3} | J_{2}, J_{4} |

L_{4} | J_{3}, J_{5}, J_{6} |

L_{5} | J_{5}, J_{4} |

L_{6} | J_{6}, J_{7} |

L_{7} | J_{4}, J_{7}, J_{8} |

L_{8,ground} | J_{1}, J_{5} |

### 4.4 Spherical RRPR Mechanism.

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*_{1}, *J*_{4} as the fixed joints, *J*_{2} as the perturbed joint, and *J*_{3}, *J*_{5} 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*_{5} can be plotted as shown in Fig. 3. The run-time of this simulation was 1.49 s.

### 4.5 Spherical Watt-I Linkage.

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*_{1}, *J*_{6} are the fixed joints, *J*_{2}, *J*_{3} are the perturbed joints, and *J*_{4}, *J*_{5}, *J*_{7}, *J*_{8} 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*_{8} as shown in Fig. 7. The run-time of this simulation was 3.33 s.

Joint | Coordinates |
---|---|

J_{1,input} | 0, 0, 1 |

J_{2} | 0.93, 0, 0.37 |

J_{3} | 0.85, −0.17, 0.51 |

J_{4} | 0.70, 0.70, 0.14 |

J_{5} | 0.73, 0.49, 0.49 |

J_{6} | 0.81, 0.41, −0.41 |

J_{7} | 0.48, −0.10, 0.87 |

J_{8} | 0.49, 0.49, 0.73 |

Joint | Coordinates |
---|---|

J_{1,input} | 0, 0, 1 |

J_{2} | 0.93, 0, 0.37 |

J_{3} | 0.85, −0.17, 0.51 |

J_{4} | 0.70, 0.70, 0.14 |

J_{5} | 0.73, 0.49, 0.49 |

J_{6} | 0.81, 0.41, −0.41 |

J_{7} | 0.48, −0.10, 0.87 |

J_{8} | 0.49, 0.49, 0.73 |

Link | Constituent joints |
---|---|

L_{1} | J_{1}, J_{2}, J_{3} |

L_{2} | J_{2}, J_{4}, J_{5} |

L_{3} | J_{4}, J_{6} |

L_{4} | J_{3}, J_{7} |

L_{5} | J_{5}, J_{7}, J_{8} |

L_{6,ground} | J_{1}, J_{6} |

Link | Constituent joints |
---|---|

L_{1} | J_{1}, J_{2}, J_{3} |

L_{2} | J_{2}, J_{4}, J_{5} |

L_{3} | J_{4}, J_{6} |

L_{4} | J_{3}, J_{7} |

L_{5} | J_{5}, J_{7}, J_{8} |

L_{6,ground} | J_{1}, J_{6} |

### 4.6 Spatial 5-SS Mechanism.

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*_{1} and *J*_{7} to actuate the mechanism [30]. This results in the mechanism having *J*_{1}, *J*_{2}, *J*_{3}, *J*_{4}, *J*_{5} as fixed joints and *J*_{6}, *J*_{7}, *J*_{8}, *J*_{9}, *J*_{10}, *J*_{11} 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*_{11}. The input prismatic link is perturbed by 0.01 units and run-time of this simulation was 0.28 s.

Joint | Coordinates |
---|---|

J_{1} | −7.2, −5.29, 7.87 |

J_{2} | 5.72, −0.71, −8.26 |

J_{3} | −8.75, −4.60, 5.49 |

J_{4} | −10.00, −8.72, 6.09 |

J_{5} | −3.88, 6.61, 8.91 |

J_{6} | −6.61, −9.99, 3.93 |

J_{7} | 8.19, −9.05, 6.07 |

J_{8} | 7.84, −5.73, −0.62 |

J_{9} | 4.15, −0.48, −2.04 |

J_{10} | −0.97, −9.04, 1.92 |

J_{11} | 2.20, −7.17, 5.36 |

Joint | Coordinates |
---|---|

J_{1} | −7.2, −5.29, 7.87 |

J_{2} | 5.72, −0.71, −8.26 |

J_{3} | −8.75, −4.60, 5.49 |

J_{4} | −10.00, −8.72, 6.09 |

J_{5} | −3.88, 6.61, 8.91 |

J_{6} | −6.61, −9.99, 3.93 |

J_{7} | 8.19, −9.05, 6.07 |

J_{8} | 7.84, −5.73, −0.62 |

J_{9} | 4.15, −0.48, −2.04 |

J_{10} | −0.97, −9.04, 1.92 |

J_{11} | 2.20, −7.17, 5.36 |

Link | Constituent joints |
---|---|

L_{1} | J_{1}, J_{6} |

L_{2} | J_{2}, J_{7} |

L_{3} | J_{3}, J_{8} |

L_{4} | J_{4}, J_{9} |

L_{5} | J_{5}, J_{10} |

L_{6} | J_{6}, J_{7}, J_{8}, J_{9}, J_{10}, J_{11} |

L_{7,ground} | J_{1}, J_{2}, J_{3}, J_{4}, J_{5} |

Link | Constituent joints |
---|---|

L_{1} | J_{1}, J_{6} |

L_{2} | J_{2}, J_{7} |

L_{3} | J_{3}, J_{8} |

L_{4} | J_{4}, J_{9} |

L_{5} | J_{5}, J_{10} |

L_{6} | J_{6}, J_{7}, J_{8}, J_{9}, J_{10}, J_{11} |

L_{7,ground} | J_{1}, J_{2}, J_{3}, J_{4}, J_{5} |

## 5 Conclusion

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.

## Footnotes

## Acknowledgment

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.