3.1. Kinematics positive solution calculation

Based on the D-H method, the coordinate system of the five-bar mechanism of the man-machine integration system is established. The thigh is defined as the link 1, and the other parts are recursively recursive, as shown in Fig. 3. The hip joint flexion as $\theta _1$ and it is between the link i and the horizontal direction, the knee joint extension as $\theta _2$ and the ankle joint extension as $\theta _3$ , and the foot handle rotation extension as $\theta _4$ thanks to the angle between the link i-1 and the link i extension line. The angle between the rehabilitation prototype shaft and the foot is defined as $\theta _4$ called the angle of the foot handle, and the angle between the prototype shaft and the frame is defined as $\theta _5$ called the rotation angle of the prototype. $\alpha _{i-1}$ is the angle between the link i-1 and the link i along the z-axis direction x-axis, $d_i$ is the distance between the link i-1 and the link i along the z-axis direction. The D-H parameters of Table II are obtained and the forward and inverse kinematics analysis is carried out.

According to the parameters in Table II, the transformation matrices ${ }_{1}^{0}\mathrm{T}$ , ${ }_{2}^{1}\mathrm{T}$ , ${ }_{3}^{2}\mathrm{T}$ , ${ }_{4}^{3}\mathrm{T}$ , ${ }_{5}^{4}\mathrm{T}$ of each coordinate system are obtained, and they are substituted into the terminal pose matrix to obtain the positive kinematics solution as follows.

(1) \begin{equation}{ }_{5}^{0} \mathrm{\,T}={ }_{1}^{0} \mathrm{\,T} \cdot{ }_{2}^{1} \mathrm{\,T} \cdot{ }_{3}^{2} \mathrm{\,T} \cdot{ }_{4}^{3} \mathrm{\,T} \cdot{ }_{5}^{4} \mathrm{\,T}=\left [\begin{array}{c@{\quad}c@{\quad}c@{\quad}c}n_{x} & o_{x} & a_{x} & p_{x} \\[5pt] n_{y} & o_{y} & a_{y} & p_{y} \\[5pt] n_{z} & o_{z} & a_{z} & p_{z} \\[5pt] 0 & 0 & 0 & 1\end{array}\right ] \end{equation}

In the Eq. (1)

\begin{align*} n_x&=\cos \!\left (\theta _1+\theta _2+\theta _3+\theta _4+\theta _5\right ) \nonumber\\[5pt] n_y&=\sin \!\left (\theta _1+\theta _2+\theta _3+\theta _4+\theta _5\right ) \nonumber\\[5pt] n_z&=0 \nonumber\\[5pt] o_x&=-\sin \!\left (\theta _1+\theta _2+\theta _3+\theta _4+\theta _5\right ) \nonumber\\[5pt] o_y&=\cos \!\left (\theta _1+\theta _2+\theta _3+\theta _4+\theta _5\right ) \nonumber\\[5pt] o_z&=0 \nonumber \\[5pt] a_x&=0 \nonumber \end{align*}

Table I. P95 percentile adult male corrected size.

Figure 2. Man-machine integration model of active training in lower limb.

Figure 3. D-H coordinate system of the man-machine integrated system for lower limb active training.

(2) \begin{align} a_y&=0 \nonumber\\[5pt] a_z&=1 \nonumber\\[5pt] p_x&=a_4 \cos \!\left (\theta _1+\theta _2+\theta _3+\theta _3\right .=1 \nonumber\\[5pt] p_y&=a_4 \sin \!\left (\theta _1+\theta _2+\theta _3+\theta _4\right )+a_2 \cos \!\left (\theta _1+\theta _2\right )+a_1 \sin \!\left (\theta _1+\theta _2\right )+a_1 \sin \theta _1+a_3 \sin \!\left (\theta _1+\theta _2+\theta _3\right ) \nonumber\\[5pt] p_z&=d_1+d_2+d_3+d_4+d_5 \end{align}

Through the inverse kinematics calculation of the lower limb active training based on the five-bar mechanism of the man-machine integration system, the expression of each joint angle is obtained, which lays a theoretical foundation for the human-computer interaction of the lower limb active rehabilitation training.

3.2. Inverse kinematics solution calculation

Inverse kinematics based on pose parameters is known. Multiplied by ${ }_{1}^{0}\mathrm{T}^{-1}$ both sides of the Eq. (1) to the left respectively, that is ${ }_{1}^{0}\mathrm{T}^{-1}{ }_{5}^{0}\mathrm{T}={ }_{2}^{1}\mathrm{T}{ }_{3}^{2}\mathrm{T}{ }_{4}^{3}\mathrm{T}{ }_{5}^{4}\mathrm{T}$ . The equality [1, 1], [1, 2], [2, 1] and [2, 2] of the matrix elements are equal on both sides respectively, which can get

(3) \begin{equation} \left \{\begin{array}{l} n_x \cos \theta _1+n_y \sin \theta _1=o_y \cos \theta _1-o_x \sin \theta _1 \\[5pt] n_y \cos \theta _1-n_x \sin \theta _1=-o_x \cos \theta _1-o_y \sin \theta _1 \end{array}\right. \end{equation}

From Eq. (3)

(4) \begin{equation} \theta _1=\arctan 2\!\left (n_y+o_x, n_x-o_y\right ) \end{equation}

Multiplied by ${ }_{2}^{1}\mathrm{T}^{-1}{ }_{1}^{0}\mathrm{T}^{-1}$ on both sides of the Equation sign of Eq. (1) to the left at the same time and take [1, 3] and [2, 3] to be Equation respectively to get

(5) \begin{equation} \theta _2=\arctan 2\!\left (a_y, a_x\right )-\arctan 2\!\left (n_x-o_y,-o_x-n_y\right ) \end{equation}

Table II. D-H parameter table of the man-machine integrated system for lower limb active training.

Taking [1, 4] and [2, 4] in Eq. (1) to corresponding Equation can be obtained

(6) \begin{equation} 2 k_1 a_3 \cos \!\left (\theta _1+\theta _2+\theta _3\right )+2 k_2 a_3 \sin \!\left (\theta _1+\theta _2+\theta _3\right )=k_1^2+k_2^2+a_3^2-a_4^2 \end{equation}

In the Eq. (6)

(7) \begin{equation} \left \{\begin{array}{l} k_1=p_x-a_2 \cos \!\left (\theta _1+\theta _2\right )-a_1 \cos \theta _1 \\[5pt] k_2=p_y-a_2 \sin \!\left (\theta _1+\theta _2\right )-a_1 \sin \theta _1^* \end{array}\right. \end{equation}

(8) \begin{equation} \left \{\begin{array}{l} a_4 \cos \!\left (\theta _1+\theta _2+\theta _3+\theta _4\right )=p_x-a_2 \cos \!\left (\theta _1+\theta _2\right )-a_1 \cos \theta _1-a_3 \cos \!\left (\theta _1+\theta _2+\theta _3\right ) \\[5pt] a_4 \sin \!\left (\theta _1+\theta _2+\theta _3+\theta _4\right )=p_y-a_2 \sin \!\left (\theta _1+\theta _2\right )-a_1 \sin \theta _1-a_3 \sin \!\left (\theta _1+\theta _2+\theta _3\right ) \end{array}\right. \end{equation}

By the simultaneous solution of Eq. (6) and Eq. (8)

(9) \begin{equation} \theta _3=\arctan 2\!\left (\mathrm{\,A} \pm \sqrt{\mathrm{A}^2+\mathrm{B}^2-\mathrm{C}^2}, \mathrm{B}-\mathrm{C}\right )-\arctan 2\!\left (a_y, a_x\right ) \end{equation}

In the Eq. (9)

(10) \begin{equation} \left \{\begin{array}{c} \mathrm{A}=2 k_2 a_3 \\[5pt] \mathrm{\,B}=2 k_1 a_3 \\[5pt] \mathrm{C}=k_1^2+k_2^2+a_3^2-a_4^2 \end{array}\right. \end{equation}

(11) \begin{equation} \begin{aligned} \theta _4= & \arctan 2\!\left (\dot{p}_y-\dot{a}_2 \sin \!\left (\theta _1+\theta _2\right )-\dot{a}_1 \sin \theta _1-\dot{a}_3 \sin \!\left (\theta _1+\theta _2+\theta _3\right ), \dot{p}_x-\dot{a}_2 \cos \!\left (\theta _1+\theta _2\right )-\dot{a}_1 c \theta _1-\right. \\[5pt] & -\dot{a}_3 \cos \!\left (\theta _1+\theta _2+\theta _3\right )-\arctan 2\!\left (\mathrm{\,A} \pm \sqrt{\mathrm{A}^2+\mathrm{B}^2-\mathrm{C}^2}, \mathrm{\,B}-\mathrm{C}\right ) \end{aligned} \end{equation}

Taking [1, 1], [1, 2], [2, 1] and [2, 2] in Eq. (1) to be Equation respectively, then

(12) \begin{equation} \theta _3=\arctan 2\!\left (\mathrm{\,A} \pm \sqrt{\mathrm{A}^2+\mathrm{B}^2-\mathrm{C}^2}, \mathrm{B}-\mathrm{C}\right )-\arctan 2\!\left (a_y, a_x\right ) \end{equation}

The inverse kinematics solution of active training of lower limb based on the five-bar mechanism of the man-machine integration system is calculated and the expression of each joint angle is obtained, which lays a theoretical foundation for human-computer interaction in lower limb active rehabilitation training.

4.1. Establishment the dynamics model based on Lagrange method

(1) Establish Lagrange coordinate system

The dynamic coordinate system and simulation of three-dimensional model of man-machine system under active training of lower limbs is established, as shown in Fig. 4, the base coordinate system is established at the hip joint position. The thigh, calf, foot and prototype shaft are named as link 1, link 2, link 3 and link 4, respectively. The connection between the hip joint and the base coordinate system is link 5, and the lengths of links 1 $\sim$ 5 are $a_1,a_2,a_3$ and $a_4$ .

Figure 4. Dynamic coordinate system and simulation of three-dimensional model.

(2) Solving the system’s kinetic energy

The dynamic Equation of the system is established on the basis of Lagrange function. The Lagrange function is obtained by the differences between the kinetic energy and potential energy of the system, and the kinetic energy and potential energy of the system are superimposed by the kinetic energy and potential energy of each link. Therefore, according to the geometric size relationship of each link in Fig. 3, the coordinates of each particle are described, as shown in the Eq. (13).

(13) \begin{equation} \left \{\begin{array}{l} x_1=a_1 \cos \theta _1/ 2 \\[5pt] y_1=a_1 \sin \theta _1/ 2, \,z_1=0 \\[5pt] x_2=a_1 \cos \theta _1+a_2 \sin \!\left (\theta _1-\theta _2\right )/ 2 \\[5pt] y_2=a_1 \sin \theta _1-a_2 \cos \!\left (\theta _1-\theta _2\right )/ 2, \,z_2=0 \\[5pt] x_3=a_1 \cos \theta _1+a_2 \sin \!\left (\theta _1-\theta _2\right )+a_3 \cos \!\left (\theta _1+\theta _3-\theta _2\right )/ 2 \\[5pt] y_3=a_1 \sin \theta _1-a_2 \cos \!\left (\theta _1-\theta _2\right )-a_3 \sin \!\left (\theta _1+\theta _3-\theta _2\right )/ 2, \,z_3=0 \\[5pt] x_4=a_1 \cos \theta _1+a_2 \sin \!\left (\theta _1-\theta _2\right )+a_3 \cos \!\left (\theta _1+\theta _3-\theta _2\right )+a_4 \cos \!\left (\theta _1+\theta _3-\theta _2-\theta _4\right )/ 2 \\[5pt] y_4=a_1 \sin \theta _1-a_2 \cos \!\left (\theta _1-\theta _2\right )-a_3 \sin \!\left (\theta _1+\theta _3-\theta _2\right )-a_4 \sin \!\left (\theta _1+\theta _3-\theta _2-\theta _4\right )/ 2, \,z_4=0 \end{array}\right. \end{equation}

In the Eq. (13), $x$ and $y$ are the coordinates of the center of mass, all of $z$ are 0.

The velocity multiplier can be obtained by adding the multipliers of the first derivative of the $x$ , $y$ and $z$ coordinates of the particle. After obtaining the $v^2$ of the four particles, the kinetic energy of the link is solved by Eq. (16). Finally, the total kinetic energy of the system can be obtained by superimposing the kinetic energy of four particles, as shown in Eq. (17) and Eq. (18).

(14) \begin{equation} \begin{array}{l} v_i^2=\dot{x}_1^2+\dot{y}_1^2+\dot{z}_1^2 \\[5pt] \end{array} \end{equation}

(15) \begin{equation} \left \{\begin{aligned} v_1^2= & \dot{a}_1^2/ 2 \\[5pt] v_2^2= & 2 \dot{a}_1^2+\dot{a}_2^2/ 2-\dot{a}_1 \dot{a}_2 \sin \theta _2 \\[5pt] v_3^2= & \dot{a}_1^2+\dot{a}_2^2+\dot{a}_3^2/ 4-2 \dot{a}_1 \dot{a}_2 \sin \theta _2+\dot{a}_1 \dot{a}_3 \cos \!\left (2 \theta _1+\theta _3-\theta _2\right )+\dot{a}_2 \dot{a}_3 \sin \!\left (2 \theta _1+\theta _3-2 \theta _2\right ) \\[5pt] v_4^2= & \dot{a}_1^2+\dot{a}_2^2+\dot{a}_3^2+\dot{a}_4^2/ 4-2 \dot{a}_1 \dot{a}_2 \sin \theta _2+2 \dot{a}_1 \dot{a}_3 \cos \!\left (2 \theta _1+\theta _3-\theta _2\right )+2 \dot{a}_2 \dot{a}_3 \sin \!\left (2 \theta _1+\theta _3-2 \theta _2\right )+ \\[5pt] & +\dot{a}_1 \dot{a}_4 \cos \!\left (2 \theta _1+\theta _3-\theta _2-\theta _4\right )+\dot{a}_3 \dot{a}_4 \cos \theta _4+\dot{a}_2 \dot{a}_4 \cos \!\left (2 \theta _1+\theta _3-2 \theta _2-\theta _4\right ) \end{aligned}\right. \end{equation}

(16) \begin{equation} \begin{array}{l} E_{k i}=m_i v_i^2/ 2 \\[5pt] \end{array} \end{equation}

(17) \begin{equation} \left \{\begin{array}{l} \begin{aligned} E_{k 1}=m_1 v_1^2/ 2= & m_1 \dot{a}_1^2/ 4 \\[5pt] E_{k 2}=m_2 v_2^2/ 2= & m_2\!\left (2 \dot{a}_1^2+\dot{a}_2^2/ 2-\dot{a}_1 \dot{a}_2 \sin \theta _2\right )/ 2 \\[5pt] E_{k 3}=m_3 v_3^2/ 2= & m_2\!\left (\dot{a}_1^2+\dot{a}_2^2+\dot{a}_3^2/ 4-2 \dot{a}_1 \dot{a}_2 \sin \theta _2+\dot{a}_1 \dot{a}_3 \cos \!\left (2 \theta _1+\theta _3-\theta _2\right )+\right. \\[5pt] & +\left .\dot{a}_2 \dot{a}_3 \sin \!\left (2 \theta _1+\theta _3-2 \theta _2\right )\right )/ 2 \\[5pt] E_{k 4}=m_4 v_4^2/ 2= & m_4\!\left (\dot{a}_1^2+\dot{a}_2^2+\dot{a}_3^2+\dot{a}_4^2/ 4-2 \dot{a}_1 \dot{a}_2 \sin \theta _2+2 \dot{a}_1 \dot{a}_3 \cos \!\left (2 \theta _1+\theta _3-\theta _2\right )+\right. \\[5pt] & +2 \dot{a}_2 \dot{a}_3 \sin \!\left (2 \theta _1+\theta _3-2 \theta _2\right )+\dot{a}_1 \dot{a}_4 \cos \!\left (2 \theta _1+\theta _3-\theta _2-\theta _4\right )+\dot{a}_3 \dot{a}_4 \cos \theta _4+ \\[5pt] & +\left .\dot{a}_2 \dot{a}_4 \cos \!\left (2 \theta _1+\theta _3-2 \theta _2-\theta _4\right )\right )/ 2 \end{aligned} \end{array}\right. \end{equation}

(18) \begin{equation} \begin{array}{l} E_k=\sum _1^4 E_{k i}=E_{k 1}+E_{k 2}+E_{k 3}+E_{k 4} \\[5pt] \end{array} \end{equation}

where $v_i^2\!\left (i=1,2,3,4\right )$ are the velocity of the center of mass, $E_{ki}^2\!\left (i=1,2,3,4\right )$ are the kinetic energy of the center of mass, $m_i^2\!\left (i=1,2,3,4\right )$ are the mass of each link, $E_k$ is the total kinetic energy.

(3) Solving the system potential energy

The potential energy of each particle of the man-machine system is solved by the Eq. (19), and the system potential energy is further superimposed, such as the Eq. (20) and Eq. (21).

(19) \begin{equation} \begin{array}{l} E_{p i}=m_i g y_i \\[5pt] \end{array} \end{equation}

(20) \begin{equation} \left \{\begin{aligned} E_{p i}= & m_1 g y_1=m_1 g a_1 \sin \theta _1/ 2 \\[5pt] E_{p 2}= & m_2 g y_2=m_2 g\!\left (a_1 \sin \theta _1-a_2 \cos \!\left (\theta _1-\theta _2\right )/ 2\right ) \\[5pt] E_{p 3}= & m_3 g y_3=m_3 g\!\left (a_1 \sin \theta _1-a_2 \cos \!\left (\theta _1-\theta _2\right )-a_3 \sin \!\left (\theta _1+\theta _3-\theta _2\right )/ 2\right ) \\[5pt] E_{p 4}= & m_4 g y_4=m_4 g\!\left (a_1 \sin \theta _1-a_2 \cos \!\left (\theta _1-\theta _2\right )-a_3 \sin \!\left (\theta _1+\theta _3-\theta _2\right )-\right. \\[5pt] & -a_4 \sin \!\left (\theta _1+\theta _3-\theta _2-\theta _4\right )/ 2 \end{aligned}\right. \end{equation}

(21) \begin{equation} \begin{array}{l} E_p=\sum _1^4 E_{p i}=E_{p 1}+E_{p 2}+E_{p 3}+E_{p 4} \\[5pt] \end{array} \end{equation}

where $E_{pi}^2\!\left (i=1,2,3,4\right )$ are the potential energy of the center of mass, $E_p$ is the total potential energy.

(4) Establishment of kinetic Equation

The Lagrange function is obtained by the differences between the kinetic energy and the potential energy of the system. The Lagrange function is calculated as shown in the Eq. (25), and the dynamic Equation of the center of mass of the foot can be obtained. It can be observed in the form of matrix Equation, which is composed of inertial force term and gravity term.

(22) \begin{equation} \begin{array}{l} L=E_k-E_p \\[5pt] \end{array} \end{equation}

(23) \begin{equation} F_i=\frac{d}{d t} \frac{\partial L}{\dot{a}_i}-\frac{\partial L}{\mathrm{a}_i} \end{equation}

(24) \begin{equation} \left \{\begin{aligned} F_1= & \left (m_1/ 2+m_2\!\left (2+\dot{a}_2 \sin \theta _2/ 2\right )+m_3\!\left (1+\dot{a}_2 \sin \theta _2+\dot{a}_3 \cos \!\left (2 \theta _1+\theta _3-\theta _2\right )/ 2\right )+\right. \\[5pt] & +\left .m_4\left (1+\dot{a}_2 \sin \theta _2+\dot{a}_3 \cos \!\left (2 \theta _1+\theta _3-\theta _2\right )+\dot{a}_4 \cos \!\left (2 \theta _1+\theta _3-\theta _2-\theta _4\right )/ 2\right )\right ) a_1^2- \\[5pt] & -\left (m_1-m_2-m_3-m_4\right ) g \sin \theta _1/ 2 \\[5pt] F_2= & \left (m_2\!\left (1/ 2-\dot{a}_1 \sin \theta _2/ 2\right )+m_3\!\left (1-\dot{a}_1 \sin \theta _2+\dot{a}_3 \sin \!\left (2 \theta _1+\theta _3-2 \theta _2\right )/ 2\right )+\right. \\[5pt] & +\left .m_4\left (1+a_1 \sin \theta _2+a_3 \sin \!\left (2 \theta _1+\theta _3-2 \theta _2\right )+a_4 \cos \!\left (2 \theta _1+\theta _3-2 \theta _2-\theta _4\right )/ 2\right )\right ) \dot{a}_2^2- \\[5pt] & -\left (-m_2-m_3-m_4\right ) g \cos \!\left (\theta _1-\theta _2\right )/ 2 \\[5pt] F_3= & \left (m_3\!\left (1/ 4+a_1 \cos \!\left (2 \theta _1+\theta _3-\theta _2\right )/ 2+a_2 \sin \!\left (2 \theta _1+\theta _3-2 \theta _2\right )/ 2\right )+\right. \\[5pt] & +m_4\!\left (1+a_1 \cos \!\left (2 \theta _1+\theta _3-\theta _2\right )+\dot{a}_2 \sin \!\left (2 \theta _1+\theta _3-2 \theta _2\right )+\dot{a}_4 \cos \theta _4\right )/ 2 a_3^2- \\[5pt] & -\left (-m_3-m_4\right ) g \sin \!\left (\theta _1+\theta _3-\theta _2\right )/ 2 \\[5pt] F_4= & \left (m _{ 4 } \!\left (1/ 4+a_1 \cos \!\left (2 \theta _1+\theta _3-\theta _2-\theta _4\right )/ 2+\dot{a}_3 \cos \theta _4/ 2+\right .\right. \\[5pt] & +a_2 \cos \!\left (2 \theta _1+\theta _3-2 \theta _2-\theta _4\right )/ 2 a_4^2-\left (-m_4\right ) g \sin \!\left (\theta _1+\theta _3-\theta _2-\theta _4\right )/ 2 \end{aligned}\right. \end{equation}

(25) \begin{equation} F_i=\left [\begin{array}{l} F_1 \\[5pt] F_2 \\[5pt] F_3 \\[5pt] F_4 \end{array}\right ]=\left [\begin{array}{c@{\quad}c@{\quad}c@{\quad}c} G_1 & 0 & 0 & 0 \\[5pt] 0 & G_2 & 0 & 0 \\[5pt] 0 & 0 & G_3 & 0 \\[5pt] 0 & 0 & 0 & G_4 \end{array}\right ]\left [\begin{array}{l} \ddot{a}_1^2 \\[5pt] a_2^2 \\[5pt] \ddot{a}_3^2 \\[5pt] \ddot{a}_4^2 \end{array}\right ]+\left [\begin{array}{c} \left (m_1-m_2-m_3-m_4\right ) \sin \theta _1 \\[5pt] \left (-m_2-m_3-m_4\right ) \cos \!\left (\theta _1-\theta _2\right ) \\[5pt] \left (-m_3-m_4\right ) \sin \!\left (\theta _1+\theta _3-\theta _2\right ) \\[5pt] \left (-m_4\right ) \sin \!\left (\theta _1+\theta _3-\theta _2-\theta _4\right ) \end{array}\right ] g/ 2 \end{equation}

(26) \begin{equation} \begin{aligned} G_1= & \left (\left (m_1/ 2+2 m_2+m_3+m_4\right )+\left (m_2/ 2+m_3+m_4\right ) \dot{a}_2 \sin \theta _2+\left (m_3/ 2+m_4\right ) \dot{a}_3 \cos \!\left (2 \theta _1+\theta _3-\theta _2\right )+\right. \\[5pt] & +\left .\left (m_4/ 2\right ) a_4 \cos \!\left (2 \theta _1+\theta _3-\theta _2-\theta _4\right )\right ) \\[5pt] G_2= & \left (m_2/ 2+m_3+m_4\right )+\left (-m_2/ 2-m_3+m_4\right ) \dot{a}_1 \sin \theta _2+\left (m_3/ 2+m_4\right ) \dot{a}_3 \sin \!\left (2 \theta _1+\theta _3-2 \theta _2\right )+ \\[5pt] & +\left (m_4/ 2\right ) \dot{a}_4 \cos \!\left (2 \theta _1+\theta _3-2 \theta _2-\theta _4\right ) \\[5pt] G_3= & \left (m_3/ 4+m_4\right )+\left (m_3/ 2+m_4\right ) a_1 \cos \!\left (2 \theta _1+\theta _3-\theta _2\right )+\left (m_3/ 2+m_4\right ) a_2 \sin \!\left (2 \theta _1+\theta _3-2 \theta _2\right )+ \\[5pt] & +\left (m_4/ 2\right ) \dot{a}_4 \cos \theta _4 \\[5pt] G_4= & m_4/ 4+m_4 \dot{a}_1 \cos \!\left (2 \theta _1+\theta _3-\theta _2-\theta _4\right )/ 2+\left (m_4/ 2\right ) \dot{a}_2 \cos \!\left (2 \theta _1+\theta _3-2 \theta _2-\theta _4\right )+ \\[5pt] & +\left (m_4/ 2\right ) \dot{a}_3 \cos \theta _4 \end{aligned} \end{equation}

4.2. Establishment of the dynamic simulation model

In order to further analyze the motion law of hip joint, knee joint and ankle joint in active training of lower limbs, and verify the correctness of theoretical derivation of dynamic model, the dynamic simulation model of five-bar mechanism of man-machine integration system is established in ADAMS environment of 2019 version of MAS company in the United States, as shown in Fig. 4. According to the data in Table III, the model size is set. The length of link of the prototype is 131 mm, the dynamic simulation parameters are shown in Table III. Fixed constraints are set on the upper body and the rehabilitation prototype to ensure that the positions of the two are unchanged, and rotation constraints are added to each joint to ensure that each joint can rotate relatively.

Table III. Body size and key information of the six subjects.

The joint angle associated motion constraint setting plays a key role in the simulation under multi-joint angle uncertainty. In this paper, based on human physiology and the coordinate system of Fig. 3, the constraint relationship among hip joint, knee joint and ankle joint is obtained for the first time by analytical method, such as Eq. (27) and Eq. (28), so as to realize the constraint setting of each joint in the dynamic simulation model. According to Fig. 4, Since the five-bar mechanism of the man-machine integration system has two degrees of freedom, a rotation drive is added to the corresponding revolute joints of the hip joint and the ankle joint.

In order to make the thigh movement stable, the hip joint is driven by sinusoidal speed [Reference Meng, Macleod, Porr and Gollee30].At the same time, the simulation speed should be close to the actual lower limb rehabilitation training speed, so the speed is set to $n_1=10/\bigtriangleup{t}\sin \!\left (2\pi t/5\right )$ and the direction is counterclockwise. According to the constraint relationship between hip joint, knee joint and ankle joint, the motion correlation curve is created [Reference Tan and Y. K.31], according to the hip joint speed through the CUBSPL function in ADAMS environment and then the ankle joint speed is set to $n_{2}=\mathrm{CUBSPL}\!\left (t,0,\mathrm{SPLINE},0\right )^{\circ }$ /s and the direction is counterclockwise. In addition, a single flexion and extension cycle time is set to 5 s, and the number of steps is 50. The simulated measurement angles of the hip, knee joint and ankle joints are the horizontal angle between thigh and x-axis, the angle between thigh and calf, and the angle between calf and foot, the steps of active rehabilitation training of lower limbs are shown in Fig. 5.

(27) \begin{equation} \left \{\begin{array}{r} a_1 \sin \theta _1+a_2 \sin \!\left (\theta _2-\theta _1\right )-a_3 \sin \!\left (\theta _2-\theta _1-\theta _3\right )-a_4 \sin \theta _4=0 \\[5pt] a_1 \cos \theta _1+a_2 \cos \!\left (\theta _2-\theta _1\right )-a_3 \cos \!\left (\theta _2-\theta _1-\theta _3\right )-a_4 \cos \theta _4-a_5=0 \end{array}\right. \end{equation}

(28) \begin{equation} a_2 \cos \!\left (\theta _2-\theta _1\right )=a_4 \cos \theta _4+a_3 \cos \!\left (\theta _2-\theta _1-\theta _3\right )-a_1 \cos \theta _1 \end{equation}

Figure 5 shows the movement steps of lower limb rehabilitation training in a single cycle. In the first half of the cycle, the ankle joint angle is the smallest, then the knee joint angle is the smallest, and then the hip joint angle is the largest. The appearance of the three angles shows a certain lag law. In the second half of the cycle, the law shows the opposite trend.

By simulating that the angle changes curves of hip joint, knee joint and ankle joint are obtained. In order to make the results convenient for analysis, each angle is transformed into the angle between lower limb and vertical direction, and the angle change law of three joints in a single flexion and extension cycle is obtained, as shown in Fig. 6. It can be observed that the angle range of hip joint flexion and extension are 80.02 $^{\circ }\sim$ 99.98 $^{\circ }$ , and the angle range of knee joint flexion and extension are 62.62 $^{\circ }\sim$ 109.05 $^{\circ }$ . The angle range of ankle joint flexion and extension is 51.24 $^{\circ }\sim$ 101.47 $^{\circ }$ .

Figure 5. The steps of active rehabilitation training of lower limbs.

Figure 6. Angle changes of the hip, knee and ankle joints.

The results show that the motion amplitude of the hip joint is the smallest, the motion amplitude of the knee joint and the ankle joint is larger and the peak angle shows a stable lag change compared with the motion law of the hip joint, but the change rules of the two are similar. The angle change curve of each angle is continuous and smooth transition, and the angle change range of the three joints is within a reasonable range, which is in line with the active training of lower limbs.

At the same time, considering that the knee joint is more vulnerable to injury than the ankle joint in actual activities [Reference Xia, Huan and Chen32, Reference Shao, Lu, Cen, Zheng, Sun and Gu33], and the ankle joint bears most of the body’s gravity in daily actions, so the joint rehabilitation is clinically suitable for resting and drug treatment, so this paper further selects the hip joint and knee joint to analyze the movement law in lower limb rehabilitation training.

4.3. Simulation of different active training speeds of lower limbs

In order to analyze the influence of different speed training on the rehabilitation [Reference Bicer, Phillips and Modenese34], the sine speed $N=n/\bigtriangleup{t}\sin \!\left (2\pi t/5\right )^{\circ }$ /s with different amplitude is further used to simulate the effect of rehabilitation training at different speeds. The selected speed is close to the actual rehabilitation training requirements to obtain more accurate rehabilitation rules. It is found that the angle, angular velocity and angular acceleration peak value and extreme point of both hip and knee joints change with “n,” but the change rule is basically unchanged.

Figure 7. Hip joint angular velocity and angular acceleration curve. (a) Hip joint angular velocity curve; (b) Hip joint angular acceleration curve.

From Fig. 7, it can be seen that the angular velocity and angular acceleration of the hip joint increase with the increase of the rotational speed, and the speed is consistent at the limit position of 1.2 s and 3.7 s, but at this time, the acceleration is positively correlated with the rotational speed, and the curve gradually appears sharp points with the increase of the rotational speed. Therefore, the rotational speed needs to be controlled within a reasonable range to reduce the secondary damage to the lower limbs caused by the impact of the limit position.

By obtaining the knee joint angular velocity and angular acceleration curves, as shown in Fig. 8, the reasonable velocity range is further analyzed.

Figure 8. Angular velocity and acceleration curves of hip joints at a different velocity. (a) n = 6; (b) n = 7; (c) n = 8; (d) n = 9; (e) n = 10; (f) n = 11.

Figure 8 shows the angular velocity of the knee joint at different speeds showed similar changes. When the speed n is in the range of 6 $\sim$ 7, the angular acceleration increases from 0 to 0.015 $^{\circ }$ / $s^2$ and then decreases to 0.005 $^{\circ }$ / $s^2$ at 0.4 $\sim$ 0.5s at the startup moment. Therefore, the acceleration fluctuation is the smallest at the speed n = 7, and the secondary damage can be avoided by extending the startup time. The angular velocity was reduced from 0.013 to 0 at 1.1 s, indicating that the rehabilitation training action was smoother when n<7. The angular acceleration increased from 40 $^{\circ }$ / $s^2$ to 93 $^{\circ }$ / $s^2$ in the range of n = 6 $\sim$ 9, and gradually decreased to 32 $^{\circ }$ / $s^2$ in the range of n = 9 $\sim$ 11. The smaller the maximum angular acceleration, the smoother the rehabilitation action. In the 3.2 $\sim$ 4.2 s time period, the angular acceleration of the speed n = 8 $\sim$ 11 fluctuates many times, while the angular velocity n = 6 $\sim$ 8 is smooth and excessive, and there is no sudden change in the acceleration when n = 6 $\sim$ 7. Therefore, in order to reduce the secondary injury of the lower limbs caused by the acceleration mutation, the hip joint speed of no more than $7/t^{\ast}\sin (2^{\ast}\Pi/5)^{\circ }$ /s should be selected for lower limb active rehabilitation training to achieve better results.

In order to verify the feasibility of the selected reasonable speed, the speed and acceleration can reflect the speed of the lower limb movement, while the angular velocity and angular acceleration reflect the speed of the lower limb rotation, and then analyze the motion characteristics of the lower limb at multiple levels. The speed, acceleration, angular velocity and angular acceleration of the thigh, calf and foot at $N=7/\bigtriangleup{t}\sin \!\left (2\pi t/5\right )^{\circ }$ /s are comprehensively analyzed, and the speed and acceleration curves are obtained as shown in Fig. 9.

Figure 9. The motion law of thigh, calf and foot. (a) Velocity curves; (b) Acceleration curves; (c) Angular velocity curves; (d) Angular acceleration curves.

It can be seen from Fig. 9(a) and Fig. 9(b) that the speed changes of thigh, calf and foot are similar, and there is no sudden change point. The delay of movement from thigh to foot is consistent with the actual movement law, indicating that the simulation method is correct. At the same time, the maximum velocity of thigh, calf and foot is 52 $^{\circ }$ /s, 130 $^{\circ }$ /s and 150 $^{\circ }$ /s, respectively, and the acceleration is 68 $^{\circ }$ / $s^2$ , 163 $^{\circ }$ / $s^2$ and 219 $^{\circ }$ / $s^2$ , respectively. At 2.5s, the acceleration of the calf fluctuates. At this time, the knee joint moves to the highest position in the rehabilitation training. Therefore, at this time, the secondary injury can be reduced by slowing the movement. It can be seen from Fig. 9(c) and Fig. 9(d) that the angular velocity and angular acceleration curves of thigh, calf and foot are smooth, indicating that they move smoothly. However, the angular velocity of the calf rotates at 0.5 s without obvious fluctuation, so that it will not cause secondary damage to the rehabilitation training. The maximum angular acceleration is 54 $^{\circ }$ / $s^2$ , and there is no obvious fluctuation in the whole cycle. Therefore, there is no sudden change in speed and angular acceleration, indicating that the active training of lower limbs based on the selected speed meets the requirements, which proves the feasibility of the speed.

5.1. Joint angle experiment

The hip and knee joint angles were measured by the joint angle experimental system shown in Fig. 10. The BWT61CL type attitude sensor with an accuracy of 0.05 $^{\circ }$ angle produced by Shenzhen Weite Intelligent Company is selected to collect the joint angle, and the sensor is connected to the computer and the initial parameter setting is carried out through Bluetooth. In the experiment, the sensor was placed near the joint position to measure the corresponding angle [Reference Li, Xue, Yang, Han and Zhang35].

Figure 10. Joint angle experimental system.

The subject’s feet were placed at the rotating shaft of the prototype, and the initial position was kept at the level of the thighs, and then the legs simultaneously exerted force to rotate the prototype’s rotating shaft in the same direction. In order to reduce the experimental error caused by the individual differences of different subjects, all subjects maintained an unified sitting posture during the experiment. Through the pre-experiment, the data of the same training duration were recorded after the lower limbs were kept at a constant speed. Each subject repeated three experiments for a total of 60 sets of data. The dynamic simulation model was established by measuring the size of each subject’s lower limb, and the angles of hip joint and knee joint were obtained by the same simulation method in sections 3.1. The angle was compared with the experimental angle to obtain more accurate lower limb motion characteristics.

In this experiment, 20 subjects (10 males and 10 females) were recruited, who were graduate students of Beijing University of Information Science and Technology (170 $\pm$ 8 cm in height, 67 $\pm$ 8 kg in weight, 25 $\pm$ 2 years old), the body size and key information were shown in Table IV, and the corresponding body mass index (BMI) was measured. BMI is calculated by dividing weight (Kg) by the square of height (m). It is currently widely used to measure the degree of body fatness, the different body types of people to be taken into account so as to eliminate the influence of personal factors. without skeletal and limb motor dysfunction. They were familiar with the experimental process and in good mental condition during the experiment. They volunteered to experiment and signed the informed consent of the experiment.

5.2. Results and discussion

The angles of the hip and knee joints measured by the active training experiment of a single subject in a flexion and extension cycle were compared with the dynamic simulation results, as shown in Fig. 11.

It can be seen from Fig. 11 that the maximum angles of hip joint experiment and simulation are 110.06 $^{\circ }$ and 105.15 $^{\circ }$ respectively, and the minimum angles are 72.64 $^{\circ }$ and 76.10 $^{\circ }$ , respectively. The maximum angles of knee joint experiment and simulation were 116.19 $^{\circ }$ and 109.54 $^{\circ }$ , respectively, and the minimum angles were 56.30 $^{\circ }$ and 59.10 $^{\circ }$ , respectively. The results show that the experimental and simulation joint angles have the same trend, and the experimental and simulation angles are within the reasonable range of human joint motion. The maximum and minimum angles appear at the highest and lowest points of joint movement in rehabilitation training, respectively. At this time, the maximum experimental and simulation errors also indicate that the lower limb movement to this position should be slowed down. The angular range of motion of the hip joint and knee joint experiment is greater than that of the simulation. The reason is that the inertia of the subject causes the prototype shaft to disengage from the limb when it is in a vertical position, and this phenomenon will become more obvious as the speed increases.

Through simulation, the joint motion law under the actual rehabilitation training conditions is obtained, and the secondary damage caused by the impact is reduced. Further comparison with the experimental results is carried out, and the errors and causes are analyzed to improve the effectiveness of rehabilitation training. The maximum error of the hip joint appears at 1.8 s and 4.2 s, and the knee joint is at 2.1 s and 4.3 s, both of which appear near the limit position of the joint angle, corresponding to the highest and lowest positions of the rehabilitation prototype shaft in the foot during active training. The maximum error of the hip joint is 9.84%, and the maximum error of the knee joint is 4.51%, both within the allowable range of 15%, which were within the allowable range of 15% error verified in Reference [Reference Ataabadi, Sarvestan and Alaei36, Reference Khuyagbaatar, Purevsuren and Ganbat37], meeting the needs of rehabilitation training.

Table IV. Body size and key information of the six subjects.

Figure 11. Comparison of hip and knee joint experiments and simulation angles under active training.

In order to obtain a more accurate joint angle change rule of lower limb active training, the experimental data of 20 subjects were further processed, and the corresponding body mass index (BMI) was measured. BMI is calculated by dividing weight (Kg) by the square of height (m). It is currently widely used to measure the degree of body fatness, so as to take into account the different body types of people and eliminate the influence of personal factors. Due to the close height, the lower limb movement rules are similar. Therefore, a total of 18 groups of representative 6 subjects were selected with height as the reference object. In order to make the experimental data more analytical and accurate, the same experimental method is used to control the same range of hip joint activity of all subjects. Setting the range of hip joint activity can reduce the secondary injury caused by the excessive range of hip joint activity and simulate the actual movement law of different groups of people under the same active training, so as to obtain a personalized rehabilitation training program. Therefore, the range of motion of the joint angle is obtained and the maximum error absolute values( $\lvert \delta _{max}\rvert$ ) of the maximum and minimum joint angles of the simulation and experiment are analyzed respectively. The simulation and experimental data of the lower limb active training of 6 subjects are obtained, as shown in Table V.

Table V. Simulation and experimental data of lower limb active training of six subjects.

From the data in Table V, it can be seen that under the same training speed and man-machine distance, the simulated hip joint angle range of motion is consistent, and the knee joint angle range of motion varies with height, weight and BMI. The larger the size of the lower limbs, the smaller the experimental angle range of the hip joint and the knee joint, and the smaller the error of the hip joint. The knee joint is the opposite, indicating that the hip joint training is better than the knee joint to maintain a smaller distance between the human and the machine, while the knee joint is the opposite. The maximum range of hip and knee joint angles of the subjects were 72.44 $^{\circ }\sim$ 110.19 $^{\circ }$ and 56.30 $^{\circ }\sim$ 113.97 $^{\circ }$ , respectively, which were within the allowable range. The error range of hip joint is $2.13\%\sim 9.84\%$ , and the error range of knee joint is $1.70\%\sim 8.17\%$ , both of which are less than 10%, which meets the experimental requirements.

In summary, the experimental hip and knee joint angles are consistent with the simulation angle, and the errors are within the allowable range. The reasons are as follows: (1) The human body size and training speed of the experiment and simulation cannot be exactly the same, and the experiment is influenced by limb inertia; (2) During the experiment, the vibration of the lower limbs caused the measurement data of the attitude angle sensor to be too large; (3) In the simulation, the lower limb movement is a relatively simple plane motion, while in the experiment, the lower limb has a complex spatial motion, which cannot be simulated.

Combining the errors of the hip and the knee joint, the angle of the hip and the knee joint measured by the subject 4 experiment was substituted into the 2.1 kinematics model to obtain the overall working space of the lower limb, as shown in Fig. 12.

Figure 12. Workspace of lower limb in active training.

It can be seen that the working space of the lower limb is a ring with an outer ring radius of 612 mm in the xoy plane, and its inner diameter of 130 mm is close to the length of the rehabilitation prototype shaft. The width of the ring corresponds to the sum of the length of the calf and the foot. Therefore, the simulation workspace verifies the theoretical calculation in Section 3.1 and conforms to the actual situation, indicating that the kinematic model is correct, which provides a theoretical reference for the design and gait planning of lower limb rehabilitation robots.

Under the same speed and man-machine distance, the lower limb movements of subjects with different BMI, height or weight are different. Through the analysis of this difference, the personalized training rules of different lower limb sizes are explored. Therefore, the relationship between the subjects with different BMI and the range of motion of the joint angle is obtained, as shown in Fig. 13.

Figure 13. The range of joint motion in subjects with different BMI.

From Fig. 13, it can be seen that the range of hip joint angle changes of subjects with different BMI is basically the same, which reflects that the change of hip joint is mainly affected by training speed. At the same time, it can be seen that under the same training speed and man-machine distance, the angle of hip joint and knee joint is the largest at the same time when BMI is 22.5, but the angle of hip joint is the smallest at 22.4, and the knee joint is the smallest at 21.3, indicating that the change of hip joint angle has little effect with BMI, while the change of knee joint angle has significant effect.

The relationship between the angle error of the hip joint and the knee joint and the height and weight is shown in Fig. 14 and Fig. 15 respectively.

Figure 14. Relationship between hip joint angle error and height and weight. (a) Man-machine distance diagram. (b) Relationship between hip angle error and body weight. (c) Relationship between error and height. (d) Relationship between error and weight. (e) Relationship between error and BMI.

Figure 15. Relationship between knee joint angle error and height and weight. (a) Man-machine distance diagram. (b) Relationship between knee angle error and body weight. (c) Relationship between error and height. (d) Relationship between error and weight. (e) Relationship between error and BMI.

From Fig. 14(a) and (b) and Fig. 15(a) and (b), it can be seen that the angle error of hip joint experiment and simulation decreases with the increase in body height and weight. The angle error of knee joint experiment and simulation increases with the increase in height and weight. It can be seen from Fig. 14(c)–(e) and Fig. 15(c)–(e) that the error of hip joint and knee joint is closer to the change rule of height, and the error is also different under the same BMI, indicating that the height which is proportional to the size of the lower limb is the main factor of error, and the weight is the secondary factor. At the same time, according to Section 3.2, it is found that the angle of joint activity decreases with the decrease of training speed. Due to the same man-machine distance, the higher the height, the larger the lower limb size, the smaller the distance between the hip joint and the knee joint and the center of the rehabilitation prototype shaft, indicating that the smaller the distance, the smaller the hip joint error, and the knee joint is the opposite.

Therefore, for the active training of lower limbs based on hip joint, the smaller the height and weight, the lower the training speed or the smaller the man-machine position distance can be appropriately reduced to obtain a hip joint training range closer to the theory; for knee-based rehabilitation training, the training speed can be appropriately reduced or the man-machine distance can be increased to reduce the error to obtain a range of activity angles closer to the knee joint theory, so as to achieve better training results.

For the hip and knee joints, the above strategies are verified by changing the man-machine distance by taking the 5 cm reduction and increase as examples, and the error effects of height, weight and BMI on the two are obtained respectively, as shown in Fig. 16.

Figure 16. Comparison of joint errors before and after adjustment. (a) Height-hip error; (b) Weigh-hip error; (c) Height-knee error; (d) Weight-knee error; (e) BMI-hip error; (f) BMI-knee error.

It can be seen from Fig. 16 that with the increase of man-machine distance, the hip joint error decreases with different body height and weight, and the knee joint also has the same effect with the decrease of man-machine distance. The errors of hip joint and knee joint are reduced by 1.2% and 4.6% respectively, and the errors of BMI on hip joint and knee joint are reduced by 1.2% and 4.1% respectively. Therefore, it shows that the accuracy of rehabilitation training can be improved by adjusting the training method of human-machine distance planning. At the same time, the feasibility of the personalized planning method proposed in this paper is verified.