Nonlinear System Stability and Behavioral Analysis for Effective Implementation of Artificial Lower Limb

Abstract

System performance and efficiency depends on the stability criteria. The lower limb prosthetic model design requires some prerequisites such as hardware design functionality and compatibility of the building block materials. Effective implementation of mathematical model simulation symmetry towards the achievement of hardware design is the focus of the present work. Different postures of lower limb have been considered in this paper to be analyzed for artificial system design of lower limb movement. The generated polynomial equations of the sitting and standing positions of the normal limb are represented with overall system transfer function. The behavioral analysis of the lower limb model shows the nonlinear nature. The Euler-Lagrange method is utilized to describe the nonlinearity in the field of forward dynamics of the artificial system. The stability factor through phase portrait analysis is checked with respect to nonlinear system characteristics of the lower limb. The asymptotic stability has been achieved utilizing the most applicable Lyapunov method for nonlinear systems. The stability checking of the proposed artificial lower extremity is the newer approach needed to take decisions on output implementation in the system design.

1. Introduction

Generally, nature is nonlinear as the responses of the physical systems in the applied fields also show nonlinear behavior. Some electromechanical processes show nonlinear behavior. Nonlinear controller can compensate the unwanted effects on the system. An important research field for nonlinear control systems is to design controllers to deal with uncertainty, mainly due to the unavailability of parametric information of the models and external disturbances. Hard nonlinearities such as dead-zones, hysteresis and saturation do not permit linear approximation of real-world systems. The method of approximation deals with the effects of nonlinearity. This can cause inaccuracy in the system. The effects of unwanted nonlinearities in the system should be analyzed to improve the dynamic performance of any real-time system development.

The applicability of forward dynamics [1] in the nonlinear systems such as lower limb models has a great impact on motion analysis. The application of Euler-Lagrange [2,3] method has been presented in the field of robotics for the movement analysis. The mathematical explanation of the systems is the newer approach using the Euler-Lagrange method to show the artificial motion [4]. There are many interesting applications of system performance analysis through mathematical model simulation in nonlinear domains, such as walking robots, Degree of Freedom (DoF) manipulator [5] working as an artificial limb, etc. The analysis and design of nonlinear control system stability has been discussed to give an overall idea to implement on different types of applications. This idea [6] has been implemented in the lower limb system performance analysis in the present work of the paper. The improvement of system stability [7] regarding robustness and control parameters should be taken into consideration to meet the need of the system design, which is incorporated in the biomedical field for human limb performance analysis. Asymptotic stability analysis from polynomials of the system to the Lyapunov method is investigated through a computational technique. The novel approach has been utilized for the system stability analysis [8] of the lower extremity design. Basic theoretical knowledge and stability-related methodology is discussed for the nonlinear systems [9]. Linearization of differential equations and phase portrait analysis of nonlinear systems has been shown in the nonlinear domain [10]. Exponential stability analysis and design of sampled-data nonlinear systems has been shown here [11]. The nonlinear system stability has been covered using the Lyapunov stability system. The most suitable stability checking method for nonlinear systems [12,13] is involved in the lower limb stability checking purpose in this present work. Stability criteria of nonlinear systems with associated features are discussed by Z. Roberto [14]. Information about Lyapunov stability analysis of nonlinear systems has been illustrated in the paper [15]. Overview of recent research works related to nonlinear system stability in different fields have been given [16]. Jacobian matrix formation for linearization of nonlinear systems with equilibrium point [17] detection has been discussed, which is implemented in the lower limb system stability analysis. The method of linearization, Lyapunov stability and Popov criteria for stability analysis of nonlinear systems have been mentioned [18]. Information on the describing function [19] related to nonlinear system behavioral analysis using a simulation approach has been discussed as the most interesting field to involve in the lower limb system behavioral analysis in the research work of this paper. Formalization of the limit cycle prediction depending on the nonlinear element by describing function is presented with great importance in the field of the system’s nature prediction [20]. The basics of describing function is discussed to predict the limit cycle [21]. The analysis of nonlinear systems based on the describing function, which is dependent on the different types of controllers such as proportional-Integral etc., are the base of the system characteristics determination as mentioned in the current work of research [22].

This approach to nonlinear control system stability determination has been chosen to incorporate the idea of the behavioral analysis of the system. The model of the lower limb of human beings is the main concern to be built representing the nonlinear system. Through some specific features, the implementation of system identification regarding forward dynamics and model generation is possible. The determination of stability in the nonlinear domain [23,24] preferably nourishes the idea of the planning of artificial limb [25,26] design. Until now, it has not been reported in any open journal that lower limb system stability can be analyzed using a nonlinear stability analysis method [27,28] in the field of bio robotics.

2. Forward Dynamics of the Artificial Lower Limb Movement

Forward dynamics [29] is the connection between the forces applied on a robotic mechanism and the produced accelerations, which is presented in Figure 1. Robot dynamics [30] is the application of rigid-body dynamics to robots. Forward dynamics [2] gives the forces and work out of the accelerations.

• Suppose, Joint position (Knee or ankle) of the artificial lower limb = $\varnothing$
• Joint force of the lower limb = $\tau$
• Joint velocity of the lower limb = $\dot{\varnothing }$
• Joint acceleration of the lower limb = $\ddot{\varnothing }$

The Lagrangian approach delivers the immediate reduction in the characteristic size to solve many equations. This describes the equation of motion. The Lagrange formulation is used to avoid the acceleration finding as shown in Figure 2. The Lagrangian equation can be expressed as the difference between the kinetic energy and the potential energy as shown in the Equation (1).

$$L\left(\varnothing , \dot{\varnothing } \right)=K\left(\varnothing , \dot{\varnothing } \right)-P(\varnothing )$$

(1)

where, Kinetic energy = $K\left(\varnothing , \dot{\varnothing } \right)$ and Potential energy = $P(\varnothing )$.

Euler-Lagrange equations are suitably obtained for the nonlinear behavioral analysis of the artificial movement of the lower limbs.

The Euler—Lagrange equation [3] is given below in Equation (2).

$$\mathsf{\tau }= \frac{d}{dt} \frac{\partial L}{\partial \dot{\varnothing }} -\frac{\partial L}{\partial \varnothing }$$

(2)

where,

• $\frac{d}{dt} \frac{\partial L}{\partial \dot{\varnothing }}=$time derivative of partial derivative of L with respect to $\dot{\varnothing }$
• $\frac{\partial L}{\partial \varnothing }=$partial derivative of L with respect to $\varnothing$
• $\mathsf{\tau }=$ power produced by the joints related to $\left(\varnothing , \dot{\varnothing } \right)$.

In the form of vector components, it can be written as shown in Equation (3).

$$T_i= \frac{d}{dt} \frac{\partial L}{\partial \dot{\varnothing }} -\frac{\partial L}{\partial \varnothing }$$

(3)

The joint1 refers to the ankle joint and the joint2 refers to the knee joint in Figure 3. From the above diagram, as shown in the Figure 3, Equation (4) can be written for the position of joint2 ($a_2, b_2$).

$$\left[\begin{array}{c}{a}_1\\ b_1\end{array}\right] = \left[\begin{array}{c}{x}_1\cos {\varnothing }_1\\ x_1\sin {\varnothing }_1\end{array}\right]$$

(4)

where, $x_{1 }$ = length of the limb from joint1 ($a_1, b_1$) to joint2.

The derivative of $\left(a_1, b_1\right)$ is expressed in Equation (5).

$$\left[\begin{array}{c}{\dot{a}}_1\\ {\dot{b}}_1\end{array}\right] = \left[\begin{array}{c}{x}_1\cos {\varnothing }_1\\ x_1\sin {\varnothing }_1\end{array}\right] \dot{{\varnothing }_1}$$

(5)

Similarly, it can be written for the position of end-effector ($a_e, b_e$) as shown in the Equation (6)

$$\left[\begin{array}{c}{a}_2\\ b_2\end{array}\right] =\left[ \begin{array}{cc}{x}_1\cos {\varnothing }_1 +& x_1\cos \left(\begin{array}{cc}{\varnothing }_1+& {\varnothing }_2\end{array}\right)\\ x_1\sin {\varnothing }_1+& x_1\sin \left(\begin{array}{cc}{\varnothing }_1+& {\varnothing }_2\end{array}\right)\end{array}\right]$$

(6)

This equation is for the position and velocity of the model as shown in Equation (7).

$$\left[\begin{array}{c}{\dot{a}}_2\\ {\dot{b}}_2\end{array}\right]= \left[\begin{array}{cc}-x_1\sin {\varnothing }_1-x_2\sin \left(\begin{array}{cc}{\varnothing }_1+& {\varnothing }_2\end{array}\right)& -x_2\sin \left(\begin{array}{cc}{\varnothing }_1+& {\varnothing }_2\end{array}\right)\\ x_1\cos {\varnothing }_1 +x_2\cos \left(\begin{array}{cc}{\varnothing }_1+& {\varnothing }_2\end{array}\right) & x_2\cos \left(\begin{array}{cc}{\varnothing }_1+& {\varnothing }_2\end{array}\right)\end{array}\right]\left[\begin{array}{c}\dot{{\varnothing }_1}\\ \dot{{\varnothing }_2}\end{array}\right]$$

(7)

where, $x_{2 }$ = length of the limb from joint2 to end-effector.

Now, the Kinetic Energy for link1 and link2 are shown in the Equations (8) and (9) respectively.

$$K.E_1= \frac{1}{2} m_1{v}_1^2=\frac{1}{2 } m_1\left({\dot{a}}_1^2+{\dot{b}}_1^2\right)= \frac{1}{2} x_1{\dot{\varnothing }}_1^2$$

(8)

where,

• $m_1$ = mass of the joint2,
• $v_1=$ velocity of the joint2.

  $$\begin{gathered} K.E_2= \frac{1}{2} m_1{v}_1^2=\frac{1}{2 } m_1\left({\dot{a}}_1^2+{\dot{b}}_1^2\right) \\ = \frac{1}{2} m_2\left(\left(x_{1 }^2+2x_1{x}_2\cos {\varnothing }_2+ x_{2 }^2 \right) {\dot{\varnothing }}_1^2+2\left(x_{2 }^2+x_1{x}_2\cos {\varnothing }_2\right)\dot{{\varnothing }_1}\dot{{\varnothing }_2} \\ + x_{2 }^2{\dot{\varnothing }}_2^2\right) \end{gathered}$$

  (9)

Now, the potential energy for link1 and link2 are as given in Equations (10) and (11),

$$P.E_1= m_{1}g{y}_1= m_{1}g{x}_1\sin {\varnothing }_1$$

(10)

$$P.E_2= m_{2}g{y}_2= m_{2}g(x_1\sin {\varnothing }_1+ x_2\sin ({\varnothing }_1+ {\varnothing }_2))$$

(11)

Now,

$$L\left(\varnothing , \dot{\varnothing } \right)= {\displaystyle \sum }_{i=1}^2(K_i- P_i)$$

(12)

$${\tau }_i= \frac{d}{dt} \frac{\partial L}{\partial \dot{{\varnothing }_i}}- \frac{\partial L}{\partial \varnothing }$$

(13)

where, $i=1, 2$.

One component of Lagrangian is $L_{comp}$ as shown in the Equation (14).

Therefore,

$$L_{comp}= m_2 x_1{x}_2\dot{{\varnothing }_1}\dot{{\varnothing }_2} \cos {\varnothing }_2$$

(14)

The effect of $L_{comp}$ on the 2nd joint torque is as given below in Equation (15).

$${\tau }_2{}_{comp}= \frac{d}{dt}\frac{\partial L_{comp}}{\partial \dot{{\varnothing }_2}}-\frac{\partial L_{comp}}{\partial {\varnothing }_2}$$

(15)

$$\begin{array}{l}= \frac{d}{dt} \left(m_2 x_1{x}_2\dot{{\varnothing }_1} \cos {\varnothing }_2\right)\\ = m_2 x_1{x}_2\ddot{{\varnothing }_1} \cos {\varnothing }_2- m_2 x_1{x}_2{\dot{\varnothing }}_1\dot{{\varnothing }_2} \sin {\varnothing }_2+ m_2 L_1{L}_2\dot{{\varnothing }_1}\dot{{\varnothing }_2} \sin {\varnothing }_2\\ = m_2 x_1{x}_2\ddot{{\varnothing }_1} \cos {\varnothing }_2\end{array}$$

(16)

Now torque for first joint is ${\tau }_1$

$$\begin{array}{ll}{\tau }_1= (m_1{x}_1^2+ & m_2\left(x_1^2+2x_1{x}_2\cos {\varnothing }_2+ x_2^2\right)){\ddot{\varnothing }}_1+m_2\left(x_1{x}_2\cos {\varnothing }_2+ x_2^2\right){\ddot{\varnothing }}_2\\ & -m_2 x_{1 }{x}_{2 }\sin {\varnothing }_2\left(2\dot{{\varnothing }_1}\dot{{\theta }_2}+ {\dot{\varnothing }}_2^2\right)+\left(m_1+ m_2\right)x_{1}g\cos {\varnothing }_1\\ & +m_{2}g{x}_2 \cos \left({\varnothing }_1+ {\varnothing }_2\right)\end{array}$$

(17)

Now torque for second joint is ${\tau }_2$

$$\begin{gathered} {\tau }_2= m_2\left(x_1{x}_2\cos {\varnothing }_2+ x_2^2\right)){\ddot{\varnothing }}_1+ m_2{x}_2^2 {\ddot{\varnothing }}_2- m_2 x_1{x}_2 {\dot{\varnothing }}_1^2\sin {\varnothing }_2 \\ + m_{2}g{x}_2 \cos \left({\varnothing }_1+ {\varnothing }_2\right) \end{gathered}$$

(18)

Some terms are not dependent on the joint acceleration and velocity. So, the torque can be expressed as shown in Equation (19).

$$\tau =M(\varnothing )\ddot{\varnothing }+C\left(\varnothing , \dot{\varnothing }\right)+g(\varnothing )$$

(19)

where,

• $M(\varnothing )\ddot{\varnothing }=n\times n$ mass matrix shown in Equation (20)
• $C\left(\varnothing , \dot{\varnothing }\right)$ = velocity product term shown in Equation (21)
• $g(\varnothing )=$ gravity term, i.e., gravity due to springs in the mechanical model for potential energy is shown in Equation (22).

  $$M(\varnothing )\ddot{\varnothing }= \left[\begin{array}{cc}{m}_1{x}_1^2+ m_2\left(x_1^2+2x_1{x}_2\cos {\varnothing }_2+ L_2^2\right)& m_2\left(x_1{x}_2\cos {\varnothing }_2+ x_2^2\right)\\ m_2\left(x_1{x}_2\cos {\varnothing }_2+ x_2^2\right)& m_2{x}_2^2\end{array}\right]$$

  (20)

  $$C\left(\varnothing , \dot{\varnothing }\right)= \left[\begin{array}{c}- m_2 x_1{x}_2 \sin {\varnothing }_2 \left(2\dot{{\varnothing }_1}\dot{{\varnothing }_2}+ {\dot{\varnothing }}_2^2\right)\\ m_2 x_1{x}_2 {\dot{\varnothing }}_1^2 \sin {\varnothing }_2\end{array}\right]$$

  (21)

  $$g(\varnothing )=\left[\begin{array}{c}(m_1+ m_2) x_{1}g\cos {\varnothing }_1+m_{2}g{x}_2 \cos \left({\varnothing }_1+ {\varnothing }_2\right)\\ m_{2}g{x}_2 \cos \left({\varnothing }_1+ {\varnothing }_2\right)\end{array}\right]$$

  (22)

The overall acceleration depends not only on the joint accelerations but also on the products of the joint velocities. Joint coordinates are non-inertial. These velocity product terms appear because joint coordinates are non-inertial.

3. Experimental Methods

3.1. Stability Analysis of the Non-Linear System

The dynamic behavior of the artificial lower extremity model shows nonlinearity. The postures such as sitting and standing conditions of the lower part of human body exhibits the nonlinear characteristics. The control system [31,32,33] involves the mathematical relationship related to the nonlinear properties such as equilibrium point, saturation, dead zone, hysteresis, describing function, etc. Nonlinear systems may not produce a specific solution for all times but gives a decision maker feature [34,35] of the analyzed system. Through phase portrait analysis and the Lyapunov method, the system stability [36,37] of the designed system definition is feasible. Thus, system design method for nonlinear control systems are application-specific [38]. The describing function is the well-established process to analyze the nonlinear function [39] with an approximation, and to determine the stability of a nonlinear system. It is basically the approximated extended version of frequency domain analysis. For a higher degree of nonlinearity, the describing function is used but no information about time domain such as steady state error is given by this. Phase plane analysis can be used, but the order of the system should be less than second order as it is a graphical system. Sometimes it is called a sinusoidal describing function.

3.2. Workflow of the System Stability Analysis

The workflow of the analysis for lower limb system [40,41] stability with respect to linear and nonlinear systems is shown in Figure 4, Figure 5 and Figure 6, respectively.

The polynomial equations generated from two different postures of the lower limb, such as sitting and standing, are considered as the primary conditions of the artificial system stability analysis. The transfer function is generated from the polynomials to achieve the mathematical model [42] of the mentioned system. The closed loop transfer function is generated with unity feedback to gather information about the controller, performance and robustness parameters of the system, which shows stability as shown in Figure 4.

The characteristics equation is generated from the polynomials of the system at two different positions of the lower limb. The roots are determined using a computational technique. In the nonlinear domain, using Eigen vectors and the generated matrix, the phase portrait analysis has been incorporated. The system shows asymptotic stability in the nonlinear domain for both the roots of the characteristics equation as shown in Figure 5. For better and stronger decision making on the stability of the system, the Lyapunov stability [43] analysis method has been performed with the equilibrium point detection. This method also shows asymptotic stability where the system is stable and roots are in a complex domain as shown in Figure 6.

4. Results and Discussions

In this present research work, a 3 DoF lower limb has been analyzed for the artificial electro-mechanical system development. The conditions of movement of the human body such as sitting and standing are the main concern in this discussion. Non applicability of the superposition and homogeneity property gives a greater emphasis on the analysis and design of the control systems in the nonlinear domain [15]. The controllability and observability [44] cannot be determined simply based on rank tests. The concept of system stability [45] and checking using Phase portrait and Lyapunov stability [17] has been incorporated to make it a complete evaluation. Polynomials for the two positions of the lower limb (sitting and standing) as shown in Figure 7 and Figure 8 are presented below in Equations (23) and (24), respectively [5]:

$${\mathrm{f}}_1\left(\mathrm{t}\right){ = \mathrm{t}}^2+25.6354\mathrm{t}+92.4473$$

(23)

See Appendix A.

$${\mathrm{f}}_2\left(\mathrm{t}\right){ = \mathrm{t}}^2+19\mathrm{t}+18$$

(24)

See Appendix A.

From Equations (23) and (24), the open loop transfer function is obtained, as given below in the Equation (25):

$$G\left(s\right)= \frac{2.815\times {10}^{14} s^2+7.216\times {10}^{15} s+2.602\times {10}^{16} }{2.815\times {10}^{14} s^2+5.348\times {10}^{15} s+5.067\times {10}^{15} }$$

(25)

The output graphical presentation of the step response of the open loop transfer function of the system has been shown in the Figure 9.

The closed loop transfer function is given below in Equation (26) derived from Equation (25):

$$G\left(s\right)= \frac{2.815\times {10}^{14} s^2+7.216\times {10}^{15} s+2.602\times {10}^{16}}{5.629\times {10}^{14} s^2+1.256\times {10}^{15 }s+3.109\times {10}^{16}}$$

(26)

The output graphical presentation of the closed loop transfer function of the system has been shown in Figure 10. The PID (Proportional-Integral-Derivative) controller has been used to tune the closed loop transfer function of the lower limb system. The controller parameters of the closed loop transfer function of the system are shown in

Table 1. Controller parameters of the closed loop transfer function of the system.

Proportional Constant (kp)	Integral Constant (ki)	Derivative Constant (kd)	Stability Status
0.92214	2.8732	0	Stable

and the performance and robustness parameters of the closed loop transfer function of the system are shown in

Table 2. Performance and Robustness parameters of the closed loop transfer function of the system.

Rise Time (s)	Settling Time (s)	Overshoot (%)	Peak (s)	Stability Status
0.626	7.77	10.3	1.1	Stable

. It is seen that the system is underdamped in Figure 10. The status related to the stability analysis of the mentioned system is prominent from the graphical presentation in Figure 10. The red line is representing the lower limb closed loop transfer function presented in Equation (26) in this paper. Comparing the reference signal represented through the blue line in Figure 10 and the lower limb system output presented in Equation (26), it is observed that the red line, the lower limb system output signal, has reached the stable condition with respect to the blue line at the settling point of the graph. The representations of the red and the blue lines are shown in the legends in Figure 10.

In Figure 10, Table 1 and Table 2 are representing the tuned graphical analysis of the lower limb transfer function presented in Equation (26), using the traditional PID controller. The output performance is not very satisfactory here. So, an advanced controller, the “Control System Designer”, is preferred in order to show the better-tuned output presented in Figure 11 of the lower limb transfer function.

Therefore, the decision making regarding the system stability is possible observing the whole performance of the advanced controller output. In

Table 3. Performance and robustness parameters of the closed loop transfer function of the system.

Rise Time (s)	Settling Time (s)	Peak Amplitude (s)	Stability Status
0.00191	5.05	0.965	Stable

, the related characteristic parameters are shown. This is clearly shown in that the peak value and the settling time is reduced to a certain level.

Modern technology requires advanced control laws to meet the design specifications for humanoid limb design, thus highlighting the position of nonlinear control systems with recent advancements in the artificial lower limb design. The application area in system design and analysis means that many technologies are developed based on this advanced nonlinear control system stability analysis method. The non-linear system stability analysis is processed as discussed below.

The generated matrix from Equations (23) and (24) is as given below:

$$A=\left[\begin{array}{cc}1& 25.6354\\ 1& 19\end{array}\right]$$

(27)

Then

$$\begin{gathered} A- \lambda I= \left[\begin{array}{cc}1& 25.6354\\ 1& 19\end{array}\right]- \left[\begin{array}{cc}\lambda & 0\\ 0& \lambda \end{array}\right] \\ = \left[\begin{array}{cc}1-\lambda & 25.6354\\ 1& 19-\lambda \end{array}\right] \end{gathered}$$

(28)

where, $\lambda =eigen value, I=Identity matrix$.

Suppose $X= A- \lambda I$

Now,

$$\det \left(X\right)= {\lambda }^2-20\lambda -6.6354$$

(29)

$\det \left(X\right)$ is called the characteristic equation of the matrix A. Since the matrix $A- \lambda I$ is singular, $\det \left(X\right) = 0$.

Therefore,

$${\lambda }^2-20\lambda -6.6354=0$$

(30)

The roots of the equation are $\lambda =10 \pm 10.3264$.

Therefore, the two roots are ${\lambda }_1=20.3264 \mathrm{and} {\lambda }_2=-0.3264$.

For ${\lambda }_1=20.3264$ from matrix X, the obtained value is given below:

$$\left[\begin{array}{cc}-19.3264& 25.6354\\ 1& -1.3264\end{array}\right]$$

Suppose, $B= \left[\begin{array}{cc}-19.3264& 25.6354\\ 1& -1.3264\end{array}\right]$.

Now solving for $B\overline{X}=0$,

$$\left[\begin{array}{cc}-19.3264& 25.6354\\ 1& -1.3264\end{array}\right]\left[\begin{array}{c}{X}_1\\ X_2\end{array}\right] =\left[\begin{array}{c}0\\ 0\end{array}\right]$$

(31)

Using the Gauss elimination method, the obtained simplified equation is

$$25.6354X_2-19.3264X_1=0$$

(32)

Or,

$$\frac{X_2}{X_1} = \frac{19.3264}{25.6354}$$

For ${\lambda }_2= -0.3264$ from matrix X, obtained value is given below:

$$\left[\begin{array}{cc}1.3264& 25.6354\\ 1& 19.3264\end{array}\right]$$

$$\mathrm{Suppose} C= \left[\begin{array}{cc}1.3264& 25.6354\\ 1& 19.3264\end{array}\right]$$

Now solving for $C\overline{X}=0$,

$$\left[\begin{array}{cc}1.3264& 25.6354\\ 1& 19.3264\end{array}\right]\left[\begin{array}{c}{X}_1\\ X_2\end{array}\right] =\left[\begin{array}{c}0\\ 0\end{array}\right]$$

(33)

Using the Gauss elimination method, the obtained simplified equation is

$$19.3264X_2+X_1=0$$

Or,

$$\frac{X_2}{X_1} = \frac{1}{-19.3264}$$

(34)

Now from Equation (27) and using the roots ${\lambda }_1=20.3264 \mathrm{and} {\lambda }_2=-0.3264$.

The Eigen vectors are given by

$$\stackrel{\rightharpoonup }{V_1}= \left(\begin{array}{c}25.6354\\ 19.3264\end{array}\right) and \stackrel{\rightharpoonup }{V_2}= \left(\begin{array}{c}-19.3264\\ 1\end{array}\right)$$

(34a)

The produced graphical phase portrait [10] plot is given below in Figure 11.

This output shows asymptotic stability because for Eigen vectors

$$\stackrel{\rightharpoonup }{V_1}= \left(\begin{array}{c}25.6354\\ 19.3264\end{array}\right) and \stackrel{\rightharpoonup }{V_2}= \left(\begin{array}{c}-19.3264\\ 1\end{array}\right)$$

(34b)

the first phase plot is going outward or away from the origin of the plane and the second phase plot is coming towards the origin of the plane as shown in Figure 12. These conditions of phase portrait plot proves an asymptotic stability. The phase plot lines are presented through red lines in Figure 12.

The theory of Lyapunov stability is a standard theory and one of the most important mathematical tools in the analysis of non-linear systems. Now, for further stability checking, the Lyapunov stability method is strongly recommended for nonlinear systems. The followed steps are given below for stability checking:

• Step 1: Representation of polynomials in characteristics equation form.
• Step 2: Determination of the equilibrium point of the system.
• Step 3: Generation of Jacobian matrix to implement in the characteristics equation.
• Step 4: Determination of roots from the characteristics equation.
• Step 5: Verifying the stability conditions from the roots.
• Step 6: Declaration of Asymptotic stability condition.

Now the process of Lyapunov stability method analysis is as shown below:

Let

$$y= x_1, \dot{y}= \dot{x_1}= x_2, \ddot{y}= \dot{x_2}$$

(35)

From Equation (23) it can be written as

$$y^2+25.6354y+92.4473=0$$

(36)

From Equation (24) it can be written as

$$y^2+19y+18=0$$

(37)

After applying a derivative process in the Equations (36) and (37) respectively, the obtained Equations are as given below:

$$2y+25.6354=0$$

(38)

$$2y+19=0$$

(39)

Now, putting in values from Equation (35), the obtained Equations are as given below:

$$2x_1+25.6354=0$$

(40)

$$2x_1+19=0$$

(41)

After applying the derivative process again in the Equations (40) and (41) respectively, the obtained Equations are as given below:

$$2\dot{x_1}=0$$

(42)

$\dot{x_1}=0$ or it can be written as

$$x_2=0$$

(43)

Now, to find the equilibrium point, it is similar as $\dot{x_2}=0$ or it can be written as

$$x_1=0$$

(44)

Therefore, the equilibrium point is at (0, 0) from Equations (43) and (44).

Now, the Jacobian matrix can be formed using the formula as given below in the Equation (45):

$$J= \left[\begin{array}{cc}\frac{\partial f}{\partial x_1}& \frac{\partial f}{\partial x_2}\\ \frac{\partial g}{\partial x_1}& \frac{\partial g}{\partial x_2}\end{array}\right]= \left[\begin{array}{cc}0& 1\\ -1& -1\end{array}\right]$$

(45)

Since,$\lambda I-A=0$ where, $A=J$

Or, $\left|\left[\begin{array}{cc}\lambda & 0\\ 0& \lambda \end{array}\right]- \left[\begin{array}{cc}0& 1\\ -1& -1\end{array}\right]\right|=0$

Or, ${\lambda }^2+ \lambda +1=0$

Therefore, the roots are ${\lambda }_1, {\lambda }_2= -\frac{1}{2} \pm j\frac{1}{\sqrt{3}}$.

The abovementioned system has roots containing real and imaginary parts. Therefore, the system can be declared as asymptotically stable. The system performance efficiency is checked through this method.

This method is involved to ensure stability of the electro-mechanical system in the field of behavioral quality improvement of the designed artificial model. According to the above results shown through the mathematical models, the nonlinearity nature present in the lower limb system is a dead zone combined with a saturation condition. Figure 13 represents the lower limb system characteristics of overall nonlinearity.

The nonlinear simulation of the lower limb system mentioned in Equation (26) is simulated using the dead zone and saturation functions to show satisfactory output in this paper for the overall simulation of the system. The graphical presentations are presented in Figure 14 and Figure 15, respectively.

The graphical scenario depicts the lower limb system performance regarding nonlinearity. In Figure 14, initially the signal has started to show the activities but after a certain time it has reached the dead zone or the “no action” condition. In Figure 15, initially the signal has started to show output but after a certain period it has reached to saturation condition, being stable at a fixed level. These dead zone and saturation conditions are the specified nonlinear characteristics of the lower limb function only.

The artificial limb design needs to be developed on the basis of knowledge-based system development. This process demands some expertise in the domain of the embedded system. This also requires the intelligent incorporation of an input–output database for the trained functional implementation of humanoid movement with an artificial lower extremity. Modern approaches for nonlinear control systems can be represented as an intelligent control. Various methods can be utilized to develop an automated system using neural networks, machine learning, and fuzzy logic, etc.

5. Conclusions

The nonlinear control systems have gained a challenging position in the technological development of the biomedical application field. In this recent work, the forward dynamics of the robotic system with the DoF concept has been explained using the Euler-Lagrange method. The idea of the mathematical model representation has been demonstrated with the polynomial equation generation from the experimental setup of the lower limb model. The comparative analysis has been incorporated for the system performance analysis using the traditional PID controller and an advanced control system designer. In future, there are lots of opportunities to use a more advanced application of controllers. The system performance stability is checked in two ways. It is found to be asymptotically stable using phase portrait analysis and the most important Lyapunov stability method. The lower limb shows nonlinear behavior. The system remains in the dead zone initially and after a certain period it starts operating. Then it reaches to the saturation level. The overall nonlinear behavior of the lower limb system has been shown here. The role of the describing function is the most important factor to determine the system nature selection. In future, the control system simulation approach can be implemented to get a more informative system identification view of the artificial system design.