Angelica Ginnante
email: aginnante@nimbl-bot.com

Enrico Simetti
email: enrico.simetti@unige.it

François Stéphane Caro
email: stephane.caro@ls2n.fr

François Leborne
email: fleborne@nimbl-bot.com

Task Priority Based Design Optimization of a Kinematic Redundant Robot

Keywords: Kinematic redundancy, Design optimization, Kinetostatic performance, Task priority kinematic control

come

Task priority based design optimization of a kinematic redundant robot

Introduction

Redundant manipulators are increasingly being used in many applications thanks to their ability to perform secondary tasks, which improves the robot performance [START_REF] Kumar | Design optimization of serial link redundant manipulator: an approach using global performance metric[END_REF]. This paper focuses on the design optimization of such robots. When building a new robot, an important step is optimizing its design with respect to relevant performance indices. Angeles, in [START_REF] Angeles | The design of isotropic manipulator architectures in the presence of redundancies[END_REF], proposed an approach to design isotropic redundant manipulators by minimizing the condition number of the robot kinematic Jacobian matrix. However, a complete way to optimize redundant robots is to consider global indices. In [START_REF] Kumar | Design optimization of serial link redundant manipulator: an approach using global performance metric[END_REF], the authors optimized a redundant serial manipulator using the global conditioning index. The design can also be optimized through both kinematic and dynamic global indices to obtain a robot with a high global dexterity through its workspace, ensuring high dynamic performance and energy efficiency [START_REF] Hwang | Design optimization method for 7 dof robot manipulator using performance indices[END_REF].

It is important to consider the main robot tasks during the robot design optimization process. An example limited to non-redundant manipulators is presented in [START_REF] Paredis | Kinematic design of serial link manipulators from task specifications[END_REF]. In [START_REF] Kim | A formulation for task based design of robot manipulators[END_REF], the authors worked with re-configurable modular manipulator systems to address the problem of task-based robot design. The work presented in [START_REF] Whitman | Task-specific manipulator design and trajectory synthesis[END_REF] explored a method to optimize the design parameters and the desired path together, considering the desired task with interesting results.

A critical issue in these task-oriented robot optimization processes is the complexity caused by the problem non-linearity. There are many ways to reduce the complexity, for example, breaking down the problem into multiple easier steps [START_REF] Kim | A formulation for task based design of robot manipulators[END_REF]. In [START_REF] Chikhaoui | Toward motion coordination control and design optimization for dual-arm concentric tube continuum robots[END_REF], the authors adopted a particle swarm optimization algorithm to determine the design parameter values of two collaborative robots. Other solutions employed a grid method [START_REF] Park | Task-oriented design of robot kinematics using the grid method[END_REF] or a complex direct search method [START_REF] Al-Dois | Task-based design optimization of serial robot manipulators[END_REF]. However, none of these techniques benefits from the robot redundancy to perform an optimization as a function of multiple tasks. In [START_REF] Maaroof | A robot arm design optimization method by using a kinematic redundancy resolution technique[END_REF], the authors developed a novel method for optimizing manipulator design using kinematic redundancy resolution. The authors use a combination of kinematic redundancy resolution and a multi-objective optimization algorithm. Kinematic redundancy resolution involves adding additional degrees of freedom in the place of the design parameters to optimize in the robot architecture. This allows performing the same task with different robot configurations. Through the Jacobian null-space projection, the optimization algorithm modifies the chosen design parameters together with the robot configuration solving the main task and some additional performance optimization sub-tasks. This algorithm gave promising results being able to identify some optimal parameter values. However, the algorithm was tested only on a two degrees of freedom non-redundant robot to which other two degrees of freedom were added for redundancy resolution optimization.

This paper stands from the idea of [START_REF] Maaroof | A robot arm design optimization method by using a kinematic redundancy resolution technique[END_REF] and extends the design optimization concept. As in [START_REF] Maaroof | A robot arm design optimization method by using a kinematic redundancy resolution technique[END_REF], virtual prismatic or revolute joints replace the design parameters to be optimized by a kinematic redundancy resolution algorithm. The first main novelty consists in the chosen kinematic control framework that can solve both equality and inequality tasks simultaneously. For example, an inequality task employed in the optimization process could be the joint limit avoidance. The joint limit task can be applied on the both real and virtual joints resulting fundamental for respecting the robot and design parameter constraints. Another important feature of the kinematic control framework is the task activation/deactivation ability to prevent over-constraining the system without causing discontinuities in the velocity generation. The second significant difference between the two design algorithms lies in the optimization methodology implementation. In [START_REF] Maaroof | A robot arm design optimization method by using a kinematic redundancy resolution technique[END_REF], the arm design is optimized for singular poses in its workspace. Here, the proposed optimization algorithm is based on a trajectory planning task for machining application. The robot end-effector moves along some trajectories that characterize a workspace area. The results are guidelines for building optimal robots for that desired workspace area. The design optimization algorithm here presented is then applied on a 21 degrees of freedom manipulator tracing a set of squared paths. The performance of the robot are improved using three kinetostatic indices.

This paper is outlined as follows. Section 2 describes the employed task-priority-based kinematic control algorithm and illustrates the tasks used to optimize the manipulator kinetostatic performance. Section 3 explains the new method for design optimization of robotic manipulators, which is the central topic. Section 4 discusses and compares the obtained results. Section 5 presents the conclusions and future work.

Background and methods

In this section, the kinematic control algorithm used to perform the task-oriented design optimization and the tasks employed to improve the robot kinetostatic performance are described.

Kinematic control algorithm

In the literature, many kinematic control algorithms exist to deal with redundant robots and exploit the redundancy to solve simultaneous tasks. In [START_REF] Escande | Hierarchical quadratic programming: Fast online humanoidrobot motion generation[END_REF], the authors present a hierarchical quadratic programming control algorithm used to find a solution to multiple and antagonistic objectives for humanoid robot motion generation. Another interesting multiple tasks control framework is presented in [START_REF] Di Lillo | Handling robot constraints within a set-based multi-task priority inverse kinematics framework[END_REF], called Set-Based Multi-Task Priority Inverse Kinematics Framework. In [START_REF] Flacco | Prioritized multi-task motion control of redundant robots under hard joint constraints[END_REF], the so-called Saturation in the Null Space (SNS) algorithm implements a predictive prioritizing technique for multiple tasks. However, these task priority algorithms generate discontinuities in the null space projector when activating or deactivating one or more tasks. Here, the new optimization method employs a taskpriority based control algorithm named Task Priority Inverse Kinematic (TPIK), presented in [START_REF] Simetti | A novel practical technique to integrate inequality control objectives and task transitions in priority based control[END_REF][START_REF] Simetti | Task priority control of underwater intervention systems: Theory and applications[END_REF][START_REF] Simetti | A task priority approach to cooperative mobile manipulation: Theory and experiments[END_REF][START_REF] Ginnante | Kinetostatic optimization for kinematic redundancy planning of nimbl'bot robot[END_REF]. The distinctive aspect of this algorithm is its ability to activate and deactivate control tasks without causing discontinuities in the control variables. In [START_REF] Ginnante | Kinetostatic optimization for kinematic redundancy planning of nimbl'bot robot[END_REF], the TPIK algorithm was employed in the control of a highly redundant robot demonstrating its ability of optimizing the robot configuration. Some definitions are taken from [START_REF] Simetti | Task priority control of underwater intervention systems: Theory and applications[END_REF] that are necessary for the introduction of the kinematic control algorithm. The vector q ∈ R n is the joint variable vector, describing the arm configuration, where n is the number of joints. The joint velocities are collected in the vector q ∈ R n .

The control objectives are used to define the robot goals. A control objective corresponds to a scalar variable x(q) computed as a function of the robot configuration vector q and represents the state of one task. There exists two types of control objectives, equality and inequality. Equality control objectives aim to satisfy the relationship x(q) = x 0 . Inequality control objectives take the form x(q) ≤ x M , or x(q) ≥ x m , or both simultaneously, where x m and x M are the lower and upper bounds of the variable x(q). A conceptual division can be done to divide the control objectives into categories depending on their scope, simplifying the identification of their priority levels. These categories are system safety objectives, e.g. joint limits or obstacle avoidance, action oriented objectives, e.g. reaching a pose or following a trajectory, and optimization objectives, e.g. minimizing the joint velocities or optimizing the kinetostatic performance indices. Then, each scalar control objective is associated with a feedback reference rate ẋ. A closed-loop rate control law drives the actual variable x(q) to the desired point x * with the associated feed-forward changing rate ẋ * and is defined as

ẋ = λ(x * -x(q)) + ẋ * , (1)
where λ is a positive gain related to the target convergence rate. The actual derivative of x is defined as a function of the joint velocity vector q as follows:

ẋ(q, q) = J task (q) q = ∂x ∂q 1 . . . ∂x ∂qn q,

where q = [q 1 . . . q n ].

An activation function a i (x) ∈ [0, 1] is associated with each control objective x(q) indicating whether or not the objective has to be taken into account at a specific time instant. The tasks associated with inequality control objectives are relevant only when the scalar variable x(q) is close or outside its validity region bounded by x m and x M . Thus, the activation function is equal to zero within the validity region of its associated inequality objective. When it is near or out of the validity region, its value is set to one with a smooth transition. For tasks associated with equality control objectives, the activation function is always equal to a i (x) 1 because these control objectives never become inactive.

A specific priority level is assigned to each task based on its objective importance. The highest priority tasks are solved as first using the available robot degrees of freedom and are not affected by the lower priority ones. Hence, lower priority tasks are handled with the remaining robot degrees of freedom. When two or more tasks have the same priority level, they are grouped into a multidimensional control task. A specific list of prioritized tasks is called control action A . A control action A collects the vectors and matrices associated with each priority level These vectors and matrices are defined in [START_REF] Ginnante | Kinetostatic optimization for kinematic redundancy planning of nimbl'bot robot[END_REF] as: ẋk = [ ẋ1,k , ẋ2,k , . . . , ẋm k ,k ] is the vector containing all the reference rates of the scalar control tasks, where m k is the number of scalar tasks for the priority level k.

J k is the Jacobian matrix associated with the k th task vector [ ẋ1,k , . . . , ẋm k ,k ] with respect to the joint velocity vector q.

A k = diag(a 1,k , . . . , a m k ,k ) is a diagonal matrix of the activation functions.

The TPIK algorithm solves a sequence of nested minimization problems to compute the system velocity reference vector q that meets all the priority requirements inside a given control action A . The solution S k to the k th priority level can be written as

S k = arg R -min q∈S k-1 ||A k ( ẋk -J k q)|| 2 , (3)
where S k-1 is the manifold of all the previous priority level solutions. The notation R -min highlights that each minimization is performed through specific regularized space projections to implement priorities among the tasks defined in [START_REF] Simetti | A novel practical technique to integrate inequality control objectives and task transitions in priority based control[END_REF]. In addition to Eq. ( 3), other regularization costs are included. These regularization costs avoid discontinuities in the system velocity vector due to kinematic and algorithmic singularities. In [START_REF] Simetti | A novel practical technique to integrate inequality control objectives and task transitions in priority based control[END_REF], the authors fully describe these regularization costs that are not analyzed in this paper. A significant advantage of the TPIK algorithm is the use of the activation functions to handle inequality control objectives without over-constraining the system. Both equality and inequality objectives require a certain amount of robot degree-of-freedom to be handled. When an inequality task is inside its validity region, the activation function is set to zero to avoid locking any degrees of freedom. So, safety tasks, like joint limits, can be placed at the top of the hierarchy without over-constraining the optimization problem.

Finally, the TPIK algorithm adopts another continuous sigmoidal function a p (p) that includes the previous and current executed actions and the time elapsed in the current step to perform a smooth transition between two actions. This function a p (p) is employed in the algorithm execution together with a i (x). More details are presented in [START_REF] Simetti | A novel practical technique to integrate inequality control objectives and task transitions in priority based control[END_REF].

Tasks for optimization of robot kinematic performance

The TPIK algorithm employs different kinetostatic performance tasks to optimize the robot design, namely dexterity, manipulability and robot transmission ratio, to optimize the robot design. Some of these tasks were used with the TPIK algorithm in [START_REF] Ginnante | Kinetostatic optimization for kinematic redundancy planning of nimbl'bot robot[END_REF]. To introduce the definition of these indices, let us consider the kinematic Jacobian matrix J e related to the end-effector velocity:

t = ṗ ω = J e (q) q = J l (q) J a (q) q, (4)
where t = ṗ , ω ∈ R 6 is the robot end-effector twist, with ṗ ∈ R 3 and ω ∈ R 3 the linear and angular velocity vectors of the end-effector, respectively. Since the kinematic Jacobian matrix J e contains non-homogeneous terms, namely linear and angular, it needs to be weighted to compute the kinetostatic performance indices correctly. The weighting of J e employs the characteristic length L that was introduced in [START_REF] Angeles | The design of isotropic manipulator architectures in the presence of redundancies[END_REF] to solve the absence of dimensional homogeneity in the kinematic Jacobian matrix entries and is computed in [START_REF] Khan | The kinetostatic optimization of robotic manipulators: The inverse and the direct problems[END_REF]. To weight J e , the revolute joint columns of the linear kinematic Jacobian matrix part are divided by L. The weighted kinematic Jacobian matrix is written as J w . The weighting is a critical issue when analyzing the kinetostatic performance of J e [START_REF] Khan | Jacobian matrix normalization -A comparison of different approaches in the context of multi-objective optimization of 6-dof haptic devices[END_REF][START_REF] Zargarbashi | Posture optimization in robot-assisted machining operations[END_REF].

Manipulability

The manipulability is an index that measures the kinematic abilities of the robotic system through its Jacobian matrix J w [START_REF] Yoshikawa | Manipulability of robotic mechanisms[END_REF]. The manipulability of a manipulator is defined as

µ = det(J w J w ), (5)
and amounts to the product of all the singular values of J w . The higher the manipulability value, the larger the manipulability hyper-ellipsoid and the better the kinematic performance of the mechanism [START_REF] Angeles | Fundamentals of robotic mechanical systems: theory, methods, and algorithms[END_REF]. It should be noted that the manipulator reaches a kinematic singularity when µ vanishes.

The derivative of the manipulability as a function of the joint variables is explained in [START_REF] Park | Analysis and control of kinematically redundant manipulators: An approach based on kinematically decoupled joint space decomposition[END_REF] and used in [START_REF] Marani | A real-time approach for singularity avoidance in resolved motion rate control of robotic manipulators[END_REF]:

∂µ ∂q i = µ trace ∂J w ∂q i J + w , (6)
where the matrix J + w is the pseudo-inverse of the weighted kinematic Jacobian matrix. Hence, the manipulability Jacobian matrix J µ as a function of the joint variables is:

J µ = ∂µ ∂q 1 . . . ∂µ ∂qn , (7)
where n, which represents the number of columns of J w , is the dimension of joint space.

Dexterity

The dexterity η(J w ) characterizes the kinematic performance of a manipulator in a given configuration and is defined as the inverse of the conditioning number κ(J w ) of its Jacobian matrix [START_REF] Angeles | Kinematic isotropy and the conditioning index of serial robotic manipulators[END_REF]:

κ(J w ) = ||J w || ||J -1 w || and η(J w ) = 1/κ(J w ). ( 8
)
The index η is bounded by 0 and 1. The higher η, the better the manipulator dexterity.

The manipulator reaches an isotropic posture when η = 1. The smaller η, the worse the manipulator dexterity and the closer to a singularity. Moreover, η can be defined as the ratio between the smallest and largest singular values of J w indicating how close the manipulability hyper-ellipsoid is to being a hyper-sphere [START_REF] Pond | Formulating jacobian matrices for the dexterity analysis of parallel manipulators[END_REF].

In order to obtain an analytical expression of η, the Frobenius norm of J w can be used [START_REF] Rakotomanga | Kinetostatic performance of a planar parallel mechanism with variable actuation[END_REF]:

η(J w ) = m trace(J w J w ) trace[(J w J w ) -1 ] , (9)
where m, which represents the number of rows of J w , is the dimension of the task space.

The following definitions are introduced to increase the readability of the equations:

γ 1 trace(J w J w ) and γ 2 trace[(J w J w ) -1 ], (10)
where is the define operator. With these definition, it follows that

η(J w ) = m γ 1 (J w ) γ 2 (J w ) . (11)
Then, the dexterity Jacobian matrix is determined to relate the velocity rate of η with respect to the joint velocity vector q. The Frobenius formula used in Eq. ( 9) expresses η as a function of joint position vector q in an analytical way allowing its derivation. So, the derivative of Eq. ( 9) with respect to each joint position q i ∈ q is

∂η ∂q i = -η ∂γ 1 ∂q i 1 γ 1 + 1 γ 2 ∂γ 2 ∂q i , (12)
where

∂γ 1 ∂q i = 1 γ 1 trace J w ∂J w ∂q i and ∂γ 2 ∂q i = 1 γ 2 trace -J w ∂J w ∂q i (J w J w ) 2 . ( 13
)
In conclusion, the dexterity Jacobian matrix J η as a function of the joint variables is:

J η = ∂η ∂q 1 . . . ∂η ∂qn , (14)
where n, which represents the number of columns of J w , is the dimension of joint space.

Robot transmission ratio

The robot transmission ratio (RTR) ρ(J w ) quantifies the effectiveness of the actuator force in producing a prescribed robot motion [START_REF] Zargarbashi | Posture optimization in robot-assisted machining operations[END_REF]. It corresponds to the angle between the joint velocity q and torque τ vectors in the joint space and is defined as

ρ = |τ q| ||τ || || q|| = | cos ∠(τ , q)|. ( 15
)
This metric is bounded between 0 and 1. In case of kinetostatic redundancy, ρ can be expressed in terms of the end-effector twist t and the wrench w applied to it, leading to

ρ = |w t| ||J w w|| ||J + w t|| , ( 16
)
where t is the robot end-effector twist defined in Eq. ( 4) and w = [f , m ] is the wrench that collects the forces f and moments m exerted by the environment on the end-effector.

To ensure that ρ is dimensionless, the linear part ṗ in t and the moment m in w are divided by the characteristic length L.

The RTR Jacobian matrix is obtained upon differentiation of Eq. ( 16) with respect to each joint position q i ∈ q:

∂ρ ∂q i = ρ w J w ∂J w ∂q i w||J + w t|| 2 -||J w w|| 2 t J + w ∂J + w ∂q i t (||J w w|| ||J + w t||) 2 , ( 17
)
where the values of the end-effector twist t and wrench w are known from the trajectory planning. The derivative of the pseudo-inverse weighted kinematic Jacobian matrix ∂J + w /∂q i is defined in [START_REF] Golub | The differentiation of pseudo-inverses and nonlinear least squares problems whose variables separate[END_REF].

The RTR Jacobian matrix J ρ as a function of the joint variables is:

J ρ = ∂ρ ∂q 1 . . . ∂ρ ∂qn , (18)
where n = columns(J w ) is the dimension of joint space.

Design and optimization method

This paper proposes a new optimization method intended for redundant robots. The main idea is to optimize the design of a robot exploiting the same algorithm used during its real-time control. In particular, the idea is to benefit from the fact that the TPIK algorithm was specifically designed to control highly redundant robot manipulators. This section presents the problem formulation and all the phases that compose this optimization method. In the considered case study, the robot main application requires following a trajectory inside a defined workspace while maintaining satisfactory kinetostatic performance. However, this optimization methodology can be applied to any robot and task.

Problem formulation

The application considered is the tracking of a set of trajectories describing a workspace area while maintaining high kinetostatic performance. So, two principal inputs are defined and given to the optimization algorithm. The first one is a series of p trajectories with desired orientations and velocities. These trajectories describe the workspace area where the robot will work. The second input is the list of tasks that the robot should achieve while tracking the trajectories. The main tasks are related to the end-effector pose and velocity to track the desired trajectories. A safety task is added for joint limit compliance. Three tasks based on dexterity, manipulability and RTR evolution are used to improve the robot kinetostatic performance. While the dexterity and the manipulability are kinematic performance indices, maximizing them forces the manipulability hyper-ellipsoid to be as big as possible and close to a hyper-sphere. So, ideally, the robot will be able to move with the same high velocity amplification factor in all directions while reducing actuator velocity limits. For this reason, the optimization method uses both dexterity and manipulability in design optimization. A final task based on the robot center of mass position is included. This task aims to to maintain the robot center of mass position as close as possible to the base and, as a result, to maximize the robot compactness.

Table 1 reports all the information about the tasks used during the optimization. The dexterity, manipulability and RTR tasks have the same priority level since they have the same relevance in the design optimization process. The optimization algorithm uses a convex combination (η, ν, ρ) of the three kinetostatic performance indices to rate and compare the kinetostatic performance of the obtained designs. It should be noted that η and ρ are bounded between [0, 1] whereas µ is not bounded, [0, ∞). So, it is necessary to bound µ in the range [0, 1] before writing the convex combination. A new index called bounded manipulability ν is defined as:

ν = 1 - 1 1 + µ . ( 19
)
When µ = 0 then ν = 0, and when µ = ∞ then ν = 1. Now, η, ν and ρ are all bounded between 0 and 1 and can then be used in a convex combination: where λ 1 , λ 2 and λ 3 are scaling factors. Since all the kinetostatic performance indices are valid in the same range, the weighting factors are selected as λ 1 = λ 2 = λ 3 = 1/3.

(η, ν, ρ) = λ 1 η + λ 2 ν + λ 3 ρ, (20)

Preliminary phase

The optimization algorithm is divided into two main phases. The first is the candidate generation phase, where several design candidates are generated. The second phase is the candidate selection one, where the robot designs obtained from the optimization process are evaluated and compared. Before these phases, some preliminary steps are necessary. The robot design parameters ζ to be optimized need to be selected, for example, the link length or the angular offset amplitude between two consecutive joints. Then, controllable virtual joints are inserted in the robot architecture to substitute these design parameters. A link is substituted by a prismatic joint to control its length, while an angular offset by a revolute joint to modify its angle. So, the vector of the design parameters ζ is converted in the virtual joint vector q v . During the candidate generation phase, the vector of the joint variables contains both the position of the real and virtual joints q = [q r , q v ] ∈ R n , with q r ∈ R r , q v ∈ R v and the robot degrees of freedom is n = r + v. The joint limit task constrains the virtual joints q v value between the desired limits and so the resulting value of the design parameters ζ is limited.

Candidate generation phase

This phase employs the robot with real and virtual joints q = [q r , q v ]. The robot configuration vector q is randomly initialized. From here, the robot is moved to track all the p trajectories in a random order to ensure more general results. The robot reaches the starting pose of the trajectories and tracks it entirely under the kinematic control of the TPIK algorithm. When the robot finishes tracking a trajectory, it is moved to another one until it follows all the p trajectories. At equidistant time steps on each trajectory, the optimization algorithm saves the values of the virtual joints q v , which will be used as design parameters ζ, and the value in that configuration. Once the robot has tracked Random initialization of q = [q r , q v ].

3:

Random shuffle of p.

4:

for k := 1 → p do 5:

Reach starting pose of p(k) using A 1 .

6:

while Trajectory p(k) not finished do 7:

Move to next step on p(k) using A 2 .

8:

if t i = equidistant time step then 9:

Save q v as ζ t i . Compute weighted average ζ for all [t 0 , . . . , t f ]. 15: end for all the p trajectories, the optimization algorithm makes the weighted average ζ of ζ using the values already stored at equidistant steps along each trajectory and the associated as weighting factor:

ζ = t 0 ζ t 0 + • • • + t f ζ t f t 0 + • • • + t f , (21)
where [t 0 , . . . , t f ] are the equidistant time steps in which the design parameter vector ζ was saved. The average result ζ is also stored. The process described above is repeated several times. At the end of this phase, the optimization algorithm has collected a list of candidate designs described by ζ. Algorithm 1 sums up all the steps of the candidate generation phase. During this phase, the optimization algorithm computes the kinetostatic performance indices η, ν and ρ of the real robot architecture without virtual joints. To correctly compute these metrics and their derivatives, the columns corresponding to the virtual joints q v in J w and its derivative ∂J w /∂q i , ∀q i ∈ q are set to zero.

Candidate selection phase

The candidate selection phase tests the designs obtained from the candidate generation phase to identify the best one. Here, the virtual joints are removed from the robot architecture vector and replaced by constant links and offsets whose value was stored inside ζ. Hence, from now on q = [q r ] and v = 0, and the total amount of degrees of freedom n = r. The robot configuration vector q is randomly initialized for each optimized design. Again, the robot is moved to the starting pose of one trajectory and tracks it entirely under the kinematic control of the TPIK algorithm. When the robot finishes tracking one trajectory, it is moved to the next one until it has followed all the p trajectories. At each step on each trajectory, the optimization algorithm stores the value of for that robot configuration. for r := 1 → number of repetitions do 3:

Random initialization of q = [q r ].

4:

for k := 1 → p do 5:

Reach starting pose of p(k) using A 1 .

6:

while Trajectory p(k) not finished do 7:

Move to next step on p(k) using A 2 .

8:

Save t i at the time step t i . Compute ˆ = 0.5 min + 0.25 ( + max ).

12:

end for 13: end for 14: Select the best design as max(ˆ ).

When the robot has tracked all the trajectories, the optimization algorithm computes ˆ as ˆ 0.5 min + 0.25 ( + max ), [START_REF] Angeles | Fundamentals of robotic mechanical systems: theory, methods, and algorithms[END_REF] where min , and max are respectively the minimum, mean and maximum of values along all the p trajectories. The ˆ value is used to compare the designs and identify the best one. The process is repeated for all the designs obtained from the previous phase. The design that has the highest ˆ is identified as the best one. Algorithm 2 sums up all the steps of the candidate selection phase.

Application and results

The following section introduces the robot used to test the proposed optimization method and discusses the obtained results.

Robot under study

A new highly redundant robot actuated by a mechanism recently developed and patented by the company Nimbl'Bot [START_REF] Dufau | Articulated robot arm[END_REF] is used as a case study in this section. Both the mechanism and the robot are fully introduced in [START_REF] Ginnante | Kinetostatic optimization for kinematic redundancy planning of nimbl'bot robot[END_REF][START_REF] Ginnante | Design and kinematic analysis of a novel 2-dof closed-loop mechanism for the actuation of machining robots[END_REF]. The mechanism that actuates the robot has two degrees of freedom, and ten mechanisms are arranged in a serial way assembling the Nimbl'Bot robot, shown in Fig. 1. The robot is divided into three regions attached by two links: a shoulder composed of three mechanisms, an elbow consisting of four mechanisms, and a wrist composed of three mechanisms. At the end of the robot, a revolute joint is attached to adjust the wrist orientation. In total, the robot has 21 degrees of freedom, i.e., r = 21. Thus, this robot is hyper-redundant, which means that the end-effector can reach a given pose in an infinite number of configurations [START_REF] Chirikjian | Kinematically optimal hyper-redundant manipulator configurations[END_REF]. In this case, the dimensions of the mechanism are constant and not included in the optimization. The optimized robot design parameters are the link lengths l 1 and l 2 and the amplitude of three angular offsets β 1 , β 2 and β 3 , shown in Fig. 2. The angular offsets are respectively inserted between the first link and first mechanism of the elbow, the second link and first mechanism of the wrist, and between the second-to-last and last mechanisms of the wrist. The angles β 1 , β 2 and β 3 rotate about the axis x, depicted in red in Fig. 2. An ending tool is mounted on the last revolute joint of the robot. The design parameters of the tool are inserted in the optimized variables, namely its length l t and orientation β t , shown in Fig. 3. So, the number of virtual joints is v = 7 and the robot with virtual joints has 28 degrees of freedom, i.e. r + v = 28. Table 2 gives the main parameter values and limits of the robot.

The design parameters under optimization are limited, but the limits are large enough not to over-constrain the optimization process. The link lengths l 1 and l 2 are constrained in the range [0, 3] m, and the angular offsets amplitude β 1 in [-3π/4, π/4] and β 2 and β 3 in [-3π/4, 3π/4]. The tool length l t is valid in the range [0, 0.5] m, and the tool orientation offsets β t in [0, 3π/4]. These ranges are used by the joint limits task to constrain the virtual joints. The Nimbl'Bot mechanism actuators can rotate infinitely and have no limits.

Workspace area

The robot application requires tracking horizontal and vertical paths in a cube whose sides are 0.7 m × 0.7 m long and centered in (x, y, z) = (0.0, 0.75, 0.65) m. The orientation of the end-effector is expressed in terms of roll φ, pitch θ and yaw ψ. The angle Figure 4 shows the p = 4 paths during the candidate generation and selection phases of the optimization process. Two paths are horizontal, i.e. the green and the magenta ones, and the other two are vertical, i.e. the blue and the black ones. The horizontal paths have the corners placed in x = (-0.35, 0.35) m and y = (0.4, 1.1) m at the height of z = 0.3 m, green path, and z = 1.0 m, magenta path. The vertical paths have the corners placed in x = (-0.35, 0.35) m and z = (0.3, 1.0) m at the depth of y = 0.4 m, blue path, and y = -1.1 m, black path. The small arrows along the paths express the z-axis orientation of the end-effector frame while following the paths. These four paths were chosen for the design optimization because they describe the cube where the robot is supposed to work in the real world. Then, the kinetostatic performance of the best and worst designs are compared on two new paths, one horizontal and one vertical, placed inside the cubic workspace area to confirm the optimization process ability to identify the best design for kinetostatic performance. Figure 5 shows the two paths used for testing the design obtained from the optimization process. The horizontal path has the corners placed in x = (-0.35, 0.35) m and y = (0.4, 1.1) m at the height of z = 0.65 m, red path. The vertical path has the corners placed in x = (-0.35, 0.35) m and z = (0.3, 1.0) m at the depth of y = 0.75 m, cyan path. Table 3 gives the trajectory details.

Computational time evaluation and comparison

Before optimizing all the design variables of Figs. 2 and3, a brute force optimization process is run to have a baseline for the comparison with the proposed methodology. In this case, only the design variables l 2 ∈ [0, 1] m and β 2 ∈ [-3π/4, 3π/4] are optimized as a function of the green and blue paths shown in Fig. 4. The other design variables are set to l 1 = 0.2 m, β 2 = β 3 = 0, l t = 0.1 m and β t = π/4. Several robot designs are generated combining the discretization of l 2 and β 2 over their allowed ranges. Namely, l 2 is bounded between 0 and 1 with a step equal to 0.1 m, leading to 11 possible values. The parameter β 2 is bounded between -3π/4 and 3π/4 with a step equal to 5 • , leading to 55 values. Combining l 2 and β 2 generates 605 designs that are tested on the green and blue trajectories. The goal is to identify the design with the highest performance in terms of ˆ . Each design is tested ten times on each trajectory to obtain general results. If a design is not able to follow one or more trajectories, it is discarded.

Then, the design optimization process presented in this paper is run optimizing l 2 and β 2 as a function of the green and blue trajectories and the kinetostatic performance. The candidate generation phase is repeated four times and produces only 76 designs. So, the Each design is tested ten times along each path to obtain general results. If one design is not able to follow one or more trajectories, it is discarded. Table 4 shows the designs obtained from the two tests described above and their kinetostatic performance ˆ , η, ν and ρ. The values of η, ν and ρ are computed in the same way of Eq. ( 22) using the η, ν and ρ values collected at path step. The best designs identified by the two methods are highly similar. However, the optimization process reached a higher ˆ because the identified values for l 2 and β 2 were not included in the discretized values. Moreover, the time required by the optimization process, 20 minutes, is lower than testing all the designs obtained from the discretization, more than 120 minutes. Therefore, the optimization method is six times faster in this case. In this case, only two design variables were considered. Testing all the possible combinations for seven design parameters, which correspond to more than two billion designs, would be highly time consuming.

Optimization test and results

Here, the optimization process is applied to all the design variables of Figs. 2 and3. The optimization is performed with respect to the four desired trajectories of Fig. 4, and the kinetostatic performance indices, dexterity η, bounded manipulability ν and RTR ρ. The candidate generation phase is run ten times generating 370 designs. Then, the candidate selection phase tests these designs on the four trajectories collecting , η, ν and ρ at each step. In total, 282 designs passed the candidate selection phase. The others could not properly track all the trajectories and were discarded. During the candidate selection phase, each design follows each trajectory ten times starting from a random configuration to obtain general results. Table 5 shows the machine and implementation details.

Table 6 reports the correlation coefficients between the kinetostatic performance indices ˆ , η, ν and ρ and the optimized design parameters for all the solutions obtained from the optimization process. The higher the absolute value of the correlation coefficient, the darker the shade of blue associated. A positive correlation coefficient indicates a directly proportional relation between two parameters. A negative correlation coefficient means a indirectly proportional relationship. The value of η, ν and ρ is computed in the same way of Eq. ( 22) using the η, ν and ρ values collected at each step on each trajectory. The diagonal of Table 6 collects the standard deviation of the kinetostatic performance indices and design parameters between all the solutions identified by the optimization process. Figures 6 and7 show the box plots of all the values assumed by the kinetostatic indices ˆ , η, ν and ρ and the design parameters for the selected candidates, respectively. From Table 6, it is clear that ν is directly linked to the size of l 1 and l 2 . The longer the links, the higher ν. Moreover, a longer link l 1 tends to be associated with a longer l 2 . Another interesting aspect concerns η. The higher η, the higher both ν and ρ. This means that the same design parameters influence An analysis of the designs that performed the best and the worst in terms of ˆ , η, ν and ρ is accomplished. The right side of Table 7 shows the parameter values of the designs that performed the best and worst for each index, i.e. ˆ , η, ν and ρ, along the four trajectories used by the optimization process, shown in Fig. 4. In some cases, one design had the best or worst performance for more than one index. Figures 8 and9 present the link lengths and offset amplitudes collected in Table 7 on graphs. The middle of Table 7 reports the values of ˆ , η, ν and ρ obtained during the candidate selection phase. The green and red colors highlight the highest and lowest values for ˆ , η, ν and ρ during the candidate selection phase. One design has the worst performance in terms of ˆ , η and especially for ν. Small links characterize this design. In fact, the design that is characterized by the best ν has the largest links. The value of ρ is generally higher for all the designs, best and worst. This means that the Nimbl'Bot robot generally has high performance in terms of ρ, independently of the design parameters. Finally, β 3 is generally close to 0 • while β t is generally high and take the place of β 3 . So, β 3 could be removed from the optimized design parameters. Figures 10 and11 show the designs that reached the best and worst ˆ . 7.

Selected designs performance inside desired workspace

The designs collected in Table 7 are tested on two trajectories inside the desired workspace area, shown in Fig. 5, which are different from the ones used during the candidate generation and selection phases. This testing phase is performed to show that the best and worst designs maintain similar performance in tracking trajectories inside the volume defined by the four ones used in the optimization process. The left side of Table 7 reports the values of ˆ , η, ν and ρ obtained during the testing phase. The green and red colors highlight the highest and lowest values for ˆ , η, ν and ρ during the testing phase. The results show the same trend during the candidate selection and testing phases meaning that each design has uniform performance in the desired workspace area. To conclude, this optimization process identified several possible optimal designs for a specific set of tasks in the desired workspace area. There is no unique global solution, but more than one design can reach high values for ˆ . For designs that share the same global index ˆ , choosing the final one is a trade-off between optimizing one of the three indices that contribute to ˆ .

Conclusions

This paper presented and analyzed a new method for designing and optimizing kinematic redundant manipulators as a function of their application and kinetostatic performance. The process employs a task priority kinematic control algorithm to solve several tasks simultaneously exploiting the robot high kinematic redundancy. First, the process was tested on a simple case with two design variables and two paths. The obtained best design was compared with the results obtained from testing some discretized combinations of the two design variables. The best designs in terms of the global kinetostatic performance index ˆ obtained from the two tests are close and the optimization process converged toward a better design. This means that the optimization process can identify the best solution for a specific problem without testing all possible combinations of design parameter values and reduce the computational time.

Then, this new method worked on more design parameters identifying several possible Nimbl'Bot manipulators with high kinetostatic performance for the desired workspace area. It also pointed out some general guidelines for building a kinetostatic efficient Nimbl'Bot robot. However, there is no specific design trend to obtain high kinetostatic performance. This is due to the high redundancy of the robot that admits infinite possible robot configurations to solve the same problem. Moreover, the kinetostatic performance index ˆ is a linear combination of dexterity η, bounded manipulability ν and RTR ρ. Choosing the best design is a trade-off between these metrics.

Future work will involve more analysis to identify if some global guidelines exist for building a Nimbl'Bot robot design with high kinetostatic performance, giving more importance to one index than others or searching for specific relationships between the design parameters. Other design parameters will be included to study how they affect the kinetostatic performance. New tasks based on self-collision and obstacle avoidance will be added to consider more complex problems.

Algorithm 1

 1 Candidate Generation PhaseRequire: Actions A 1 = Reach Pose and A 2 = Follow Path and trajectory vector p.1: for r := 1 → number of repetitions do 2:

Algorithm 2

 2 Candidate Selection Phase Require: Actions A 1 = Reach Pose and A 2 = Follow Path and trajectory vector p. 1: for d := 1 → number of optimized designs do 2:

Figure 1 :

 1 Figure 1: The 21 degrees of freedom robot under optimization; the shoulder and the wrist are composed of three mechanisms, and the elbow is made of four mechanisms. In the bottom left corner, the Nimbl'Bot mechanism is shown with its design parameters half height h and slope α.

Figure 2 :

 2 Figure 2: Robot design parameters under optimization: link lengths l 1 and l 2 , angular offsets β 1 , β 2 and β 3 .

Figure 3 :

 3 Figure 3: Ending tool design parameters under optimization: tool length l t and orientation β t .Table 2: Robot details. Mechanism half height h 0.07 m Mechanism slope α 15 • Max/min revolute joint velocities ±1.0 rad/s Max module revolute joint accelerations/decelerations 2.5 rad/s 2 Max/min prismatic joint velocities ±1.0 m/s Max module prismatic joint accelerations/decelerations 2.5 m/s 2

Figure 4 :

 4 Figure 4: The four paths used in the optimization process and a version of Nimbl'bot robot of size: l 1 = 0.5 m, l 2 = 0.5 m, β 1 = π/4, β 2 = π/4, β 3 = π/4, l t = 0.1 m, , β t = π/4.

Figure 5 :

 5 Figure 5: The two paths used to validate the designs obtained from the optimization process and a version of Nimbl'bot robot of size: l 1 = 0.5 m, l 2 = 0.5 m, β 1 = π/4, β 2 = π/4, β 3 = π/4.

Figure 6 :

 6 Figure 6: Box plot of the values assumed by the kinetostatic indices ˆ , η, ν and ρ for the selected candidates. Symbol × represents the mean and • indicates the outliers. The kinetostatic index range is displayed next to the plots.

Figure 8 :

 8 Figure 8: Link lengths of best and worst design for indices ˆ , η, ν and ρ during candidate selection phase.

Figure 9 :

 9 Figure 9: Offset amplitudes of best and worst design for indices ˆ , η, ν and ρ during candidate selection phase.

Figure 10 :

 10 Figure 10: Representation of the design with the best ˆ together with the four trajectories used for the optimization, its design parameter values are in the first line of Table7.

Figure 11 :

 11 Figure 11: Representation of the design with the worst ˆ together with the four trajectories used for the optimization, its design parameter values are in the fourth line of Table7.

Table 1 :

 1 Details about the task names, control objective types, and hierarchy levels. Symbol (E) represents the equality control objective tasks and (I) the inequality ones. The last two columns list the task hierarchies for actions A 1 (Reach Pose) and A 2 (Follow Trajectory). The symbol "/" means that a task is not present in an action.

	Hierarchy level

10 :

 10 Save t i .

	11:	end if
	12:	end while
	13:	end for
	14:

Table 3 :

 3 Trajectory details.

	Steps	560
	Magnitude velocity vector || v||	0.002 m/s
	Magnitude tangential force vector || f t ||	60 N
	Magnitude radial force vector || f r ||	20 N
	Time	1400 s

Table 4 :

 4 Design parameter and kinetostatic performance values ˆ , η, ν and ρ for the best designs obtained from discretization and optimization, respectively.

	Method	l 2	β 2	ˆ	η	ν	ρ	Time
	Discretization 0.1 m -105 • 0.269 0.222 0.087 0.558 120 min
	Optimization 0.06 m -102 • 0.281 0.221 0.086 0.595 20 min

Table 5 :

 5 Machine details.

	Operating System	Linux
	Distribution	Ubuntu 20.04
	CPUs number	4
	CPU model	Intel Core i7 10th Gen, 1.30GHz
	Language	C++
	Control frequency	10Hz
	Time to track one trajectory	5 s
	candidate selection phase takes less time in testing the 76 designs on the desired trajectories.

Table 6 :

 6 Correlation coefficients between kinetostatic performance indices ˆ , η, ν and ρ and design parameters, bottom left, and their standard deviation along the diagonal, gray. The higher the correlation coefficient absolute value, the darker the blue shade. Correlation coefficients ∈ [-1, 1].

		ˆ	η	ν	ρ	l 1	l 2	β 1	β 2	β 3	l t	β t
	ˆ	0.06
	η 0.947	0.09
	ν 0.654 0.579 ρ 0.810 0.770 0.176 0.06 l 1 0.394 0.310 0.823 -0.026 0.41 m 0.05 S t a n d a r d d e v ia t io n l 2 0.326 0.206 0.753 -0.035 0.743 0.4 m
	β 1 0.090 0.048 -0.081 0.280 -0.104 0.100 27.2 •
	β 2 -0.326 -0.336 -0.225 -0.328 0.049 -0.008 -0.261 116.7 •
	β 3 -0.141 -0.148 -0.132 -0.098 0.076 -0.089 -0.081 0.426 52.5 •
	l t 0.155 0.187 0.182 0.082 0.105 0.105 -0.118 0.111 -0.083 0.08 m
	β t -0.016 -0.018 0.092 -0.165 -0.012 -0.079 -0.107 -0.069 -0.213 -0.098 13.7 •

Table 7 :

 7 Design parameters and kinetostatic performance indices ˆ , η, ν and ρ for best and worst designs obtained during candidate selection phase. Green and red colors indicate the highest and lowest values for ˆ , η, ν and ρ obtained during the candidate selection and testing phases.

		l 1	l 2	β 1	β 2	β 3	l t	β t		Candidate selection	Testing
									ˆ	η	ν	ρ	ˆ	η	ν	ρ
	Best ˆ ρ	2.19 m 2.33 m -61 • -131 • -3 • 0.44 m 131 •	0.555 0.695 0.214 0.698	0.556 0.671 0.209 0.718
	Best η	2.05 m 1.99 m -108 • -131 • 17 • 0.44 m 97 •	0.537 0.705 0.259 0.631	0.538 0.703 0.242 0.639
	Best ν	2.94 m 2.78 m -131 • -108 • 38 • 0.44 m 105 •	0.509 0.599 0.372 0.564	0.449 0.526 0.320 0.527
	Worst ˆ η ν 0.07 m 0.31 m 11 • -105 • 1 • 0.44 m 73 •	0.247 0.280 0.018 0.458	0.292 0.350 0.024 0.522
	Worst ρ	1.54 m 2.32 m -131 • 131 • -14 • 0.44 m 105 •	0.309 0.363 0.137 0.435	0.435 0.552 0.176 0.540

Acknowledgments

This paper was supported by the ANRT (Agence Nationale Recherche Technologie) grant CIFRE no 2020/1051 and Nimbl'bot company.

them. On the contrary, the correlation coefficient between ν and ρ is low meaning that different design variables influence these two metrics. However, the correlation coefficients between η, ν and ρ are always positive meaning that these performance indices are not antagonistic. A consideration can be done on the correlation coefficients between β 2 and the indices ˆ , η, ν and ρ, which are negative. This means that the lower β 2 the larger ˆ , η, ν and ρ. Moreover, there is a relevant positive correlation between β 2 and β 3 , meaning that a bigger β 2 is associated to a bigger β 3 . The other design parameters do not present significant correlation coefficients. However, this fact does not mean that any value can be used for these design parameters. For example, the standard deviation of the kinetostatic indices ˆ , η, ν and ρ is small meaning that the obtained results are closer to each other. The plots in Fig. 6 show how each index is concentrated in a smaller area compared to the available one. Considering the design parameters, the standard deviation of l t is much lower than l 1 and l 2 . This means that almost all the obtained designs have a similar value for l t . So, the optimization algorithm identified a small optimal range of values for l t . This can also be noticed by comparing the interquartile size, i.e. the distance between the lower and upper quartile, of the box plots in Fig. 7a. Similarly, β 1 , β 3 and β t have a standard deviation lower than β 2 meaning that the optimization algorithm found a smaller optimal range of values for these parameters. Moreover, the β 2 interquartile size in Fig. 7b is much bigger than the other ones. This can justify why the correlation coefficients between β 1 , β 3 and β t and ˆ , η, ν and ρ are low as well. There is no design that can be identified globally as the best one, but several combinations can lead to a good kinetostatic design for the Nimbl'Bot robot. This is mainly due to the high kinetostatic redundancy of the robot that allows infinite possible solutions for the same problem. Moreover, more than one kinetostatic performance index is considered in the optimization process. So, the designs identified by the optimization process can perform better in dexterity, manipulability or RTR. It is necessary to make a trade-off for choosing the best design.