5.1. Position analysis

We start with the position analysis and attach reference frame $Oxyz$ to the base link according to Figure 7. Next, we can consider the vector-loop equation of the mechanism and examine its projections onto the $Oy$ and $Oz$ axes:

(7) \begin{align} \begin{cases} l_1 \cos (\theta _1) + l_2 \cos (\theta _2) - l_4 \cos (\theta _4) = l_5, \\[3pt] l_1 \sin (\theta _1) + l_2 \sin (\theta _2) - l_4 \sin (\theta _4) = 0, \end{cases} \end{align}

where $\theta _2$ is an angle between the screw axis and the $Oy$ axis; $l_1$ , $l_4$ , and $l_5$ are constant lengths of the crank, the rocker, and the base, respectively; $l_2$ is a variable length of the coupler.

In addition, we have the following relation between angles $\theta _2$ and $\theta _4$ (Figure 7):

(8) \begin{align} \theta _4 = \theta _2 + \pi / 2. \end{align}

Using Eq. (8), we can eliminate angle $\theta _2$ and rewrite Eq. (7) as follows:

(9a) \begin{align} l_1 \cos (\theta _1) + l_2 \sin (\theta _4) - l_4 \cos (\theta _4) = l_5, \\[-25pt]\nonumber\end{align}

(9b) \begin{align} l_1 \sin (\theta _1) - l_2 \cos (\theta _4) - l_4 \sin (\theta _4) = 0. \end{align}

Equations (9) includes two unknown parameters: length $l_2$ and angle $\theta _4$ . To find the latter, we multiply both sides of Eq. (9a) by $\cos (\theta _4)$ and both sides of Eq. (9b) by $\sin (\theta _4)$ . Summing up the resulting equations, we get:

(10) \begin{align} (l_1 \cos (\theta _1) - l_5) \cos (\theta _4) + l_1 \sin (\theta _1) \sin (\theta _4) - l_4 = 0. \end{align}

Equation (10) is a trigonometric equation with respect to single variable $\theta _4$ whose solutions have the following form [Reference Waldron, Schmiedeler, Siciliano and Khatib32, sec. 2.7]:

(11) \begin{align} \theta _4 = 2{\textrm{atan}}\left (\dfrac{C_2 \pm \sqrt{C_1^2 + C_2^2 - C_3^2}}{C_1 - C_3}\right ), \end{align}

with coefficients $C_1$ , $C_2$ , and $C_3$ deduced from Eq. (10):

\begin{equation*} C_1 = l_1 \cos (\theta _1) - l_5, \qquad C_2 = l_1 \sin (\theta _1), \qquad C_3 = -l_4. \end{equation*}

The “ $\pm$ ” sign in Eq. (11) corresponds to two different solutions and two different assembly modes of the mechanism. We select a suitable solution according to the mechanism geometry and operating conditions. To find remaining angle $\theta _3$ , we have to compute length $l_2$ first. This time, we multiply both sides of Eq. (9a) by $\sin (\theta _4)$ and both sides of Eq. (9b) by $\cos (\theta _4)$ . Subtracting the resulting equations from each other, we get:

(12) \begin{align} l_2 = l_1 \sin (\theta _1 - \theta _4) + l_5 \sin (\theta _4). \end{align}

Finally, we calculate angle $\theta _3$ using the expression below:

(13) \begin{align} \theta _3 = 2\pi (l_2 - l_2^0) / p, \end{align}

where $l_2^0$ is a value of $l_2$ for $\theta _1 = 0$ , which we can determine using Eqs. (11) and (12). For the considered mechanism (Figure 7), we have $l_2^0 = ((l_5 - l_1)^2 - l_4^2)^{1/2}$ .

We have found all the parameters that describe the configuration of the output link. In the dynamic analysis (Section 6), we will also need the value of angle $\theta _2$ , which we can determine from Eq. (8). This completes the position analysis of the mechanism.

5.2. Velocity analysis

Velocity analysis of the mechanism is the next step. We can find angular speed $\omega _4$ by taking the time derivative of Eq. (10). After rearranging the terms in the resulting equation, we get:

(14) \begin{align} l_2 \omega _4 - l_1 \omega _1 \sin (\theta _1 - \theta _4) = 0, \end{align}

from where we can express angular speed $\omega _4$ .

To find angular speed $\omega _3$ , we differentiate Eq. (13) with respect to time:

(15) \begin{align} \omega _3 = 2\pi \upsilon _2 / p, \end{align}

where $\upsilon _2$ is a time derivative of length $l_2$ , which comes from Eq. (12):

(16) \begin{align} \upsilon _2 = l_1 (\omega _1 - \omega _4) \cos (\theta _1 - \theta _4) + l_5 \omega _4 \cos (\theta _4). \end{align}

During the dynamic analysis, we will need angular speed $\omega _2$ of the coupler. From Eq. (8), we get $\omega _2 = \omega _4$ . This completes the velocity analysis of the mechanism.

5.3. Acceleration analysis

We finish the kinematic analysis of the mechanism by computing accelerations of its links. Differentiating Eq. (14) with respect to time, we obtain angular acceleration $\varepsilon _4$ of the rocker:

(17) \begin{align} \varepsilon _4 = \left(l_1 \left(\varepsilon _1 \sin (\theta _1 - \theta _4) + \omega _1 (\omega _1 - \omega _4) \cos (\theta _1 - \theta _4)\right) - \upsilon _2 \omega _4\right) / l_2. \end{align}

Next, we compute angular acceleration $\varepsilon _3$ by taking a time derivative of Eq. (15):

(18) \begin{align} \varepsilon _3 = 2\pi a_2 / p, \end{align}

where $a_2$ is a time derivative of speed $\upsilon _2$ , which we find from Eq. (16):

(19) \begin{align} a_2 = l_1 \left((\varepsilon _1 - \varepsilon _4) \cos (\theta _1 - \theta _4) - (\omega _1 - \omega _4)^2 \sin (\theta _1 - \theta _4)\right) + l_5 \left(\varepsilon _4 \cos (\theta _4) - \omega _4^2 \sin (\theta _4)\right). \end{align}

Finally, we compute angular acceleration $\varepsilon _2$ of the coupler, which is equal to angular acceleration $\varepsilon _4$ of the rocker, i.e., $\varepsilon _2 = \varepsilon _4$ . Thus, we have expressed positions, angular speeds, and angular accelerations of all the links in terms of the crank motion defined by parameters $\theta _1$ , $\omega _1$ , and $\varepsilon _1$ . This completes the kinematic analysis of the 4R1H mechanism.

6.1. Problem formulation

Kane’s method requires selecting generalized speeds of the mechanism. The considered mechanism has a single DOF, and it is natural to choose angular speed $\omega _1$ of the crank as the generalized speed. Thus, there will be only one equation of motion, which has the following form:

(20) \begin{align} F + F^\star = 0, \end{align}

where $F$ and $F^\star$ are the generalized active and inertia forces, respectively, corresponding to generalized speed $\omega _1$ .

Let $i = 1, \dots, 4$ be an index number of a moving link according to Figure 5. Suppose the $i$ -th link is affected by a set of active forces and torques. We can reduce these forces and torques to resultant force $\mathbf{F}_i$ , whose axis passes through center of mass (COM) $C_i$ , and resultant torque $\mathbf{T}_i$ . Given the above, we can calculate generalized forces $F$ and $F^\star$ in the following way:

(21a) \begin{align} F &= \sum _{i = 1}^{4} \left (\boldsymbol{\unicode{x03C5}}_{\omega 1}^{Ci} \cdot \mathbf{F}_i + \boldsymbol{\unicode{x03C9}}_{\omega 1}^i \cdot \mathbf{T}_i \right ), \\[-27pt] \nonumber \end{align}

(21b) \begin{align} F^\star &= \sum _{i = 1}^{4} \left (\boldsymbol{\unicode{x03C5}}_{\omega 1}^{Ci} \cdot \mathbf{F}_i^\star + \boldsymbol{\unicode{x03C9}}_{\omega 1}^i \cdot \mathbf{T}_i^\star \right ), \end{align}

where $\mathbf{F}_i^\star$ is an inertia force acting on the $i$ -th link in its COM; $\mathbf{T}_i^\star$ is an inertia torque acting on the $i$ -th link; $\boldsymbol{\unicode{x03C5}}_{\omega 1}^{Ci}$ is a partial velocity of point $C_i$ ; $\boldsymbol{\unicode{x03C9}}_{\omega 1}^i$ is a partial angular velocity of the $i$ -th link.Footnote ² We suppose all vector variables are written with respect to the same reference frame.

We compute inertia force $\mathbf{F}_i^\star$ and inertia torque $\mathbf{T}_i^\star$ acting on the $i$ -th link from the following relations:

(22a) \begin{align} \mathbf{F}_i^\star &= -m_i \mathbf{a}_{Ci}, \\[-27pt] \nonumber \end{align}

(22b) \begin{align} \mathbf{T}_i^\star &= -\mathbf{I}_{Ci} \boldsymbol{\unicode{x025B}}_i - \boldsymbol{\unicode{x03C9}}_i \times \mathbf{I}_{Ci}\boldsymbol{\unicode{x03C9}}_i, \end{align}

where $\mathbf{a}_{Ci}$ is an acceleration of the $i$ -th link COM; $\boldsymbol{\unicode{x03C9}}_i$ and $\boldsymbol{\unicode{x025B}}_i$ are angular velocity and acceleration of the $i$ -th link; $m_i$ is a mass of the $i$ -th link; $\mathbf{I}_{Ci}$ is an inertia matrix of the $i$ -th link about its COM.

Partial velocity $\boldsymbol{\unicode{x03C5}}_{\omega 1}^{Ci}$ and partial angular velocity $\boldsymbol{\unicode{x03C9}}_{\omega 1}^i$ are partial derivatives of velocity $\boldsymbol{\unicode{x03C5}}_{Ci}$ of point $C_i$ and angular velocity $\boldsymbol{\unicode{x03C9}}_i$ with respect to generalized speed $\omega _1$ :

(23) \begin{align} \boldsymbol{\unicode{x03C5}}_{\omega 1}^{Ci} = \frac{\partial \boldsymbol{\unicode{x03C5}}_{Ci}}{\partial \omega _1}, \quad \boldsymbol{\unicode{x03C9}}_{\omega 1}^{i} = \frac{\partial \boldsymbol{\unicode{x03C9}}_i}{\partial \omega _1}. \end{align}

To summarize, we should compute the following parameters to derive Eq. (20):

1. 1. Angular velocities and angular accelerations of the links.

2. 2. Accelerations of their COMs.

3. 3. Partial velocities of the COMs and partial angular velocities of the links.

4. 4. Inertia matrices.

5. 5. Active forces and torques.

The next subsection will show that we can compute the parameters listed above with almost no effort using the kinematic expressions from Section 5.

6.2. Kinematic and dynamic parameters

As mentioned in the previous subsection, we consider all the vector parameters are expressed in the same reference frame. Without loss of generality, we select base frame $Oxyz$ as this reference frame. To make the subsequent analysis more convenient, we also attach reference frame $O_ix_iy_iz_i$ to the $i$ -th link, $i = 1, \dots, 4$ , according to Figure 7. Frame $O_3x_3y_3z_3$ of the output link is attached in such a way that the $O_3x_3$ axis is parallel to the $Ox$ axis of the base frame when $\theta _1 = 0$ . Rotation matrices $\mathbf{R}_i$ describe orientation of frames $O_ix_iy_iz_i$ relative to $Oxyz$ , and we calculate these matrices as follows:

(24) \begin{align} \mathbf{R}_i = \mathbf{R}_x(\theta _i), \quad i = 1, 2, 4, \qquad \mathbf{R}_3 = \mathbf{R}_2 \mathbf{R}_y(\theta _3), \end{align}

where $\mathbf{R}_x(\ast )$ and $\mathbf{R}_y(\ast )$ are the matrices of elementary rotations about the corresponding axes [Reference Waldron, Schmiedeler, Siciliano and Khatib32, sec. 2.2]. Angles $\theta _i$ , $i = 2, 3, 4$ , are determined according to Eqs. (8), (11), and (13) we derived in Subsection 5.1. We are now ready to compute all the parameters required for the dynamic analysis.

6.2.1. Angular velocities and angular accelerations

We start with calculating angular velocities $\boldsymbol{\unicode{x03C9}}_i$ , $i = 1, \dots, 4$ . Let $\hat{\mathbf{x}} = [1\,0\,0]^{\mathrm{T}}$ be a unit vector parallel to the $Ox$ axis of frame $Oxyz$ . Using this notation, we can write:

(25) \begin{align} \boldsymbol{\unicode{x03C9}}_i = \omega _i \hat{\mathbf{x}}, \quad i = 1, 2, 4, \qquad \boldsymbol{\unicode{x03C9}}_3 = \boldsymbol{\unicode{x03C9}}_2 + \omega _3 \hat{\mathbf{y}}_2, \end{align}

where $\hat{\mathbf{y}}_2$ is a unit vector equal to the second column of rotation matrix $\mathbf{R}_2$ . Angular speeds $\omega _i$ , $i = 2, 3, 4$ , can be found from Eqs. (14) and (15) given in Subsection 5.2.

Next, we compute angular accelerations $\boldsymbol{\unicode{x025B}}_i$ using Eq. (25):

(26) \begin{align} \boldsymbol{\unicode{x025B}}_i = \varepsilon _i \hat{\mathbf{x}}, \quad i = 1, 2, 4, \qquad \boldsymbol{\unicode{x025B}}_3 = \boldsymbol{\unicode{x025B}}_2 + \varepsilon _3 \hat{\mathbf{y}}_2 + \omega _2 \omega _3 \hat{\mathbf{z}}_2, \end{align}

where $\hat{\mathbf{z}}_2$ is a unit vector equal to the third column of rotation matrix $\mathbf{R}_2$ . To get Eq. (26), we have used condition $\dot{\hat{\mathbf{y}}}_2 = \mathbf{\boldsymbol{\unicode{x03C9}}}_2 \times \hat{\mathbf{y}}_2 = \omega _2 \hat{\mathbf{z}}_2$ . Angular accelerations $\varepsilon _i$ , $i = 2, 3, 4$ , come from Eqs. (17) and (18) derived in Subsection 5.3.

6.2.2. Accelerations of COMs

Let vector $\mathbf{r}_{Ci}$ define the COM coordinates of the $i$ -th link in its reference frame $O_ix_iy_iz_i$ . This is a constant parameter that depends on the link geometry and its mass distribution. Vector $\mathbf{p}_{OiCi} = \mathbf{R}_i \mathbf{r}_{Ci}$ represents these coordinates relative to base frame $Oxyz$ . Given these notations, we compute acceleration $\mathbf{a}_{Ci}$ of the COM as follows:

(27) \begin{align} \mathbf{a}_{Ci} = \mathbf{a}_{Oi} + \boldsymbol{\unicode{x025B}}_i \times \mathbf{p}_{OiCi} + \boldsymbol{\unicode{x03C9}}_i \times (\boldsymbol{\unicode{x03C9}}_i \times \mathbf{p}_{OiCi}), \end{align}

where $\mathbf{a}_{Oi}$ is an acceleration of point $O_i$ :

(28) \begin{align} \begin{aligned} \mathbf{a}_{O1} &= \mathbf{a}_{O4} = \mathbf{0}, \\[3pt] \mathbf{a}_{O2} &= \boldsymbol{\unicode{x025B}}_1 \times \mathbf{p}_{O1O2} + \boldsymbol{\unicode{x03C9}}_1 \times (\boldsymbol{\unicode{x03C9}}_1 \times \mathbf{p}_{O1O2}), \\[3pt] \mathbf{a}_{O3} &= \boldsymbol{\unicode{x025B}}_4 \times \mathbf{p}_{O4O3} + \boldsymbol{\unicode{x03C9}}_4 \times (\boldsymbol{\unicode{x03C9}}_4 \times \mathbf{p}_{O4O3}), \end{aligned} \end{align}

with $\mathbf{p}_{O1O2}$ and $\mathbf{p}_{O4O3}$ being the vectors directed from point $O_1$ to $O_2$ and from point $O_4$ to $O_3$ , respectively. We determine these vectors according to Figure 7:

(29) \begin{align} \mathbf{p}_{O1O2} = \mathbf{R}_1 \left[\begin{array}{l@{\quad}l@{\quad}l} 0 & l_1 & 0 \end{array}\right]^{\mathrm{T}}, \qquad \mathbf{p}_{O4O3} = \mathbf{R}_4 \left[\begin{array}{l@{\quad}l@{\quad}l} 0 & l_4 & 0 \end{array} \right]^{\mathrm{T}}. \end{align}

6.2.3. Partial velocities of COMs and partial angular velocities

We begin with computing partial angular velocities $\boldsymbol{\unicode{x03C9}}_{\omega 1}^i$ , $i = 1, \dots, 4$ . Using Eq. (25) together with Eqs. (14)–(16), we get:

(30) \begin{align} \begin{aligned} \boldsymbol{\unicode{x03C9}}_{\omega 1}^1 &= \hat{\mathbf{x}}, \\[3pt] \boldsymbol{\unicode{x03C9}}_{\omega 1}^2 &= \boldsymbol{\unicode{x03C9}}_{\omega 1}^4 = \dfrac{l_1 \sin (\theta _1 - \theta _4)}{l_2} \hat{\mathbf{x}}, \\[3pt] \boldsymbol{\unicode{x03C9}}_{\omega 1}^3 &= \boldsymbol{\unicode{x03C9}}_{\omega 1}^2 + \dfrac{2\pi l_1}{p} \left (\cos (\theta _1 - \theta _4) - \dfrac{\left(l_1 \cos (\theta _1 - \theta _4) - l_5 \cos (\theta _4)\right) \sin (\theta _1 - \theta _4)}{l_2}\right ) \hat{\mathbf{y}}_2. \end{aligned} \end{align}

To find partial velocity $\boldsymbol{\unicode{x03C5}}_{\omega 1}^{Ci}$ of the $i$ -th link COM, we can first write velocity $\boldsymbol{\unicode{x03C5}}_{Ci}$ of point $C_i$ :

(31) \begin{align} \boldsymbol{\unicode{x03C5}}_{Ci} = \boldsymbol{\unicode{x03C5}}_{Oi} + \boldsymbol{\unicode{x03C9}}_i \times \mathbf{p}_{OiCi}, \end{align}

where $\boldsymbol{\unicode{x03C5}}_{Oi}$ is a velocity of point $O_i$ :

(32) \begin{align} \boldsymbol{\unicode{x03C5}}_{O1} = \boldsymbol{\unicode{x03C5}}_{O4} = \mathbf{0}, \qquad \boldsymbol{\unicode{x03C5}}_{O2} = \boldsymbol{\unicode{x03C9}}_1 \times \mathbf{p}_{O1O2}, \qquad \boldsymbol{\unicode{x03C5}}_{O3} = \boldsymbol{\unicode{x03C9}}_4 \times \mathbf{p}_{O4O3}. \end{align}

The equations above allow us to compute the partial velocities:

(33) \begin{align} \begin{aligned} \boldsymbol{\unicode{x03C5}}_{\omega 1}^{C1} &= \boldsymbol{\unicode{x03C9}}_{\omega 1}^1 \times \mathbf{p}_{O1C1}, \\[3pt] \boldsymbol{\unicode{x03C5}}_{\omega 1}^{C2} &= \boldsymbol{\unicode{x03C9}}_{\omega 1}^1 \times \mathbf{p}_{O1O2} + \boldsymbol{\unicode{x03C9}}_{\omega 1}^2 \times \mathbf{p}_{O2C2}, \\[3pt] \boldsymbol{\unicode{x03C5}}_{\omega 1}^{C3} &= \boldsymbol{\unicode{x03C9}}_{\omega 1}^4 \times \mathbf{p}_{O4O3} + \boldsymbol{\unicode{x03C9}}_{\omega 1}^3 \times \mathbf{p}_{O3C3}, \\[3pt] \boldsymbol{\unicode{x03C5}}_{\omega 1}^{C4} &= \boldsymbol{\unicode{x03C9}}_{\omega 1}^4 \times \mathbf{p}_{O4C4}. \end{aligned} \end{align}

6.2.4. Inertia matrices

Let matrix $\mathbf{I}_{Ci}'$ be an inertia matrix of the $i$ -th link computed in its COM with respect to the reference frame whose axes are parallel to the corresponding axes of link frame $O_ix_iy_iz_i$ . This matrix is a constant parameter that depends on link mass distribution of the link. Now, we can compute variable inertia matrix $\mathbf{I}_{Ci}$ , which depends on the link configuration and appears in Eq. (22b), as follows [Reference Roithmayr and Hodges34, sec. 4.5]:

(34) \begin{align} \mathbf{I}_{Ci} = \mathbf{R}_i \mathbf{I}_{Ci}' \mathbf{R}_i^{\mathrm{T}}, \quad i = 1, \dots, 4. \end{align}

6.2.5. Active forces and torques

Concerning the assumptions given at the beginning of this section, we suppose the only active forces are the weights of the mechanism links. There is only one active torque—torque $\tau$ generated by the motor. We ignore other active forces and torques (such as payload) because they will depend on a specific application, but we can always include them in our analysis using the same techniques. Thus, we obtain the following expressions of resultant forces $\mathbf{F}_i$ and torques $\mathbf{T}_i$ :

(35) \begin{align} \begin{aligned} \mathbf{F}_i &= m_i \mathbf{g}, & i &= 1, \dots, 4, \\[3pt] \mathbf{T}_i &= \mathbf{0}, & i &= 2, 3, 4, \qquad \mathbf{T}_1 = \tau \hat{\mathbf{x}}, \end{aligned} \end{align}

where $\mathbf{g}$ is a vector of gravitational acceleration. Thus, we have found all the kinematic and dynamic parameters required to compose the motion equation of the mechanism.

6.3. Inverse dynamics

The subsequent dynamic analysis goes as follows. For each link, we compute inertia force $\mathbf{F}_i^\star$ and inertia torque $\mathbf{T}_i^\star$ using Eq. (22). Next, we determine generalized inertia force $F^\star$ according to Eq. (21b). Let $F'$ denote the first term of generalized active force $F$ in Eq. (21a) that considers the effect of active forces $\mathbf{F}_i$ . Given this notation and concerning Eqs. (30) and (35), we can rewrite generalized active force $F$ :

(36) \begin{align} F = F' + \sum _{i = 1}^{4} \left (\boldsymbol{\unicode{x03C9}}_{\omega 1}^i \cdot \mathbf{T}_i\right ) = F' + \boldsymbol{\unicode{x03C9}}_{\omega 1}^1 \cdot \mathbf{T}_1 = F' + \tau . \end{align}

Substituting the results into Eq. (20), we obtain the solution to the inverse dynamic problem:

(37) \begin{align} \tau = -F' - F^\star . \end{align}

This concludes the dynamic analysis of the 4R1H mechanism.

7.1. Virtual prototyping

Numerical simulations of the developed 4R1H mechanism require specifying the values of its geometrical and inertial parameters. For this purpose, we have developed a CAD model (virtual prototype) of the mechanism using the Autodesk Inventor software (Figure 8). This model allows us to set the material of the links and accurately determine their lengths, locations of COMs, masses, and inertia matrices. The virtual prototype can also be used to perform kinematic and dynamic simulations in the CAD software and validate the proposed theoretical formulations. In addition, we can specify friction coefficients and carry out dynamic simulations with friction effects, which are challenging to consider in analytical computations. Subsection 7.2 will show that friction can significantly affect the computed motor torque.

Figure 8.

Virtual prototype (CAD model) of the 4R1H mechanism.

In the prototype, the rotation axes of the crank and the rocker do not lie in the horizontal plane (Figure 8). This means the $Oy$ axis of base frame $Oxyz$ is not horizontal either. This axis is inclined relative to the horizontal line for angle $\alpha ={\textrm{atan2}}(l_5^z, l_5^y)$ , where $l_5^y$ and $l_5^z$ are the horizontal and vertical displacements of the rocker axis relative to the crank axis. Such a design affects the computation of forces $\mathbf{F}_i$ in Eq. (35), where we get:

\begin{equation*} \mathbf {g} = g {\begin {bmatrix} 0 & -\sin\! (\alpha ) & -\cos\! (\alpha ) \end {bmatrix}}^{\mathrm {T}}, \end{equation*}

with $g = 9.81\text{ m/s$^2$}$ .

According to the CAD model, the mechanism has the following geometrical parameters: $p = 25\text{ mm}$ , $l_1 = 20\text{ mm}$ , $l_4 = 40\text{ mm}$ , $l_5^y = 160\text{ mm}$ , $l_5^z = -20\text{ mm}$ , and we compute $l_5 = ((l_5^y)^2 + (l_5^z)^2)^{1/2}$ . Table I lists the inertial parameters of the links. Inertia matrices $\mathbf{I}_{Ci}'$ , $i = 1, \dots, 4$ , are computed with respect to the $C_ix_{Ci}y_{Ci}z_{Ci}$ reference frames located at the COMs according to Figure 9. Thus, we have determined all kinematic and dynamic parameters, which are necessary for numerical simulations.

Table I.

Inertial parameters of the mechanism links.

Figure 9.

Mechanism links and their reference frames: (a) crank; (b) coupler; (c) nut; (d) rocker. Blue frames correspond to Figure 7. Red frames have the same orientation as the blue ones, and they are located at the COMs of the links. The COMs are determined using the virtual prototype (Figure 8).

7.2. Numerical simulations

This subsection will look at two examples of the kinematic and dynamic analysis of the 4R1H mechanism. In the first example, we suppose the crank rotates at a constant speed and consider three different values of this speed: $5$ rpm, $10$ rpm, and $20$ rpm. To perform the kinematic and dynamic analysis, we have implemented the algorithms from Sections 5 and 6 in Julia.Footnote ³ Figure 10 illustrates the results of the kinematic analysis for one revolution of the crank. As expected, the shape of the plots for angles $\theta _3$ and $\theta _4$ does not depend on crank speed $\omega _1$ : the plots coincide for all three considered values of the crank speed. The absolute values of angular speeds $\omega _3$ and $\omega _4$ increase linearly as the crank speed increases: this behavior matches Eqs. (14)–(16). On the other hand, the absolute values of angular accelerations $\varepsilon _3$ and $\varepsilon _4$ grow more rapidly because they are quadratic functions of crank speed $\omega _1$ , as we see from Eqs. (17)–(19). We have also performed similar simulations in Autodesk Inventor and got identical results, which we omit here for clarity.

Figure 10.

Results of the kinematic analysis for one revolution of the crank at various angular speeds.

The results show angle $\theta _3$ increases and reaches the highest value of about $600^\circ$ when the crank makes a half-turn. This agrees with the mobility analysis we did in Section 3: the output link changes the direction of its rotation about the screw axis when the crank becomes aligned with the base link. Angular speeds $\omega _3$ and $\omega _4$ vary smoothly—the mechanism does not meet the singular configuration discussed at the end of Section 4.

Table II.

Comparison of the torque computation results between the analytical and CAD methods.

Figure 11.

Computed motor torque for one revolution of the crank at various angular speeds.

Figure 11 shows the solution to the inverse dynamic problem computed with various methods for three different values of the crank speed. First, we applied the analytical method proposed in Section 6 (the red lines in Figure 11). Table II also presents maximum torque $\tau _{\text{max}}$ for each considered value of $\omega _1$ . The motor torque does not change much when the crank speed increases, and its highest value is about $7.5$ N $\cdot$ mm. Even if the crank speed increases four times, the torque increases by only $1$ %. These results indicate that the active forces (the weight of the links) dominate over the inertia forces, although the values of the latter should increase quadratically according to Figure 10 and the previous discussion. Next, we simulated the same motion in the Autodesk Inventor software using the developed virtual prototype (the green lines in Figure 11). Table II lists maximum torque $\tau _{\text{max}}^{\text{CAD}}$ , maximum torque deviation $\Delta _{\text{max}}$ between the two methods along the whole rotation, and the root-mean-square deviation for each value of $\omega _1$ . We see the torque deviations between the analytical and CAD methods increase with the increase in the crank speed. This behavior can be explained as follows. Kane’s method, implemented in the analytical approach, allows computing the motor torque without the need to determine the joint reactions. On the other hand, dynamic simulations in the CAD software always determine the motor torque together with these reactions. Since the considered 4R1H mechanism is overconstrained, there is an infinite number of solutions to the inverse dynamic problem. The CAD software finds the unique solution numerically, for example, using the minimum norm least squares method [Reference Callejo, Gholami, Enzenhöfer and Kövecses39], and the computed torque will differ from the real one. Thus, the torque calculations performed by the analytical method are more accurate.

Finally, we performed CAD simulations considering friction effects. We adopted a simple Coulomb model [Reference Farhat, Mata, Page and Díaz-Rodríguez40] with the following friction coefficients: $0.002$ for the bearings (in the base–crank and rocker–nut couplings) and $0.460$ for the joints without bearings (in other couplings). For the latter, we considered a plastic–plastic contact with no lubrication. The blue lines in Figure 11 show the peak torque is about $14$ N $\cdot$ mm for each considered value of $\omega _1$ , which is $87$ % higher than in the frictionless simulations. We can also observe a slight discontinuity in the computed torque when $\theta _1 = 180^\circ$ . In this configuration, the output link changes its direction of rotation about the screw axis (Figure 10), and friction forces in the H pair change their direction as well. Thus, the simulations show the friction effects can significantly affect the results. At the same time, there exists the same problem of solving the inverse dynamic problem for the overconstrained mechanism, as we considered in the previous paragraph, so the obtained results should be treated with caution. In Section 8, we will also discuss why it is difficult to consider friction within the analytical dynamic model.

In the second example, we add acceleration and deceleration phases to the crank rotation and consider three revolutions of the crank. During the first revolution, crank angular speed $\omega _1$ increases linearly from $0$ to $5$ rpm. For the second turn, the crank rotates at a constant speed. In the third and final revolution, the crank angular speed decreases linearly to zero.

Figure 12.

Results of the kinematic and dynamic analysis for three revolutions of the crank with acceleration and deceleration phases and a middle phase with constant angular speed $\omega _1 = 5$ rpm. The dashed lines separate acceleration, constant speed, and deceleration phases of the motion. Results of the kinematic analysis coincide for the analytical and CAD computations, and the latter are omitted for clarity.

Figure 13.

The actuated full-scale physical prototype of the 4R1H mechanism.

Figure 12 shows the results of the kinematic and dynamic analysis. The plots for angles $\theta _3$ and $\theta _4$ have three identical portions, which coincide with Figure 10. The amplitudes of angular speeds $\omega _3$ and $\omega _4$ and angular accelerations $\varepsilon _3$ and $\varepsilon _4$ increase during the acceleration phase, than they have the same shape as in Figure 10. After that, these amplitudes decrease as the crank decelerates. Unlike the first example, angular speed $\omega _4$ starts and finishes with a zero value because of the added transient modes. On the other hand, angular acceleration $\varepsilon _4$ still starts and finishes with a nonzero value and also has a slight discontinuity after the second and third revolutions of the crank. This is caused by the discontinuity of crank angular acceleration $\varepsilon _1$ . In practice, this behavior can be undesirable, and one can consider more complex transient regimes (e.g., using high-order polynomials or trigonometric functions [Reference Biagiotti and Melchiorri41]) to avoid any acceleration discontinuities.

Figure 12 also indicates that the motor torque behaves similarly in each of the three motion phases, and the shape of the torque curves match Figure 11. This result agrees with the first example where we found that the inertia effects did not contribute much to the motor torque. Finally, numerical simulations in Autodesk Inventor verify the analytical computations in both the kinematic and dynamic analysis. The small torque deviations between the CAD and analytical results are caused by the same reasons, as we discussed in the first example.Footnote ⁴

7.3. Physical prototyping

Using the CAD model, we have designed a physical prototype of the mechanism (Figure 13). The thin rod attached to the nut serves for visualization purposes and helps distinguishing the nut rotations. Dimensions of the links correspond to Subsection 7.1. Most prototype components are 3D-printed plastic parts made of polylactic acid (PLA). Other components represent easy-to-access standard elements, including aluminum profiles, ball bearings (in the base–crank and rocker–nut couplings), and fasteners.

We have designed the prototype mainly for illustrative purposes: it is supposed to operate at a crank speed of $5$ rpm without any payload. Based on the numerical simulations discussed in the previous subsection, we have selected a $12$ V ”Zheng“ ZGA37RG DC motor with a $1\,:\,627$ gearbox to actuate the mechanism. The output rated speed and rated torque of the motor are $5$ rpm and $1$ N $\cdot$ m, respectively, which matches the desired motion.

Figure 14 shows the physical prototype during operation for one revolution of the crank performed at a crank speed of $5$ rpm. The observed motion corresponds to the kinematic simulations. The mechanism does not meet any singularities during the motion cycle and has no jamming in the helical pair, which verifies its performance ability.

Figure 14.

The physical prototype during operation for one revolution of the crank.