Modeling, Identification, and Control of a Discrete Variable Stiffness Actuator (DVSA)

A branch of robotics, variable impedance actuation, along with one of its subfields variable stiffness actuation (VSA) targets the realization of complaint robotic manipulators. In this paper, we present the modeling, identification, and control of a discrete variable stiffness actuator (DVSA), which will be developed for complaint manipulators in the future. The working principle of the actuator depends on the involvement of series and parallel springs. We firstly report the conceptual design of a stiffness varying mechanism, and later the details of the dynamic model, system identification, and control techniques are presented. The dynamic parameters of the system are identified by using the logarithmic decrement algorithm, while the control schemes are based on linear quadratic control (LQR) and computed torque control (CTC), respectively. The numerical simulations are performed for the evaluation of each method, and results showed the good potentialities for the system. Future work includes the implementation of the presented approach on the hardware.


Introduction
Safety is the most important element during human-robot interactions in smart manufacturing and domestic scenarios [1][2][3].Bio-inspired robotics is inspired by soft manipulation with improvement towards human-robot interaction [4][5][6].Remarkable progress has been made towards complaint robots, and human-robot interaction (HRI) has become the most important feature in the robot design and development.Although far from reaching the performance of their biological counterpart, robots have been offering a wide range of applications in many different fields, from medicine to industry [7].Human-robot interaction (HRI) aims at developing hardware (actuation) and software (algorithms) to allow more natural and effective communication and interaction between humans and robots [2, 8,9].The complaint actuators are growing in the scientific community because of the performance limitations of traditional stiff actuation approaches in terms of safety, energy efficiency, and the ability to interact with the environment [10].
The possible advantages that can be obtained with complaint actuators are very well known, but most of the robots in industry still lack these features.Towards this direction, many complaint actuation systems have been proposed in the last two decades.The research on complaint manipulators was driven by the need for more safety and robustness in human-robot interactions [1,11].Robustness and soft interaction are mainly achieved by either actively regulating the compliance of robotic joints [1,2] or tuning the intrinsic softness, acting on the passive characteristics of the robot body ware [12][13][14].Incorporating complaint elements in actuators comes with many challenges from the mechatronic and control points of view.For this reason, the use of complaint actuators in robotic manipulators is still limited and challenging [11].
Many designs and prototypes of variable stiffness actuators are presented by the research community.Some remarkable examples are VSA-I [15], VSA-II [16], AMASC [17], and the biologically inspired joint stiffness control mechanism [6].An alternative approach of varying the stiffness was achieved through the principle of a lever mechanism [6,[18][19][20][21][22][23].Another new method of stiffness variation was proposed through the addition and subtraction of elastic element involvement in actuator mechanism.The arrangements of elastic elements can be done in series [24] or in parallel, which are also called series-parallel elastic actuators (SPEAs).The approach of series-parallel arrangement is efficient in terms of energy and peak torque [25][26][27].Based on the concept of SPEAs, we firstly proposed a passive version of the discrete variable stiffness actuator (DVSA) (see Figure 1) [28].In our all previous works, we mainly focused on the passive versions of the actuators.In this regard, we proposed various design topologies and optimized their design specifications according to the intended applications ranging from industrial to haptics and medical.Following this motivation, in [24], we presented an approach of realizing the arm exoskeleton (TELEXOS-I) for applications such as haptic-teleoperation and rehabilitation.The key motivation behind our design topology is the need for instantaneous switching between stiffness levels when applied for virtual reality (VR) or remote environment (RE) stiffness mapping applications, in addition to the need of low-energy consumption.Moreover, in [29], we presented a case scenario of using another developed passive discrete variable stiffness joint (pDVSJ-II) in a remote palpation application.Palpation is defined as the application of the fingers with light pressure to the surface of the human body for the purpose of determining the condition of the parts beneath during physical diagnosis.Palpation with a robotic device, or the artificial tactile exploration of soft viscoelastic and nonhomogeneous objects, is an important area of study for various fields including medical, VR, or teleoperation applications.
joints [1,2] or tuning the intrinsic softness, acting on the passive characteristics of the robot body ware [12][13][14].Incorporating complaint elements in actuators comes with many challenges from the mechatronic and control points of view.For this reason, the use of complaint actuators in robotic manipulators is still limited and challenging [11].
Many designs and prototypes of variable stiffness actuators are presented by the research community.Some remarkable examples are VSA-I [15], VSA-II [16], AMASC [17], and the biologically inspired joint stiffness control mechanism [6].An alternative approach of varying the stiffness was achieved through the principle of a lever mechanism [6,[18][19][20][21][22][23].Another new method of stiffness variation was proposed through the addition and subtraction of elastic element involvement in actuator mechanism.The arrangements of elastic elements can be done in series [24] or in parallel, which are also called series-parallel elastic actuators (SPEAs).The approach of series-parallel arrangement is efficient in terms of energy and peak torque [25][26][27].Based on the concept of SPEAs, we firstly proposed a passive version of the discrete variable stiffness actuator (DVSA) (see Figure 1) [28].In our all previous works, we mainly focused on the passive versions of the actuators.In this regard, we proposed various design topologies and optimized their design specifications according to the intended applications ranging from industrial to haptics and medical.Following this motivation, in [24], we presented an approach of realizing the arm exoskeleton (TELEXOS-I) for applications such as haptic-teleoperation and rehabilitation.The key motivation behind our design topology is the need for instantaneous switching between stiffness levels when applied for virtual reality (VR) or remote environment (RE) stiffness mapping applications, in addition to the need of low-energy consumption.Moreover, in [29], we presented a case scenario of using another developed passive discrete variable stiffness joint (pDVSJ-II) in a remote palpation application.Palpation is defined as the application of the fingers with light pressure to the surface of the human body for the purpose of determining the condition of the parts beneath during physical diagnosis.Palpation with a robotic device, or the artificial tactile exploration of soft viscoelastic and nonhomogeneous objects, is an important area of study for various fields including medical, VR, or teleoperation applications.Now that the ability of passive actuators has been successfully demonstrated to change their stiffness at different levels using the proposed topologies, we are in the process of developing a compact and active version (with a motor).In [30], we introduced the concept of the active binary-controlled variable stiffness actuator (BcVSA), which is an example that falls under the category of DVSA.One of the benefits of our approach is the selection of a desired level of stiffness without the need for succession involvement, which improves the response time in switching the level of stiffness.In order to control the actuator, we presented the model predictive control (MPC) technique.The control design and identification of the dynamic parameters are quite challenging tasks and play a great role in the efficiency of the actuator.In fact, system identification is an important and crucial problem in a control system.Any uncertainty in the dynamic parameters will eventually affect the performance of the system if the controller does not address it.In the state of the art, several algorithms for identifying Now that the ability of passive actuators has been successfully demonstrated to change their stiffness at different levels using the proposed topologies, we are in the process of developing a compact and active version (with a motor).In [30], we introduced the concept of the active binary-controlled variable stiffness actuator (BcVSA), which is an example that falls under the category of DVSA.One of the benefits of our approach is the selection of a desired level of stiffness without the need for succession involvement, which improves the response time in switching the level of stiffness.In order to control the actuator, we presented the model predictive control (MPC) technique.The control design and identification of the dynamic parameters are quite challenging tasks and play a great role in the efficiency of the actuator.In fact, system identification is an important and crucial problem in a control system.Any uncertainty in the dynamic parameters will eventually affect the performance of the system if the controller does not address it.In the state of the art, several algorithms for identifying the dynamic parameters have been introduced to enhance the performance of the system.One of these algorithms is the logarithmic decrement algorithm, and it is introduced with suitable measuring techniques in [31,32].Other publications have concentrated on identifying dynamic parameters for various types of robots based on the logarithmic decrement algorithm, which is presented in [33,34].
In this paper, we use the same concept of the BcVSA for varying stiffness but propose system identification and different control techniques.Thus, the contribution of this paper lies in identifying the dynamic parameters of our active actuator and implementing new suitable control techniques.The aim is to evaluate the ability of the proposed controllers using the identified dynamic parameters to control our actuator and to test the controllers robustness in different cases while changing the actuator's stiffness in two cases-offline (when the system is stopped) and online (while the system is running).Both methods are simulated using actuator dynamics under the effect of various combinations of active (engaged) and inactive (disengaged) clutches.In particular, we have used the logarithmic decrement algorithm to perform system identification and have evaluated two control schemes using the identified parameters in the model.The control techniques implemented are the linear quadratic regulator (LQR) [12,13] and the computed torque control (CTC) [14,[35][36][37].We performed numerical simulations to evaluate the effectiveness of our approach.The purpose of the numerical simulations was to choose the proper control scheme for experiments on the hardware of the system in the near future.
The rest of the paper is organized as it follows.The conceptual design and the stiffness modeling are presented in Section 2. Section 3 describes the system identification by using a logarithmic decrement algorithm.Section 4 presents the two control schemes based on linear quadratic control (LQR) and computed torque control (CTC), respectively.The simulation results are reported in Section 5.Moreover, conclusions and future work are detailed in Section 6.

Design and Stiffness Model of the Discrete Variable Stiffness Actuator
The motivation for the proposed actuator lies in the series-parallel variable stiffness mechanism, as shown in Figure 2, which shows the conceptual 3D CAD representation and schematic illustration of the proposed discrete variable stiffness actuator.A DC motor is coupled with the sun gear of the first planetary gear set, which is engaged with its plane gears.The torsional springs are serially connected to the inline clutches at one end and, to another planet, gears of the second planet gear set.The sun gear of this planet gear stage drives the load shaft (see Figure 2a).If any of the clutches are active, the actuator torque is transmitted to load through the corresponding spring.Here, a stiffness bit consists of a clutch, a torsional spring, and a planet gear, which transmit the power from actuator to load.Several stiffness bits can be connected in parallel to achieve the overall stiffness of the actuator.
The novelty of our discrete variable stiffness actuator lies in the design topology (see Figure 2a).Three springs with different stiffness values each are serially attached to clutches.The clutches are connected to the planet gears which in turn are attached to the sun gear.The other end of each spring is connected to a planet gear that contributes part of the resultant torque on the other sun gear.The load shaft is attached to the sun gear of the system.If any clutch is disengaged, the corresponding spring will be detached.The achievable levels of stiffness in binary representation are reported in Table 1, where "0" represents inactive stiffness and "1" active stiffness.The clutch reduces the stiffness switching time.Moreover, scalability can be achieved either by changing the value of the lower stiffness (K0) or by adding extra stiffness-bits.

Stiffness Bit 2 2 (4K0) Stiffness Bit 2 1 (2K0) Stiffness Bit 2 0 (K0) Stiffness Level
In order to derive the stiffness model, we can use the actuator's kinematics.If the motor's shaft is clamped (blocked from motion) and an external load is applied, the external torque will rotate the sun gear.In the scenario where no stiffness bit is active, the spring will rotate freely, while in the scenario where any of the stiffness bits are active, the load torque will be transmitted through the corresponding spring, which deflects and produces a counter torque on the human hand.
The resultant torque can be written as follows: where  ,  , and  are the torque transmitted by each planet gear to the sun gear, which can be further written as In order to derive the stiffness model, we can use the actuator's kinematics.If the motor's shaft is clamped (blocked from motion) and an external load is applied, the external torque will rotate the sun gear.In the scenario where no stiffness bit is active, the spring will rotate freely, while in the scenario where any of the stiffness bits are active, the load torque will be transmitted through the corresponding spring, which deflects and produces a counter torque on the human hand.
The resultant torque can be written as follows: where τ 0 , τ 1 , and τ 2 are the torque transmitted by each planet gear to the sun gear, which can be further written as where: • N: the gear ratio; The total torque (τ Σ ) and the total torsional stiffness (K Σ ) can be obtained as As can be seen in Equation ( 5), the involvement of a spring does not depend on other springs, which means that the level of stiffness can be uniquely selected at any position without the need of returning to the equilibrium position.Moreover, the system is scalable.

System Identification of Discrete Variable Stiffness Actuator
This section focuses on identifying the unknown dynamical parameters of the discrete variable stiffness actuator model.If we assume that there is no gears backlash, the dynamic equations of the system can be written as where: θ m , and .. θ m are the angular position, the angular velocity, and the angular acceleration of the motor, respectively; θ L , and .. θ L are the angular position, the angular velocity, and the angular acceleration of the motor, respectively; • T in is the input torque of the motor; • T ext is the torque needed to lift the external load.
Note that some of the dynamic parameters shown in the model are fixed in the design phase.In particular, the stiffness of each spring (hence K eq ) and the gear ratio (N).However, our system may have eight possible values for K eq depending upon the status of clutches, as shown in Table 1.The goal of this section is to identify the remaining unknown model parameters-notably, the load inertia (I L ), the motor inertia (I m ), the load damping (B L ), and the motor damping (B m ).
The parameters identification is carried out in two steps.At both steps, the two-Degrees-of-Freedom dynamic model is converted into a single-DOF one by constraining the other.Firstly, we clamp θ L and set T in equal to zero to identify I m and B m .Thus, the differential equation for the system can be written as ..
The system is excited by any initial condition (IC), e.g., θ m (0) = 0.1 rad, and the response of θ m is recorded.The logarithmic decrement algorithm is applied to identify both parameters.For a more detailed description on the logarithmic decrement algorithm and how it works, the reader is encouraged to refer to [31][32][33][34].The system vibrates under the effect of the IC while the vibrations (θ m ) are recorded.As a result, the periodic time (τ d1 ), the damped natural frequency (ω d1 ), and the natural frequency (ω n1 ) are calculated from the response of θ m .On the other hand, the damping ratio (ζ 1 ) of the system is calculated by using the following equation: where: • θ 1 : the peak amplitude corresponding to time t 1 ; • θ m+1 : the peak amplitude corresponding to time t m+1 ; • m: number of complete cycles between θ 1 and θ m+1 ; • δ 1 : the logarithmic decrement for the motor side under the effect of initial condition (IC) θ m (0). Therefore, where ζ 1 is the damping ratio.As a result, I m and B m are calculated using the following equations: The whole process is repeated to identify I L and B L , but in this case θ m is clamped while keeping θ L free.The dynamic system as previously becomes a single DOF system, where it is excited by IC (θ L (0) = θ L0 = 0.1 rad).Hence, τ d2 , ω d2 , ω n2 , δ 2 , and ζ 2 are calculated from the response of θ L .Similarly, I L and B L are calculated using the following equations: where: • I L : the inertia load; • B L : the damping in the load side; • ζ 2 : the damping ratio for the load side under the effect of the IC θ L (0); • ω n2 : the natural frequency under the constraints of the second experiment.

Control System Design of Discrete Variable Stiffness Actuator
In this section, two control techniques are implemented on the dynamic model of the DVSA.The first control method is based on a type-one servo system design using linear quadratic regulator (LQR).This controller is classified as an optimal controller because it minimizes a quadratic cost function with respect to the dynamics of the system, and it has the capability to maximize energy saving as well as coping with the tracking of the desired response, as detailed in [12,13].On the other hand, the computed torque control (CTC) scheme is classified as a model-based control.In the following, we present both methods and numerical simulations to evaluate the performance of each for our DVSA.

Type-One Servo System Design Based on LQR
The system dynamic model is linear as shown in Equations ( 7) and (8).The idea is to construct a tracking system using LQR control to track the reference signals θ m and θ L .Thus, differential equations of our system can be written as .. . .
The corresponding state space representation for the binary controlled variable stiffness actuator (BcVSA) can be derived as follows: where q = q 1 q 2 q 3 q 4 T is called the state vector.Therefore, the state space representation is defined as follows: . where: • q = q 1 q 2 q 3 q 4 T is called the state vector and u = T in T ext T is called the input vector.
• C: the output matrix, which is an identity matrix.• D: the feedforward matrix, which is a zero matrix.
Consequently, a type-one servo system with an integral action is designed (see Figure 3).The corresponding state space representation for the binary controlled variable stiffness actuator (BcVSA) can be derived as follows: where  = [     ] is called the state vector.Therefore, the state space representation is defined as follows: where: •  = [     ] is called the state vector and  = [  ] is called the input vector.
• C: the output matrix, which is an identity matrix.
• D: the feedforward matrix, which is a zero matrix.
Consequently, a type-one servo system with an integral action is designed (see Figure 3).Thus, another state is added to the system to eliminate the steady state error.However, the following equations are used to construct the control system of the BcVSA.More details on this method can be found in [12]. where: is the dynamic error on the new state space vector.Thus, another state is added to the system to eliminate the steady state error.However, the following equations are used to construct the control system of the BcVSA.More details on this method can be found in [12]. where: • e(t) = q e (t) is the dynamic error on the new state space vector.
• u e (t) = u(t)-u(∞) is the dynamic error on the input vector.
where K is the gain matrix calculated in MATLAB by using the following command: where: • Q is a semi-positive definite weighting matrix or positive definite weighting matrix; • R is a positive definite weighting matrix.
The main features of this tracking system include its ability to track the reference commands and to regulate the other states.This controller, in reality, is classified as an optimal controller because it is minimizing a quadratic cost function to calculate the optimal response, which is achieved by solving the reduced form for the Riccati equation [12,13].

Tracking System Based on Computed Torque Control (CTC)
In general, this method is based on cancelling the nonlinear terms from the dynamic equations.However, the model of the BcVSA does not include any nonlinear term.We write the dynamic equation in the following form: where: • M = I m 0 0 I L is the mass matrix.; • B = B m 0 0 B L is the damping matrix; is the torque input vector; T is the angular position for the motor and for the load, respectively.
Based on that, the CTC provides the system with an input vector (T) to cancel the damping and the stiffness matrices and to convert the model to a unit mass model but with 'n' DOF (see Equation ( 26)).However, this is achieved by identifying the mass, the damping, and the stiffness matrices, respectively.
where M, B, and K are the estimated mass, the estimated damping, and the estimated stiffness matrices.In a case of have a perfect estimation process, it is expected that the controller will work perfectly.Otherwise, the performance of the controller will be affected based on the error in the estimation process.
with good estimates . where: where e is the error between the actual and the estimated values for the angular positions of the motor and the load, respectively.In general, it is recommended that high gains for k d and k p are used to reject the errors.

Numerical Simulations
In this section, we present the numerical simulations of our model based on the system identification and control methods detailed in previous sections.The proposed results are intended to be used on the hardware platform of our actuator in the near future.This section is divided into three major subsections.In particular, the first subsection focuses on verifying the proposed method for identifying the unknown dynamic parameters, while the second and the third subsections evaluate the performance of both control techniques, i.e., LQR and CTC, respectively.Both control techniques are tested under different possible conditions, i.e., while altering the status of each clutch in different combinations.

Evaluation of the System Identifications Algorithm
In order to verify the previously introduced parameters identification method, we built a Simulink scheme based on the dynamical equations of the BcVSA to identify the same parameters.The main idea is to assume that the parameter values are the same as previously identified and re-determine the parameter values under different conditions.Hence, the proposed algorithm will be considered correct if the estimated values confirm the assumed values or vice versa.Accordingly, in the first test, θ L is clamped while θ m is free, so the goal is to identify I m and B m .This test is repeated under various conditions, e.g., by choosing different values of K eq through the selection of different combinations of active and inactive clutches.Moreover, we also simulated this test under the effect of various values of ICs (less than 0.3 rad).The resulting parameters are detailed in Table 2.
In the second test, θ L is free while θ m is clamped to identify I L and B L .Again, the test is performed under the effect of various combinations of active and inactive clutches.The resultant relative errors for the estimated parameters are shown in Table 3.As a result, the logarithmic decrement algorithm method showed good performance in estimating the assumed parameters without any significant relative error, as detailed in Tables 2 and 3.

Evaluation of the Servo System Design Based on a Linear Quadratic Regulator (LQR)
A cubic trajectory planning algorithm was used to create the desired trajectories for θ m and θ L , respectively.The design parameters in this scheme are Q and R matrices.Nevertheless, suitable parameters were selected to perform the required task.Under these conditions, we tested our system by selecting different stiffness levels through the status (engaged/disengaged) of the clutches.5.2.1.All Clutches Are Active (K eq = 7K 0 ) the Entire Time In the first case, the status of all the clutches is active, i.e., K eq = 7K 0 .The resultant simulated response and the input signals are shown in Figures 4 and 5, respectively.As is clear from the plots, the LQR controller was able to track the desired reference trajectory.As a result, the logarithmic decrement algorithm method showed good performance in estimating the assumed parameters without any significant relative error, as detailed in Tables 2 and 3.

Evaluation of the Servo System Design Based on a Linear Quadratic Regulator (LQR)
A cubic trajectory planning algorithm was used to create the desired trajectories for  and  , respectively.The design parameters in this scheme are Q and R matrices.Nevertheless, suitable parameters were selected to perform the required task.Under these conditions, we tested our system by selecting different stiffness levels through the status (engaged/disengaged) of the clutches.

All Clutches Are Active (𝐾 = 7𝐾 ) the Entire Time
In the first case, the status of all the clutches is active, i.e.  = 7 .The resultant simulated response and the input signals are shown in Figures 4 and 5, respectively.As is clear from the plots, the LQR controller was able to track the desired reference trajectory.

Altering the Clutches' Status While the System is Running at Different Time Intervals
In this case, we used the same design parameters as chosen in the previous case, but this time we tested the controller ability to control the actuator by altering the stiffness level while the system

Altering the Clutches' Status While the System Is Running at Different Time Intervals
In this case, we used the same design parameters as chosen in the previous case, but this time we tested the controller ability to control the actuator by altering the stiffness level while the system is running.The cubic trajectory-planning algorithm was used for three seconds.Different stiffness levels (K eq ) were selected at different sub-time intervals as follows:
The simulation result for this case are shown in Figure 6.The corresponding input signal is shown in Figure 7.As a result, the controller successfully tracked the reference trajectory.The simulation result for this case are shown in Figure 6.The corresponding input signal is shown in Figure 7.As a result, the controller successfully tracked the reference trajectory.

Evaluating the Tracking System based on Computed Torque Control (CTC)
In order to evaluate the performance of the CTC, we repeated the simulations with CTC while using the same status of the clutches as previously chosen for LQR.

Evaluating the Tracking System Based on Computed Torque Control (CTC)
In order to evaluate the performance of the CTC, we repeated the simulations with CTC while using the same status of the clutches as previously chosen for LQR.5.3.1.All Clutches Are Active (K eq = 7K 0 ) the Entire Time In this test, all clutches are active i.e., K eq = 7K 0 .The simulation results are reported in Figure 8, while the corresponding input signal is shown in Figure 9.The CTC scheme also showed its ability to track the desired reference trajectory with little error.The major source of these errors are the saturation blocks on the control scheme, because the control signals must not exceed 35 Nm. when θ = 10 • (which is fixed according to actuator design criteria).However, the controller has the capability of tracking the desired reference signals without any errors in principle as long as the required input signals are lower than the saturated levels.In this test, all clutches are active i.e.  = 7 .The simulation results are reported in Figure 8, while the corresponding input signal is shown in Figure 9.The CTC scheme also showed its ability to track the desired reference trajectory with little error.The major source of these errors are the saturation blocks on the control scheme, because the control signals must not exceed 35 Nm. when  = 10° (which is fixed according to actuator design criteria).However, the controller has the capability of tracking the desired reference signals without any errors in principle as long as the required input signals are lower than the saturated levels.

Altering the Clutches Status at Different Time Intervals While the System Is Running
Similarly, in this case, the same procedure of applying different stiffness levels was achieved by altering the status of clutches at different sub-time intervals.Details are as follows:

•
In the first time interval [0 ≤ t < 1] second, the Bit 2 1 clutch is on.
The results of this case are presented in Figure 10 while the applied input signals are shown in Figure 11.The results of this case are presented in Figure 10 while the applied input signals are shown in Figure 11.In summary, both control schemes showed good performance in controlling our BcVSA model in a simulation environment.It was observed that the CTC method has comparatively fewer errors but is comparatively less efficient in terms of energy consumption.In addition, the CTC mainly depends on our knowledge about the model itself and on reducing any source of uncertainty on the dynamic parameters.However, the preceding problem can be partially solved by increasing the value of the gain matrices   and   .On the other hand, the LQR method showed a good capability to track the desired signals with relatively less energy consumption, in spite of the fact that the presented results were achieved by selecting higher values for the entries of the  matrix and smaller values for the entries of the  matrix.These choices mean that the priority is focused on decreasing the error between the desired and the simulated responses instead of focusing on minimizing the energy consumption.Overall, it was observed that the CTC method is more reliable for our system, as the error In summary, both control schemes showed good performance in controlling our BcVSA model in a simulation environment.It was observed that the CTC method has comparatively fewer errors but is less efficient in terms of energy consumption.In addition, the CTC mainly depends on our knowledge about the model itself and on reducing any source of uncertainty on the dynamic parameters.However, the preceding problem can be partially solved by increasing the value of the gain matrices k d and k p .On the other hand, the LQR method showed a good capability to track the desired signals with relatively less energy consumption, in spite of the fact that the presented results were achieved by selecting higher values for the entries of the Q matrix and smaller values for the entries of the R matrix.These choices mean that the priority is focused on decreasing the error between the desired and the simulated responses instead of focusing on minimizing the energy consumption.Overall, it was observed that the CTC method is more reliable for our system, as the error between the desired and the simulated responses are less than the LQR (see Figures 6 and 8).Table 4 shows the quantitative values of power and energy consumption for each case using both types of controllers.
We aim to study the robustness of both control approaches in the simulation and the prototype of our actuator in the future.

Conclusions
In this work, the design, modeling, system identification, and control of discrete variable stiffness actuators are presented.We used a logarithmic decrement algorithm for the system identification.In conclusion, the implemented approach for identifying the dynamic parameters is working perfectly as shown in the numerical results, but it is expected that, during the practical phase, the efficiencies of all the gears will cause an increase in errors, i.e., the backlashes may cause this problem.Nevertheless, these errors can be dramatically reduced by taking into account the efficiency of the gears into the dynamic equations.Additionally, this paper presented two control methods based on LQR and CTC along with their simulated results.Both methods were simulated under the effect of various combinations of active (engaged) and inactive (disengaged) clutches.Based on relatively fewer errors and a shorter settling time as compared to the LQR, results of the CTC are more reliable for our model.On the contrary, the LQR method was proven to be more energy-efficient as compared to the CTC, as shown in the simulated results.In future work, we plan to test the presented approaches on our hardware prototype and compare the results with the simulated ones.

Figure 1 .
Figure 1.The passive version of the discrete variable stiffness actuator (DVSA)

Figure 1 .
Figure 1.The passive version of the discrete variable stiffness actuator (DVSA).

Figure 2 .
Figure 2. (a) The CAD representation of discrete variable stiffness actuator.(b) The schematic representation of the working principle of the actuator.

Figure 2 .
Figure 2. (a) The CAD representation of discrete variable stiffness actuator.(b) The schematic representation of the working principle of the actuator.

Figure 3 .
Figure 3. Type-one servo system for controlling the binary controlled variable stiffness actuator (BcVSA).

Figure 3 .
Figure 3. Type-one servo system for controlling the binary controlled variable stiffness actuator (BcVSA).

Figure 4 .
Figure 4.The responses of the BcVSA based on the linear quadratic control (LQR) scheme under the condition when all the clutches are active ( = 7  ) during the entire simulation time.(Top) the input trajectory response; (Bottom) the output response.

Figure 4 .
Figure 4.The responses of the BcVSA based on the linear quadratic control (LQR) scheme under the condition when all the clutches are active (K eq = 7 K 0 ) during the entire simulation time.(Top) the input trajectory response; (Bottom) the output response.

Figure 4 .
Figure 4.The responses of the BcVSA based on the linear quadratic control (LQR) scheme under the condition when all the clutches are active ( = 7  ) during the entire simulation time.(Top) the input trajectory response; (Bottom) the output response.

Figure 5 .
Figure 5.The input signals of BcVSA based on the LQR scheme under the condition when all the clutches are active ( = 7  ) during the entire simulation time.

Figure 5 .
Figure 5.The input signals of BcVSA based on the LQR scheme under the condition when all the clutches are active (K eq = 7 K 0 ) during the entire simulation time.

Figure 6 .
Figure 6.The responses of DVSA based on the LQR scheme in the scenario of altering the status of clutches while the system is running.(Top) the input trajectory response; (Bottom) the output response.

Figure 6 .
Figure 6.The responses of DVSA based on the LQR scheme in the scenario of altering the status of clutches while the system is running.(Top) the input trajectory response; (Bottom) the output response.

Figure 6 .
Figure 6.The responses of DVSA based on the LQR scheme in the scenario of altering the status of clutches while the system is running.(Top) the input trajectory response; (Bottom) the output response.

Figure 7 .
Figure 7.The input signals for the LQR scheme while altering the clutches' status at different time intervals while the system is running.

Figure 7 .
Figure 7.The input signals for the LQR scheme while altering the clutches' status at different time intervals while the system is running.

Figure 8 .
Figure 8.The responses of BcVSA by using the computed torque control (CTC) scheme under the condition when all the clutches are active ( = 7 ) for the entire simulation time.(Top) the input trajectory response; (Bottom) the output response.

Figure 8 .
Figure 8.The responses of BcVSA by using the computed torque control (CTC) scheme under the condition when all the clutches are active (K eq = 7K 0 ) for the entire simulation time.(Top) the input trajectory response; (Bottom) the output response.

Figure 8 .
Figure 8.The responses of BcVSA by using the computed torque control (CTC) scheme under the condition when all the clutches are active ( = 7 ) for the entire simulation time.(Top) the input trajectory response; (Bottom) the output response.

Figure 9 .
Figure 9.The input signals for the CTC scheme when all the clutches are active ( = 7 ) for the entire simulation time.

Figure 9 .
Figure 9.The input signals for the CTC scheme when all the clutches are active (K eq = 7K 0 ) for the entire simulation time.

Figure 10 .
Figure 10.The responses of the BcVSA based on the CTC scheme in the scenario of altering the status of clutches while the system is running.(Top) the input trajectory response; (Bottom) the output response.

Figure 10 .
Figure 10.The responses of the BcVSA based on the CTC scheme in the scenario of altering the status of clutches while the system is running.(Top) the input trajectory response; (Bottom) the output response.

Figure 10 .
Figure 10.The responses of the BcVSA based on the CTC scheme in the scenario of altering the status of clutches while the system is running.(Top) the input trajectory response; (Bottom) the output response.

Figure 11 .
Figure 11.The input signals for the CTC scheme during altering status of the clutches at different time intervals while the system is running.

Figure 11 .
Figure 11.The input signals for the CTC scheme during altering status of the clutches at different time intervals while the system is running.

Table 1 .
The total possible stiffness levels of the DVSA.

Table 1 .
The total possible stiffness levels of the DVSA.

Table 2 .
Relative errors for the estimated parameters under the effect of active and inactive clutches in the case of the first test.

Table 3 .
Relative errors for the estimated parameters under the effect of active and inactive clutches in the case of the second test.

Table 4 .
The quantitative values of power and energy consumption for each case using both types of controllers.