**3.1 Kinematics**

Denavit-Hartenberg (DH) method is used to analyze the kinematics of human motion [14]. In the Denavit-Hartenberg method any point *<sup>i</sup>* **r** can be transferred to the global reference frame as <sup>0</sup> **r** in **Figure 3** and it can be presented as Eq. (1).

**Figure 3.** *Articulated chain.*

$$\mathbf{^0r} = \mathbf{^0T\_i} \mathbf{^ir} \tag{1}$$

*Qi* <sup>¼</sup> *<sup>d</sup> dt*

and *qi* is the generalized coordinate. Total kinetic energy is derived as

*<sup>K</sup> <sup>j</sup>* <sup>¼</sup> <sup>1</sup> 2 X*n i*¼1

*m <sup>j</sup>***g***<sup>T</sup>***A** *<sup>j</sup> j* **<sup>r</sup>** *<sup>j</sup>* �X*<sup>n</sup> j*¼1 **f** *T <sup>k</sup>* **A***<sup>j</sup> k*

gravity force vector, **f***<sup>k</sup>* is an external force which is defined in global reference frame and acting on the body segment expressed in *k*th reference frame, *<sup>k</sup>*

> X*n j*¼*i m <sup>j</sup> i* **T** *j j* **r** *<sup>j</sup>* � **f** *T k ∂***A***<sup>i</sup> ∂qi*

where *τ<sup>i</sup>* is the joint torque acting on the joint represented with generalized coordinate *qi*, **G***<sup>i</sup>* is external moment which is defined in the global reference frame

transformation the mechanical system of human body. Obtained joint torques are used to evaluate joint torque constraints and objective functions in later section.

Optimization based motion simulation in this chapter is performed according to

Kinematics and dynamics are covered in previous section, so we will discuss the

joint variables, objective function, and constraints evaluation in this section.

To generate weight lifting motion, our optimization variable is joint angle profiles. These joint angle profiles are approximated using B-spline function

. Ground reaction force can be calculated using global force

location of external force acting on the link expressed in the *k*th reference frame, δ*jk* is Kronecker delta. Then, the equations of motion can be derived from above

*<sup>K</sup>* <sup>¼</sup> <sup>X</sup>*<sup>n</sup> j*¼1

where **A***<sup>i</sup>* is joint angle matrix and **J***<sup>i</sup>* is inertia matrix at *i*

*j*¼1

where *mj* is the mass of the body segment represented as *j*

� **<sup>g</sup>***<sup>T</sup> <sup>∂</sup>***A***<sup>i</sup> ∂qi*

Step2: Function approximations for joint variables.

Step7: Print out if converge. Otherwise go back to step 2.

*<sup>V</sup>* ¼ �X*<sup>n</sup>*

the potential energy for system can be given as

*DOI: http://dx.doi.org/10.5772/intechopen.98391*

equations (Eqs. (7)–(9))

*∂***A***<sup>i</sup> ∂qi*

and **z**<sup>0</sup> is [0 0 1 0]T

the following process:

**4.1 Variables**

**163**

X*n j*¼*i i* **<sup>T</sup>** *<sup>j</sup>***J***j***A**€ *<sup>T</sup> j*

!

**4. Optimization formulation**

Step1: Prepare input data.

Step3: Kinematics analysis. Step4: Dynamics analysis.

Step5: Objective function evaluation. Step6: Constraints evaluation.

*τ<sup>i</sup>* ¼ *tr*

*∂L ∂q*\_ *i* � �

*Optimization Based Dynamic Human Motion Prediction with Modular Exoskeleton Robots…*

where *L* is Lagrangian which is *L=K – V* (kinetic energy – potential energy),

� *∂L ∂qi*

*tr* **<sup>A</sup>**\_ *<sup>i</sup>***J***i***A**\_ *<sup>T</sup> i*

� �*:* (8)

X*n j*¼*i i* **T** *j k*

th reference frame. Also,

th reference frame, **g** is

*rf δjk* � **G***i***A***<sup>i</sup>*�<sup>1</sup>**z**<sup>0</sup>

**r***<sup>f</sup>* is the

(10)

**r** *<sup>f</sup> δjk* (9)

(7)

Transformation matrix <sup>0</sup> **T***<sup>i</sup>* in Eq. (1) can be obtained as follows.

$$\mathbf{T}^0 \mathbf{T}\_i = {}^0 \mathbf{T}\_1 {}^1 \mathbf{T}\_2 \cdots {}^{i-1} \mathbf{T}\_i = \prod\_{n=1}^i {}^{n-1} \mathbf{T}\_n \tag{2}$$

$$\mathbf{T}\_{i} = \begin{bmatrix} \cos q\_i & -\cos a\_i \sin q\_i & \sin a\_i \sin q\_i & a\_i \cos q\_i \\ \sin q\_i & \cos a\_i \cos q\_i & -\sin a\_i \cos q\_i & a\_i \sin q\_i \\ \mathbf{0} & \sin a\_i & \cos a\_i & d\_i \\ \mathbf{0} & \mathbf{0} & \mathbf{0} & \mathbf{1} \end{bmatrix} \tag{3}$$

where the parameters *q*, *α*, *a*, *d*, are DH parameters. The *qi* is the joint angle between the *xi*-1 axis and the xi axis about the *zi*-1 axis according to the right hand rule. The *di* is the distance between the origin of the *i*-1th coordinate frame and the intersection of the *zi*-1 axis with the *xi* axis along the *zi*-1 axis. *ai* is the distance between the intersection of the *zi*-1 axis with the *xi* axis and the origin of the *i* th frame along the *xi* axis. *α<sup>i</sup>* is the angle between the *zi*-1 axis and the *zi* axis about the *xi* axis according to the right hand rule. In here, we use *qi* as our generalized coordinates. As mentioned before, the inertial reference frame is located at the point between foot. The origin of inertial reference frame is O in the above **Figure 2**. Thus, all kinematic chain starts from origin of inertial reference frame. For the efficiency of calculation time, we use recursive way for kinematic information as follows:

$$\mathbf{A}\_{i} = \mathbf{A}\_{i-1} \mathbf{T}\_{i} \tag{4}$$

$$\mathbf{B}\_{i} = \dot{\mathbf{A}}\_{i} = \mathbf{B}\_{i-1} \mathbf{T}\_{i} + \mathbf{A}\_{i-1} \frac{\partial \mathbf{T}\_{i}}{\partial q\_{i}} \dot{q}\_{i} \tag{5}$$

$$\mathbf{C}\_{i} = \dot{\mathbf{B}}\_{i} = \mathbf{C}\_{i-1}\mathbf{T}\_{i} + 2\mathbf{B}\_{i-1}\frac{\partial \mathbf{T}\_{i}}{\partial q\_{i}}\dot{q}\_{i} + \mathbf{A}\_{i-1}\frac{\partial^{2}\mathbf{T}\_{i}}{\partial q\_{i}^{2}}\dot{q}\_{i}^{2} + \mathbf{A}\_{i-1}\frac{\partial \mathbf{T}\_{i}}{\partial q\_{i}}\ddot{q}\_{i} \tag{6}$$

where **A***<sup>i</sup>* is position matrix, **B***<sup>i</sup>* is velocity matrix, and **C***<sup>i</sup>* is acceleration matrix.

#### **3.2 Dynamics**

Once we obtain the kinematic information, we can use them to calculate dynamics of motion simulation. The equations of motion are derived from Lagrange's equation. *Optimization Based Dynamic Human Motion Prediction with Modular Exoskeleton Robots… DOI: http://dx.doi.org/10.5772/intechopen.98391*

$$Q\_i = \frac{d}{dt} \left(\frac{\partial L}{\partial \dot{q}\_i}\right) - \frac{\partial L}{\partial q\_i} \tag{7}$$

where *L* is Lagrangian which is *L=K – V* (kinetic energy – potential energy), and *qi* is the generalized coordinate. Total kinetic energy is derived as

$$K = \sum\_{j=1}^{n} K\_j = \frac{1}{2} \sum\_{i=1}^{n} tr\left(\dot{\mathbf{A}}\_i \mathbf{J}\_i \dot{\mathbf{A}}\_i^T\right). \tag{8}$$

where **A***<sup>i</sup>* is joint angle matrix and **J***<sup>i</sup>* is inertia matrix at *i* th reference frame. Also, the potential energy for system can be given as

$$V = -\sum\_{j=1}^{n} m\_j \mathbf{g}^T \mathbf{A}\_j{}^j \overline{\mathbf{r}}\_j - \sum\_{j=1}^{n} \mathbf{f}\_k^T \mathbf{A}\_j{}^k \mathbf{r}\_f \delta\_{jk} \tag{9}$$

where *mj* is the mass of the body segment represented as *j* th reference frame, **g** is gravity force vector, **f***<sup>k</sup>* is an external force which is defined in global reference frame and acting on the body segment expressed in *k*th reference frame, *<sup>k</sup>* **r***<sup>f</sup>* is the location of external force acting on the link expressed in the *k*th reference frame, δ*jk* is Kronecker delta. Then, the equations of motion can be derived from above equations (Eqs. (7)–(9))

$$\tau\_i = \text{tr}\left(\frac{\partial \mathbf{A}\_i}{\partial q\_i} \sum\_{j=i}^n i \mathbf{T}\_j \mathbf{J}\_j \tilde{\mathbf{A}}\_j^T \right) - \mathbf{g}^T \frac{\partial \mathbf{A}\_i}{\partial q\_i} \sum\_{j=i}^n m\_j{}^j \mathbf{T}\_j{}^j \tilde{\mathbf{r}}\_j - \mathbf{f}\_k^T \frac{\partial \mathbf{A}\_i}{\partial q\_i} \sum\_{j=i}^n i \mathbf{T}\_j{}^k r\_f \delta\_{jk} - \mathbf{G}\_i \mathbf{A}\_{i-1} \mathbf{z}\_0 \tag{10}$$

where *τ<sup>i</sup>* is the joint torque acting on the joint represented with generalized coordinate *qi*, **G***<sup>i</sup>* is external moment which is defined in the global reference frame and **z**<sup>0</sup> is [0 0 1 0]T . Ground reaction force can be calculated using global force transformation the mechanical system of human body. Obtained joint torques are used to evaluate joint torque constraints and objective functions in later section.

## **4. Optimization formulation**

Optimization based motion simulation in this chapter is performed according to the following process:

Step1: Prepare input data.

Step2: Function approximations for joint variables.

Step3: Kinematics analysis.

Step4: Dynamics analysis.

Step5: Objective function evaluation.

Step6: Constraints evaluation.

Step7: Print out if converge. Otherwise go back to step 2.

Kinematics and dynamics are covered in previous section, so we will discuss the joint variables, objective function, and constraints evaluation in this section.

#### **4.1 Variables**

To generate weight lifting motion, our optimization variable is joint angle profiles. These joint angle profiles are approximated using B-spline function

<sup>0</sup>**<sup>r</sup>** <sup>¼</sup> <sup>0</sup>**T***<sup>i</sup> i*

where the parameters *q*, *α*, *a*, *d*, are DH parameters. The *qi* is the joint angle between the *xi*-1 axis and the xi axis about the *zi*-1 axis according to the right hand rule. The *di* is the distance between the origin of the *i*-1th coordinate frame and the intersection of the *zi*-1 axis with the *xi* axis along the *zi*-1 axis. *ai* is the distance between the intersection of the *zi*-1 axis with the *xi* axis and the origin of the *i*

along the *xi* axis. *α<sup>i</sup>* is the angle between the *zi*-1 axis and the *zi* axis about the *xi* axis according to the right hand rule. In here, we use *qi* as our generalized coordinates. As mentioned before, the inertial reference frame is located at the point between foot. The origin of inertial reference frame is O in the above **Figure 2**. Thus, all kinematic chain starts from origin of inertial reference frame. For the efficiency of calculation

<sup>B</sup>*<sup>i</sup>* <sup>¼</sup> <sup>A</sup>\_ *<sup>i</sup>* <sup>¼</sup> <sup>B</sup>*<sup>i</sup>*�1T*<sup>i</sup>* <sup>þ</sup> <sup>A</sup>*<sup>i</sup>*�<sup>1</sup>

*∂*T*<sup>i</sup> ∂qi*

*q*\_ *<sup>i</sup>* þ A*<sup>i</sup>*�<sup>1</sup>

where **A***<sup>i</sup>* is position matrix, **B***<sup>i</sup>* is velocity matrix, and **C***<sup>i</sup>* is acceleration matrix.

Once we obtain the kinematic information, we can use them to calculate dynamics of motion simulation. The equations of motion are derived from Lagrange's equation.

<sup>0</sup>**T***<sup>i</sup>* <sup>¼</sup> <sup>0</sup>**T**<sup>1</sup>

time, we use recursive way for kinematic information as follows:

<sup>C</sup>*<sup>i</sup>* <sup>¼</sup> <sup>B</sup>\_ *<sup>i</sup>* <sup>¼</sup> <sup>C</sup>*<sup>i</sup>*�1T*<sup>i</sup>* <sup>þ</sup> 2B*<sup>i</sup>*�<sup>1</sup>

1 **T**2⋯*<sup>i</sup>*�<sup>1</sup>

**T***<sup>i</sup>* in Eq. (1) can be obtained as follows.

**<sup>T</sup>***<sup>i</sup>* <sup>¼</sup> <sup>Y</sup> *i*

cos *qi* � cos *α<sup>i</sup>* sin *qi* sin *α<sup>i</sup>* sin *qi ai* cos *qi* sin *qi* cos *α<sup>i</sup>* cos *qi* � sin *α<sup>i</sup>* cos *qi ai* sin *qi* 0 sin *α<sup>i</sup>* cos *α<sup>i</sup> di* 00 0 1

*n*¼1

*n*�1

A*<sup>i</sup>* ¼ A*<sup>i</sup>*�1T*<sup>i</sup>* (4)

*q*\_ *<sup>i</sup>* (5)

€*qi* (6)

*∂*T*<sup>i</sup> ∂qi*

*∂*T*<sup>i</sup> ∂qi*

*∂*2 T*i ∂q*<sup>2</sup> *i q*\_ 2 *<sup>i</sup>* þ A*<sup>i</sup>*�<sup>1</sup>

Transformation matrix <sup>0</sup>

*Collaborative and Humanoid Robots*

**Figure 3.** *Articulated chain.*

> *i*�1 **T***<sup>i</sup>* ¼

**3.2 Dynamics**

**162**

**r** (1)

**T***<sup>n</sup>* (2)

(3)

th frame

approximation [15]. The control points in B-spline function approximation are updated in each iteration through the optimization process. Also, we used clamped B-spline which is that the starting point and end point of spline curve are matching to the control points. Following equation describes the B-spline approximation for the joint angle profiles *qi*.

$$q\_i(\mathbf{t}, \mathbf{p}) = \sum\_{j=0}^{m} N\_j(\mathbf{t}) p\_j \tag{11}$$

**4.4 Analytical gradients**

*∂q*1 *∂x*<sup>1</sup> þ *∂C ∂q*2

*∂q*1 *∂x*<sup>2</sup> þ *∂C ∂q*2

*∂q*1 *∂xm* þ *∂C ∂q*2

**5. Numerical experiments**

**5.1 Weight lifting motion**

**5.2 Simulation setting**

Commercial software MATLAB.

**165**

*∂C ∂x*<sup>1</sup>

*∂C ∂x*<sup>2</sup>

*∂C ∂xm* <sup>¼</sup> *<sup>∂</sup><sup>C</sup> ∂q*1

<sup>¼</sup> *<sup>∂</sup><sup>C</sup> ∂q*1

<sup>¼</sup> *<sup>∂</sup><sup>C</sup> ∂q*1

constraints and objective function are as follows:

þ ⋯ þ

þ ⋯ þ

þ ⋯ þ

¼ *∂τ ∂q j*

objective function is calculated by using Eqs. (12) and (18).

*∂C ∂qn*

*∂C ∂qn*

> *∂C ∂qn*

> > *∂q j ∂xk* þ *∂τ ∂q*\_ *j*

*∂qn ∂xm* þ *∂C ∂q*\_ 1

*∂qn ∂x*<sup>1</sup> þ *∂C ∂q*\_ 1

*∂qn ∂x*<sup>2</sup> þ *∂C ∂q*\_ 1

*∂q*2 *∂x*<sup>1</sup>

*DOI: http://dx.doi.org/10.5772/intechopen.98391*

*∂q*2 *∂x*<sup>2</sup>

⋮

*∂q*2 *∂xm*

> *∂τ ∂xk*

The analytical gradients of each constraints and objective function are provided to the optimization solver. These analytical gradients improve the accuracy as well as calculation time for convergence in optimization process. Analytical gradients for

*Optimization Based Dynamic Human Motion Prediction with Modular Exoskeleton Robots…*

*∂q*\_ 1 *∂x*<sup>1</sup>

*∂q*\_ 1 *∂x*<sup>2</sup>

*∂q*\_ *j ∂xk* þ *∂τ ∂*€*q j*

where *C* represents constraints, *qi* is joint angle profile, *xi* is control point of Bspline function approximation, *τ* is joint torque profile. The analytical gradients of each constraint are calculated in the form of Eq. (17). For the analytical gradient of

**Figure 4** depicts weight lifting motion. The weight is defined as W and it is applied as and external load in dynamics equilibrium equations of motion. Weight object is virtual and only hand location is guided by position constraints. Foot is located on the ground and weight is moving up vertically in 0.15 m from 0.5 m above ground and away from heal position by 0.25 m as shown in **Figure 4**. Dynamic stability constraint which is ZMP are imposed on the foot support region area. Joint angle limits and joint torque limits are imposed based on the literature review. It is assumed that the assistive moment is applied hip joint and knee joint

For the optimization, sequential quadratic programming (SQP) is adopted. The sequential quadratic programming is very effective for the large scale nonlinear constrained optimization problem. The commercial software SNOPT is used which is well known as the effectiveness for large scale nonlinear problem SQP solver [21]. Total number of variables for optimization is 330 and total number of constrains is 1,766. The weighting parameters in objective function are set to equally in current study. Numerical experiment is performed while the weights are set to 15 kg and 30 kg. Also, different assistive moments from exoskeleton robot are tested in the experiment. The postprocessing for animation and snapshot was performed using

which axes are parallel to the horizontal axis of inertial reference frame.

*∂q*\_ 1 *∂xm*

þ ⋯ þ

þ ⋯ þ

þ ⋯ þ

*∂*€*q j ∂xk*

*∂C ∂q*\_ *n*

*∂C ∂q*\_ *n* *∂q*\_ *n ∂x*<sup>1</sup>

*∂q*\_ *n ∂x*<sup>2</sup>

*∂C ∂q*\_ *n*

*∂q*\_ *n ∂xm*

þ ⋯ þ

þ ⋯ þ

*∂C ∂*€*qn*

*∂C ∂*€*qn*

þ ⋯ þ

*∂*€*qn ∂x*<sup>1</sup>

*∂*€*qn ∂x*<sup>2</sup>

*∂C ∂*€*qn*

*∂*€*qn ∂xm* (17)

(18)

where **t** = {*t*1, … ,*ts*} is knot vector, **p** = {*p*1, … ,*pm*} is the control points, *Nj* is the basis function of B-spline function approximation. Here, control points **P** becomes optimization variables.

#### **4.2 Performance measure and objective function**

We define the energy consuming (dynamic effort) for the given motion as the performance measure for weight lifting motion prediction and simulation. In this study, we use the joint torque square which is proportional to the mechanical energy. This mechanical energy is a reasonable criterion to minimize [4]. Then, it is formulated as objective function in optimization formulation as follows:

$$f = \int\_0^T \left[\sum w\_1 \tau\_i^2\right] dt\tag{12}$$

where *wi* is the weighting parameters of each joint, *τ<sup>i</sup>* is the joint torque of each joint. The joint torque is calculated from above dynamics equation (Eq. (10)).

#### **4.3 Constraints**

Constraints are the one of the motion control way in this optimization formulation and we used minimal set of constraints to generate motion in this formulation. The list of constraints are as follows:


The joint angle limits and torque limits for the human motion are determined based on the literature [16–19]. Zero Moment Point(ZMP) method is used for the stability condition constraint [20]. Each constraint is formulated as follows accordingly:

$$\mathbf{q}^{L} \le \mathbf{q}(t) \le \mathbf{q}^{U} \tag{13}$$

$$
\mathfrak{r}^L \le \mathfrak{r}(t) \le \mathfrak{r}^U \tag{14}
$$

$$\mathbf{r}\_i(t) = \mathbf{\tilde{r}}\_i \tag{15}$$

$$z\_{xmp} \in \text{FSR}, \ \varkappa\_{xmp} \in \text{FSR} \tag{16}$$

*Optimization Based Dynamic Human Motion Prediction with Modular Exoskeleton Robots… DOI: http://dx.doi.org/10.5772/intechopen.98391*

## **4.4 Analytical gradients**

approximation [15]. The control points in B-spline function approximation are updated in each iteration through the optimization process. Also, we used clamped B-spline which is that the starting point and end point of spline curve are matching to the control points. Following equation describes the B-spline approximation for

ð Þ¼ **<sup>t</sup>**, **<sup>p</sup>** <sup>X</sup>*<sup>m</sup>*

*j*¼0

where **t** = {*t*1, … ,*ts*} is knot vector, **p** = {*p*1, … ,*pm*} is the control points, *Nj* is the basis function of B-spline function approximation. Here, control points **P** becomes

We define the energy consuming (dynamic effort) for the given motion as the performance measure for weight lifting motion prediction and simulation. In this study, we use the joint torque square which is proportional to the mechanical energy. This mechanical energy is a reasonable criterion to minimize [4]. Then, it is

> <sup>X</sup>*w*1*τ<sup>i</sup>* <sup>2</sup> h i

where *wi* is the weighting parameters of each joint, *τ<sup>i</sup>* is the joint torque of each

Constraints are the one of the motion control way in this optimization formulation and we used minimal set of constraints to generate motion in this formulation.

The joint angle limits and torque limits for the human motion are determined based on the literature [16–19]. Zero Moment Point(ZMP) method is used for the stability condition constraint [20]. Each constraint is formulated as follows

> **<sup>q</sup>***<sup>L</sup>* <sup>≤</sup>**q**ð Þ*<sup>t</sup>* <sup>≤</sup>**q***<sup>U</sup>* (13) **<sup>τ</sup>***<sup>L</sup>* <sup>≤</sup>**τ**ð Þ*<sup>t</sup>* <sup>≤</sup> **<sup>τ</sup>***<sup>U</sup>* (14) **r***i*ðÞ¼ *t* ~**r***<sup>i</sup>* (15)

*zzmp* ∈ FSR, *xzmp* ∈FSR (16)

joint. The joint torque is calculated from above dynamics equation (Eq. (10)).

formulated as objective function in optimization formulation as follows:

*f* ¼ ð*T* 0

*N <sup>j</sup>*ð Þ**t** *p <sup>j</sup>* (11)

*dt* (12)

*qi*

**4.2 Performance measure and objective function**

the joint angle profiles *qi*.

*Collaborative and Humanoid Robots*

optimization variables.

**4.3 Constraints**

The list of constraints are as follows:

3.Foot position and hand position

1. Joint angle limits

2. Joint torque limits

4. Stability condition

accordingly:

**164**

The analytical gradients of each constraints and objective function are provided to the optimization solver. These analytical gradients improve the accuracy as well as calculation time for convergence in optimization process. Analytical gradients for constraints and objective function are as follows:

*∂C ∂x*<sup>1</sup> <sup>¼</sup> *<sup>∂</sup><sup>C</sup> ∂q*1 *∂q*1 *∂x*<sup>1</sup> þ *∂C ∂q*2 *∂q*2 *∂x*<sup>1</sup> þ ⋯ þ *∂C ∂qn ∂qn ∂x*<sup>1</sup> þ *∂C ∂q*\_ 1 *∂q*\_ 1 *∂x*<sup>1</sup> þ ⋯ þ *∂C ∂q*\_ *n ∂q*\_ *n ∂x*<sup>1</sup> þ ⋯ þ *∂C ∂*€*qn ∂*€*qn ∂x*<sup>1</sup> *∂C ∂x*<sup>2</sup> <sup>¼</sup> *<sup>∂</sup><sup>C</sup> ∂q*1 *∂q*1 *∂x*<sup>2</sup> þ *∂C ∂q*2 *∂q*2 *∂x*<sup>2</sup> þ ⋯ þ *∂C ∂qn ∂qn ∂x*<sup>2</sup> þ *∂C ∂q*\_ 1 *∂q*\_ 1 *∂x*<sup>2</sup> þ ⋯ þ *∂C ∂q*\_ *n ∂q*\_ *n ∂x*<sup>2</sup> þ ⋯ þ *∂C ∂*€*qn ∂*€*qn ∂x*<sup>2</sup> ⋮ *∂C ∂xm* <sup>¼</sup> *<sup>∂</sup><sup>C</sup> ∂q*1 *∂q*1 *∂xm* þ *∂C ∂q*2 *∂q*2 *∂xm* þ ⋯ þ *∂C ∂qn ∂qn ∂xm* þ *∂C ∂q*\_ 1 *∂q*\_ 1 *∂xm* þ ⋯ þ *∂C ∂q*\_ *n ∂q*\_ *n ∂xm* þ ⋯ þ *∂C ∂*€*qn ∂*€*qn ∂xm* (17)

$$\frac{\partial \mathfrak{r}}{\partial \mathfrak{x}\_{k}} = \frac{\partial \mathfrak{r}}{\partial q\_{j}} \frac{\partial q\_{j}}{\partial \mathfrak{x}\_{k}} + \frac{\partial \mathfrak{r}}{\partial \dot{q}\_{j}} \frac{\partial \dot{q}\_{j}}{\partial \mathfrak{x}\_{k}} + \frac{\partial \mathfrak{r}}{\partial \ddot{q}\_{j}} \frac{\partial \ddot{q}\_{j}}{\partial \mathfrak{x}\_{k}} \tag{18}$$

where *C* represents constraints, *qi* is joint angle profile, *xi* is control point of Bspline function approximation, *τ* is joint torque profile. The analytical gradients of each constraint are calculated in the form of Eq. (17). For the analytical gradient of objective function is calculated by using Eqs. (12) and (18).

#### **5. Numerical experiments**

#### **5.1 Weight lifting motion**

**Figure 4** depicts weight lifting motion. The weight is defined as W and it is applied as and external load in dynamics equilibrium equations of motion. Weight object is virtual and only hand location is guided by position constraints. Foot is located on the ground and weight is moving up vertically in 0.15 m from 0.5 m above ground and away from heal position by 0.25 m as shown in **Figure 4**. Dynamic stability constraint which is ZMP are imposed on the foot support region area. Joint angle limits and joint torque limits are imposed based on the literature review. It is assumed that the assistive moment is applied hip joint and knee joint which axes are parallel to the horizontal axis of inertial reference frame.

#### **5.2 Simulation setting**

For the optimization, sequential quadratic programming (SQP) is adopted. The sequential quadratic programming is very effective for the large scale nonlinear constrained optimization problem. The commercial software SNOPT is used which is well known as the effectiveness for large scale nonlinear problem SQP solver [21]. Total number of variables for optimization is 330 and total number of constrains is 1,766. The weighting parameters in objective function are set to equally in current study. Numerical experiment is performed while the weights are set to 15 kg and 30 kg. Also, different assistive moments from exoskeleton robot are tested in the experiment. The postprocessing for animation and snapshot was performed using Commercial software MATLAB.

**Lifting weight**

**Table 1.**

**Figure 6.**

**167**

*Numerical experiment results for weight lifting motion.*

*The contour plot of numerical experiment results.*

**Exo assistive force (N/m)**

*DOI: http://dx.doi.org/10.5772/intechopen.98391*

**Objective function value (energy consumption)**

**Lifting weight Exo assistive**

**consumption) waist knee waist knee**

15 kg 0 0 3.3260043054–002 30 kg 0 0 6.6957529363–002

*Optimization Based Dynamic Human Motion Prediction with Modular Exoskeleton Robots…*

**force (N/m)**

**Objective function value (energy**

**Figure 4.** *Weight lifting motion with waist and knee exoskeleton robots.*
