Concentric Tube Robots (CTR) are becoming more and more popular in medical robotics community. In this thesis, a general literature survey on existing works covering the research topics of CTR is first presented. The kinematics of CTR is more specifically detailed since it is the basics of the main contributions of this thesis. The first contribution is a concept study of exploiting CTR for resection of deep brain tumors located at the frontal lobe. Grid searching has been used as the optimization method for the CTR tubes design. This method allows to avoid the crucial problem of weights preselection which is required in all scalarization methods existing in literature. Instead, the grid searching method used in this work allows to choose the optimal parameters with the help of graphical illustration of calculation results distribution with respect to the selection criteria. The elastic stability dues to the bending and torsion interaction between tubes is considered and evaluated with a new approach introduced in this work. The second contribution then is to deal with the kinematic uncertainties in motion control of CTR. The proposed control method designed at the actuator level shows that the control design of actuator input with task-space feedback and approximate Jacobian matrix provides robustness in handling inaccuracy in kinematic model and maintains good control performance at the same time.

Résumé

Les robots à tubes concentriques deviennent de plus en plus populaires dans la communauté de la robotique médicale. Dans cette thèse, un état de l'art général des travaux existants et qui couvre les thématiques de recherche en robots à tubes concentriques (RTC) est présenté dans un premier temps. Les modélisations géométrique (directe et inverse) et cinématique des RTC sont détaillées car elles servent de base pour les contributions de cette thèse. La première contribution consiste en une étude de concept d'utilisation des RTCs pour la résection des tumeurs profondes situées au niveau du lobe frontal du cerveau. 'Grid searching' a été utilisée comme méthode d'optimisation pour la conception des tubes des RTCs. Cette méthode permet d'éviter le problème crucial de présélection des coefficients de pondération, cette pondération étant nécessaire dans toutes les méthodes de scalarisation existantes dans la litérature. La méthode de 'grid searching' utilisée dans ce travail permet la sélection des paramères optimaux avec l'aide d'une illustration graphique de la distribution des résultats de calcul concernant les critères de séléction. La stabilité élastique due aux interactions des tubes en flexion et en torsion est incluse dans les critères de séléction et est évaluée avec une nouvelle approche introduite dans ce travail. La deuxième contribution de cette thèse repose sur la synthèse d'une loi de commande qui permet de faire face aux incertitudes cinématiques dans le contrôle de mouvement des RTCs. L'étude réalisée a montré qu'un contrôle au niveau des couples moteurs avec un retour dans l'espace opérationnel et une matrice Jacobienne The evolution of surgery from traditional open surgery to minimally invasive surgery (MIS) has witnessed an improvement of the surgical outcomes (reduced patient trauma, less blood-loss, shorter recovery duration, etc). This improvement is due to the miniature incisions compared to large openings on patient body. However, remarkable limitations are associated to this technique. The accessibility to the surgical target passing through complex 3D paths and the control of the rigid instruments used is more difficult (constrained manipulability, dexterity, degrees of freedom, etc) [START_REF] Vitiello | Emerging robotic platforms for minimally invasive surgery[END_REF]. All these have motivated the use of robotic assistance for better benefits of MIS. Robotic assisted surgery started in 2000, with the objective of enhancing the surgical procedure while maintaining the advantages of open surgery in MIS [START_REF] Satava | Future applications of robotics[END_REF]. Robotic systems offer many advantages to MIS procedures: facilitate access to the surgical target with more dexterity and degrees of freedom; enhance the vision of the anatomical structures; enhance the precision of the instrument positioning; eliminate the hand tremor of the surgeon manipulating a teleoperation arm and generate smoother trajectory of the instrument than the trajectory of the hand manipulation ; reduce the fatigue of surgeon thanks to ergonomic studies (comfortable seat position, 3D visual feedback, etc), etc [START_REF] Rosen | Surgical Robotics[END_REF].

CHAPTER 1. GENERAL INTRODUCTION AND MOTIVATION

In MIS interventions, the often used instruments can be classified into three main classes [Dupont et al., 2010a]:

-The first class consists of flexible needles with a straight shaft. They are used for percutaneous insertions in tissues. The needle relays on the lateral forces to steer and change its direction thanks to the guidance forces applied at its base and thanks to its bevel tip.

The limitation of this class of instrument is its dependence on tissues, i.e. it can not steer in curved tubular cavities. Moreover, these instruments are not able to produce lateral tip motion without relying on solid tissues penetration.

-The second class of instruments includes straight and stiff shaft with an articulated tool mounted on the tip. The insertion of the instrument from the body entry to the target should be straight and the manipulation of the tool during the operation is ensured by pivoting around the body entry which on the other hand can cause tissues deformation.

The straight and stiff shaft is difficult to pass through complex cavity as well.

-The last class is composed of elongated steerable devices such as catheters and multistage microrobot devices typically mounted at the distal extremity of a rigid shaft and includes a flexible backbone.

Flexible robots with small size are more suitable for MIS as they can reach surgical sites of high access difficulty with straight rigid instruments [START_REF] Jessica | Continuum Robots for Medical Applications: A Survey[END_REF]. Concentric tube robots belong to this category of flexible robots with the advantage of good steerability compared to rigid devices and controllability of the shape and force compared to catheters with small size. The concept of CTRs is briefly presented in the next section.

Concentric tube robots

CTR is composed of pre-shaped flexible tubes, mostly made of Nitinol. The tubes are inserted concentrically one into another forming a flexible and controllable robot shaft. Usually the tubes used are composed of two parts: a straight part followed by a curved one (cf. Figure

Motivations of the thesis 1.3.1 Medical motivation

It is noted that in 2012, % of the over  patients with brain tumor died in the world, this number of deaths is predicted to reach  in 2020 [IAR]. Anterior deep brain tumors are difficult to access with MIS, and its resection requires a cranium access for good accessibility to the surgical target currently. A part of the cranium is severed to ensure an open surgery then fixed once the tumor is resected. This trans-cranial procedure is very traumatic

CHAPTER 1. GENERAL INTRODUCTION AND MOTIVATION

as the surgeon needs to lift the brain in order to reach the deep tumor located under it. This manipulation causes a swelling of the brain after the operation which creates a high pressure at the brain level.

The concept of an endonasal approach (see Figure 1.2) is studied in this thesis by using CTR which is flexible and adaptable with the anatomical complexity. This operation requires bone drilling in order to access the surgical site. But the insertion path is from the nostril to the anterior wall of the frontal sinus (location of the drilling point) and is surrounded by structures of high sensitivity. The optimization algorithm used is dedicated to this surgical application and aims to obtain the optimal CTR design.

Figure 1.2 -Nasal access to the frontal lobe using CTR.

Technical motivation

For a CTR robot, the shape is defined by the interaction between its tubes. This interaction is created by actuating the tubes in axial rotations and linear translations. Once the tubes interact between each other, they are deformed and result in bending and torsion. The deformation depends on the pre-curvatures of the tubes, their length, their stiffness, and also their translation and rotation. The existing kinematic models which allow to define the robot shape with respect to the rotations and translations of the tubes and taking into account the bending and the torsion between the tubes are complicated, time consuming and not very ac-1.4. ORGANIZATION OF THE THESIS curate. Moreover, most of position control strategies developed for CTR are considered in the joint space, which keeps the position control performances highly dependent on the kinematic modelling accuracy. A task-space controller at the actuator level is developed in this thesis in order to reduce the effect of kinematic inaccuracies and guarantee motion control accuracy at the same time.

Another issue in CTR position control is the bifurcation phenomenon. The interaction between the tubes creates an elastic energy which can be stored and suddenly dissipated for some specific configurations. This leads to tip jump which is not safe for medical interventions. In general, the jumps occur for tubes with high curvatures and large lengths and when the difference of rotation angles of the tubes is in the neighborhood of π. In literature, two solutions are proposed to avoid this problem: optimize the tubes dimensions which reduce the bifurcation, and consider the critical configurations at the level of motion planning to avoid them. In this thesis, a new approach to evaluate the bifurcation effect based on the tip trajectory smoothness is proposed. This bifurcation evaluation index and the reachability of the desired surgical volume are used as criteria of a Pareto grid-searching optimization algorithm for CTR design.

The objective is to calculate the optimal tubes curvatures which maximize the reachability and minimize the proneness to bifurcation. The Pareto grid-searching optimization method used allows an interpretation of the results without any weights pre-selection comparing to the scalarization methods in literature which require a pre-selection of the weights for the different cost sub-functions composing the global cost function.

Organization of the thesis

The thesis is organized as follows:

-In Chapter 2, the existing works of CTR in literature are briefly presented under different categories: the design of CTR; the kinematic modelling of CTR; the motion planing; and the control of CTR. The objective of this chapter is to provide a general framework of CTR literature although this thesis mainly focuses on design optimization and control.

-In Chapter 3, the forward, differential, and inverse kinematic models of CTR are presented in detail since they serve as the basis for CTR design and control.

CHAPTER 1. GENERAL INTRODUCTION AND MOTIVATION

-The Pareto grid-searching optimization method for the design of a CTR robot dedicated to frontal lobe tumor removal and a new method to evaluate the tubes elastic stability are introduced in Chapter 4.

-The chapter 5 presents a new control strategy from the actuator level to compensate the inaccuracies due to the kinematic modelling. Simulation study results are shown to justify the effectiveness of the proposed control method comparing to traditional inverse kinematic method.

-The last chapter concludes the work presented in this thesis and provides some perspectives on future works.

-The prototype of a CTR robot together with its CAD design, motors and sensors selection are presented in the Appendix. The motivation to design concentric-tube robot (CTR) is to fulfill clinical tasks requirements of navigation through complex paths and high manipulation dexterity. To reach this objective, efforts have been made to improve the design of CTR, to develop their appropriate kinematic modelling and their path planers, and to improve the performance of posture control. This chapter presents a general literature survey on existing works covering the aforementioned research topics of CTR [START_REF] Boushaki | Conception, modélisation et commande des robots à tubes concentriques : vers des applications médicales[END_REF].

Introduction

Since CTRs have been introduced with the objective of more contribution and enhancement to robotic assisted surgery in general thanks to their dexterity and miniature-size, they have mainly witnessed four research axes: the design of the actuation unit and of the tubes; the kinematics modelling; the path planing; and the control. The progress in design of the actuation unit was mainly to propose MRI compatible prototypes with different solutions (pneumatic and piezoelectric actuation). The tubes design and path planing works are often presented in a medical context i.e. define the optimal tubes parameters and the optimal joint sequence to navigate through a path while ensuring the sensitive anatomical structures. Many works have been conducted in kinematics modelling to improve the accuracy by considering more realistic physical phenomena and some control strategies have been developed for CTR position control. In this chapter, existing works in literature on all these topics are briefly introduced to provide a general framework of CTR study, though this thesis mainly focuses on the tubes design and motion control problems.

Design of CTR

In terms of mechanical structure, CTRs have special characteristics. They look like serial robot but with flexible and length variable links without discrete articulations, and also they look like parallel robot because the actuators are not distributed along the robot arm. The robot then is composed of two main parts: the actuation unit which ensures the translation and the rotation of the tubes, and the second part is the robot arm composed of concentric tubes. To accomplish a specific task, the integration of other mechanical and/or electronic devices on the robot arm would be required for sensing and manipulating. In this section, existing works on design of both actuation unit and arm (tubes) are presented.

DESIGN OF CTR

Actuation unit design

Since the concept of CTR was introduced, different mechanisms of actuation unit (AU) have been designed. The AU is composed of actuation blocks which carry the tubes in translation and are equipped with rotation actuators. The number of the actuation blocks is equal to the number of the actuated tubes. In terms of repartition of the actuation blocks, CTR AUs can be classified into two categories: AU with serial blocks and AU with parallel blocks. In this part, the majority of existing AUs (serial and parallel) are presented according to:

-the number of actuated tubes;

-the coupling of the rotation and translation motions;

-the number of the robot arms;

-the actuators nature and the MRI compatibility of the robot.

Serial actuation unit

A serial AU is composed of a set of actuation blocks ensuring the rotation and the translation (which can be coupled or not) of each tube. These blocks are distributed in a one followed by another configuration (see Figure 2.1). Here, some typical serial AU are introduced and their actuation mechanism are explained.

For the robot illustrated in Figure 2.2, the AU is composed of three actuation blocks to drive three tubes in rotation and translation. The mechanism of actuation is a nut/screw system with a guidance rail which guides the block during translation. The nut and the screw are rotated independently using one motor for each. For the screw, the motion transmission is ensured by a belt, and the nut is fixed on the motor shaft. The extremity of each tube is fixed on its corresponding nut. Holes are drilled on the nuts to allow the concentric insertion of the tubes.

The actuation of a screw creates a block sliding on the guidance rail, which produces a translation of the tube. The rotation of the nut produces two motions: linear translation and axial rotation of the tube. To create a pure rotation of the tube, both nut and screw actuation are required, but with opposite directions of rotation. This actuation system results in the coupled rotation and translation of a tube, i.e. the actuation of one motor creates rotation and translation motions. The rotation of one tube requires the actuation of two motors. The quantities of translations and rotations of the tubes are indispensable to define the arm shape. They are calculated arithmetically from the motors rotation angles and other parameters (gears, pulleys diameters, screw/nut pitches). This block distribution makes the distance between the base of the innermost tube and the base of the outermost tube too long which requires a long innermost tube. As the torsion deformation of the tube is proportional to its length, the innermost tube is then more prone to this deformation.

The actuators used are DC motors which are not MRI compatible.

The prototype designed by [START_REF] Burgner | A bimanual teleoperated system for endonasal skull base surgery[END_REF] and [START_REF] Burgner | [END_REF] was the first bimanual CTR (Figure 2.3). The two arms are actuated by two parallel Figure 2.2 -CTR with serial actuation blocks and coupled rotation/translation for three tubes [START_REF] Webster | Design and Mechanics of Continuum Robot for Surgery[END_REF] axes. Each axis is composed of a screw which supports three serial blocks. Each block is composed of a motorized nut to drive the block in translation, and another motor to actuate the rotation. When the blocks are translated to their maximal extremities, the space between blocks is reduced thanks to the perpendicular positioning of the blocks with the axes, comparing to the previous prototype in which they are co-axial. Two parallel guidance rails are fixed at the base of the device to support and guide the blocks. The actuators of this MRI non-compatible device are DC motors, and for each DOF only one motor is used. In total the actuation structure includes six motors per arm, with three for the rotations and three for the translations independently without any coupling.

Figure 2.3 -RTC equipped with two three-tubes arms and an actuation serial unit without coupling [START_REF] Burgner | A bimanual teleoperated system for endonasal skull base surgery[END_REF]] [Burgner et al., 2014].

Pneumatic actuators in [START_REF] Cardona | A mri compatible concentric tube continuum robot with pneumatic actuation[END_REF] and piezoelectric actuators in [START_REF] Su | A MRI-guided concentric tube continuum robot with piezoelectric actuation: A feasibility study[END_REF] have been used to ensure the MRI compatibility of CTR. These two prototypes are illustrated The robot AU designed in [START_REF] Cardona | A mri compatible concentric tube continuum robot with pneumatic actuation[END_REF] allows to drive two tubes. The innermost tube is controlled only in translation, but, the middle tube is controlled in both translation and rotation. A passive outermost tube is added. The translations are actuated by pneumatic cylinders controlled through valves (cylinder-valve system). The rotation is also actuated pneumatically with a cylinder connected to a rack and pinion system. The rotation stage is carried by the translation one and the whole block relies on an aluminum guidance rail. This prototype is MRI compatible with no coupling rotation and/or translation. The encoders used for the joint position sensing are optical.

In [START_REF] Su | A MRI-guided concentric tube continuum robot with piezoelectric actuation: A feasibility study[END_REF], the specification of the structure is that the base of the CTR is not Figure 2.5 -Three-tubes CTR actuated pneumatically with translation coupling [START_REF] Comber | Sliding Mode Control of an MRI-Compatible Pneumatically Actuated Robot[END_REF].

Parallel actuation unit

The blocks of a parallel AU are distributed and translated in parallel with each other, as show in Figure 2.6.

The first parallel AUs have been reported by Dupont et al. [Dupont et al., 2010a] and Gosline et al. [Gosline et al., 2012a] with both prototypes having the same mechanical structure. The particularity of these prototypes compared to the serial ones is the parallelization of the actuation blocks. Three parallel nut/screw axes are used for the translations, with one motor per nut fixed at its extremity. Theses axes are used for the tubes translations, and each translated screw carries with it a motor dedicated to the rotation. The rotation motion is transmitted through a belt in [Gosline et al., 2012a] (see Figure 2.7.a) or a gearing in [Dupont et al., 2010a] (see Figure 2.7.b). The motions of the tubes are decoupled with each motor creating only one motion for one tube: a linear translation or an axial rotation. Another advantage of this mechanism is the reduction of the tubes lengths as they are actuated by parallel blocks rather than serial ones. The presented devices actuate a single arm with three tubes using DC motors and hence is not MRI compatible. In [START_REF] Xu | Position control of concentrictube continuum robots using a modified Jacobian-based approach[END_REF], two motorized linear slides are assembled perpendicularly on the base [START_REF] Xu | Position control of concentrictube continuum robots using a modified Jacobian-based approach[END_REF] As the sterilization of medical devices is an important issue, Hendrick et al. have isolated the actuators from the transmission stage with a sterile barrier [START_REF] Hendrick | A multi-arm hand-held robotic system for transurethral laser prostate surgery[END_REF]. The

CHAPTER 2. FROM DESIGN TO CONTROL OF CONCENTRIC TUBE ROBOTS: A LITERATURE SURVEY

designed robot is a two-arm CTR inserted in a common sheath (see Figure 2.9). The first arm is a two-tube robot with three DOFs. The first tube is actuated in both rotation and translation, and the other is actuated only in translation. The second arm has three tubes with six DOFs.

All the tubes are actuated in both rotation and translation. The tubes are fixed on the translated blocks, and the rotation motion is transmitted through a rotated square shaft. The transmission part is designed to be biocompatible and autoclavable. This robot is also particular in the sense that it is intended to be a hand-held robot.

a) b)

Figure 2.9 -a) Two arms CTR with b) parallel actuation unit [START_REF] Hendrick | A multi-arm hand-held robotic system for transurethral laser prostate surgery[END_REF].

The serial AU is more oftenly used than the parallel as the latter is more complicated in design. Despite the particularity of each prototype, their skeletons resemble each other. The structure for the majority of the prototypes is mainly composed of translating blocks thanks to a nut/screw system and supporting rotation motors. The rotation is usually transmitted to the tube through a belt or gearing. The actuation characteristics of most existing CTR prototypes are summarized in the Table .2.1. A new actuating technique of CTR has been introduced in [START_REF] Chikhaoui | Kinematics and performance analysis of a novel concentric tube robotic structure with embedded soft micro-actuation[END_REF] based on embedded soft micro-actuation.

Tubes design

The arm of a concentric tube robot is a set of pre-shaped super-elastic tubes combined concentrically. The tubes used are usually composed of two sections: a straight section with an extremity fixed on its corresponding rotary axis and a distal pre-curved section. Once the tubes are inserted inside each other and actuated in rotation and translation, elastic interactions stiffness, etc). This motivates the optimization works introduced in this section, which mainly aim to define the optimal lengths and curvatures of the tubes for a given clinical application.

Before design optimization algorithms were developed, a first guideline of the tube specifications selection has been proposed by Dupont et al. [Dupont et al., 2010a]. The objective of the general rules proposed in that work is to select the specifications of the tubes which satisfy two criteria: the ability to manipulate distal and proximal sections independently, and the ability to navigate through complex curved 3D paths. It is also important that the robot has the capacity to steer in tissues by exerting minimum lateral forces. Mainly, these rules explain how to select the stiffnesses and the curvatures of the tubes but in a general manner.

The order of the robot sections is also defined in term of stiffness, curvature and extension (in descending order from the proximal to the distal section for the stiffnesses and in ascending order for the curvatures). The work provides only general rules and does not provide optimal definition of the tube parameters.

After this proposed general guideline, several algorithms have been proposed to optimize the tube dimensions. In order to reduce the complexity of the optimization problem, the number of optimization variables was restricted. In most of the works, the only considered parameters of tubes are the curvatures and the lengths. The number of tubes and their diameters are pre-selected and the only tube pre-shape considered is the circular form except the work of [START_REF] Gilbert | Can concentric tube robots follow the leader[END_REF] where the helical form is considered.

The three first design optimization algorithms have been introduced in the same period [START_REF] Bedell | Design optimization of concentric tube robots based on task and anatomical constraints[END_REF][START_REF] Anor | Algorithms for design of continuum robots using the concentric tubes approach: A neurosurgical example[END_REF][START_REF] Burgner | A bimanual teleoperated system for endonasal skull base surgery[END_REF]. The work of Bedell et al.

[ [START_REF] Bedell | Design optimization of concentric tube robots based on task and anatomical constraints[END_REF] aims to define the optimal lengths and curvatures for the entire robot sections. However, the number of tubes and sections and their nature (variable or fixed curvature)

and order are all preselected. In that work, the optimization problem is decomposed into two sub-problems: navigation optimization problem and manipulation optimization problem. The navigation part of the robot is used to position the manipulation part at the surgical site. Once the surgical site is reached by the manipulation part, the navigation part is held and only the manipulation part moves to accomplish the task (resection, injection, etc). The optimization algorithm used is Generalized Pattern Search (GPS) [START_REF] Torzcon | On the Convergence of Pattern Search Algorithms[END_REF], and the cost function of the optimization problem is defined to achieve that: the tip position stays inside the de-sired workspace limited by the anatomical constraints and the length and the curvatures of the sections are minimal. This work has served later to optimize the design of a robot used in intracardiac beating heart surgery in [Gosline et al., 2012b].

The second design optimization work has been proposed in [START_REF] Anor | Algorithms for design of continuum robots using the concentric tubes approach: A neurosurgical example[END_REF] which addresses another medical application. From a given set of discrete target points, and knowing the anatomical constraints, the algorithm allows to find the simplest robot structure with the minimum required number of sections. It also calculates the length and the curvature of each section. All these criteria are formulated by a set of cost functions which are minimized using Pattern Search Algorithm (PSA) [START_REF] Audet | Analysis of Generalized Pattern Searches[END_REF]].

These two presented works require an inverse kinematics calculation to avoid the anatomical constraints in the optimization procedure. The inverse kinematic calculation for CTR is not trivial. It is complicated and time consuming as will be explained in the next section.

The third work has been introduced by Burgner et al. [START_REF] Burgner | A bimanual teleoperated system for endonasal skull base surgery[END_REF], which allows to avoid the inverse kinematics calculation. Only the forward kinematics calculation is required. The application considered is the trans-nasal skull base surgery. After geometrical description of the workspace based on computed tomography (CT) images, the workspace is discretized to allow an evaluation of the robot within reachability. The reachability is evaluated based on the percentage of the reached volume in the desired workspace. It is quantified by the convex hull of the discrete points reached in the surgical workspace which is considered as a cost function of the optimization problem. The optimization variables are the lengths and the curvatures of the tubes, and the optimization algorithm uses Nonlinear Simplex Method.

The kinematic models used in the aforementioned works are simplified ones. In [START_REF] Anor | Algorithms for design of continuum robots using the concentric tubes approach: A neurosurgical example[END_REF] the infinite stiffness model [START_REF] Furusho | Curved multi-tube systems for fetal blood sampling and treatments of organs like brain and breast[END_REF] was used which neglects the interaction bending and torsion between tubes, and in both works of [START_REF] Bedell | Design optimization of concentric tube robots based on task and anatomical constraints[END_REF] and [START_REF] Burgner | A bimanual teleoperated system for endonasal skull base surgery[END_REF] the bending only model was used which neglects the tubes torsion. The advantage of using simplified models is to reduce the computation time and the complexity of the optimization problem formulation. This simplification comes at the cost of accuracy. For this reason, more accurate kinematic models have been used later in the design optimization algorithms.

Torres et al. [START_REF] Luis G Torres | Task-oriented design of concentric tube robots using mechanics-based models[END_REF] use a more complete kinematic model which considers the torsion of the tubes [START_REF] Rucker | Mechanics of Bending, Torsion, and Variable Precurvature in Multi -Tube Active Cannulas[END_REF] The limitation of the presented works so far is the fact that they only consider a discrete set of points to be reached. As an alternative, a volume-based approach was proposed by Burgner et al. in [Burgner et al., 2013a] which uses Nelder-Mead Simplex algorithm [START_REF] Lagarias | Convergence Properties of the Nelder-Mead Simplex Method in Low Dimensions[END_REF]] to calculate the optimal tubes dimensions. This work offers more coverage of the surgical workspace compared to the previous works. The kinematic model used in this work

is the torsionally compliant model described in [Rucker et al., 2010b].

All the aforementioned works consider only the case of piecewise-constant pre-curvatures.

However, there are other options that could be explored. Junhyoung et al. [START_REF] Ha | Achieving elastic stability of concentric tube robots through optimization of tube precurvature[END_REF] have proposed to utilize a variable pre-curvature of the tube along the length. This work shows that the obtained optimal curvatures allow to avoid bifurcation, which is defined in the end of the next chapter. The optimization problem is formulated by maximizing a cost function and studying its monotony.

In order to study the possibility of steering the robot in a follow-the-leader way, i.e. the shaft of the robot keeps its form once the tip progresses forward, Gilbert et al. [START_REF] Gilbert | Can concentric tube robots follow the leader[END_REF] explored this concept for the first time and showed that this deploying manner is possible only for very specific tube pre-curvatures.

The helical form as another alternative of the tube pre-shaping has been proposed as well.

The mathematical conditions of follow-the-leader deployment and of elastic stability have been formulated in that work. The two cases of curved and helical tubes have been studied.

The extension of this work with experimental results have been presented in [START_REF] Gilbert | Concentric tube robots as steerable needles: Achieving follow-the-leader deployment[END_REF].

The latest work in design optimization [START_REF] Bergeles | Concentric tube robot design and optimization based on task and anatomical constraints[END_REF] has introduced an optimization framework based on specifications of the task and the anatomical structure to generate the robot tube parameters. The elastic stability was considered in the optimization problem and the possibility of follow-the-leader insertion was studied. The optimization algorithm used is Nelder-Mead downhill simplex method [START_REF] Nelder | Asimplex method for function minimization[END_REF], and it was tested with two different clinical applications: intra-cardiac and neuro-surgeary.

The works introduced in this section are summarized in Table .2.2. J, 2009b]. In this section, the main existing kinematic models of CTR in literature are introduced. The mathematical formulation of each model will be detailed in the next chapter.

KINEMATIC MODELS OF CTR

• The first model used considers the stiffness of the outer-most tube of the robot as infinite compared to all the inner tubes (Infinite Testiness Model) [START_REF] Furusho | Curved multi-tube systems for fetal blood sampling and treatments of organs like brain and breast[END_REF]. The limitation of this model is that if the stiffness of the outermost tube is not high enough to satisfy the assumption then the accuracy of the model will degrade. In the case of very high stiffness of the outermost tube compared to the inners, the shape variation and the dexterity of the robot are limited. otherwise, the shape variation is limited with less dexterity.

• Torsion-Free Model (TFM) takes into account the bending interaction between tubes, but the torsion of tubes is neglected. As the tubes are made of Nitinol, the torsion impact on the rotation angles deformation of the tubes can be significant in some cases which affects the estimation of the shape [START_REF] Sears | A Steerable Needle Technology Using Curved Concentric Tubes[END_REF]].

• Kinematics modelling of CTR is improved by including the transmissional torsion between actuators and the first curved link in the Transmissional Torsion Model (TTM) which reduces the error caused by torsion phenomenon in the straight part, but the torsion in the curved part of the tubes is neglected [START_REF] Webster | Kinematics and calibration of active cannulas[END_REF] [R. J. Webster et al., 2009].

• The Torsionally Compliant Model (TCM) is more accurate than the prior works, which is based on Cosserat Rod Theory (Special Cosserat Rod Equilibrium Model) [START_REF] Cosserat | Théorie des corps déformables[END_REF]Cosserat, 1909][Antman, 1995]. Two approaches have been used to calculate the forward kinematic model: the Newton equilibrium wrench in [Dupont et al., 2010a] [Rucker et al., 2010b] and the Euler-Lagrange energy minimization in [START_REF] Dupont | Torsional kinematic model for concentric tube robots[END_REF]] [Rucker et al., 2010a]. The torsions of the straight and curved parts of the tubes are both considered in the model, but the latter is described by a second order nonlinear differential equation causing problems of boundary conditions and heavy computational burden in real-time. Another limitation is that the frictions are not included in this model. For this reason, this model was extended by considering the effect of the friction between tubes in the kinemato-static formulation [START_REF] Lock | Friction Modeling in Concentric Tube Robots[END_REF]].

• In order to reduce the computation time, [START_REF] Xu | A fast torsionally compliant kinematic model of concentric-tube robots[END_REF] reformulate the torsionally compliant model and linearize it, then they propose to use a force sensor to solve the unknown initial conditions problem: the torsion efforts of the tubes at the base are 2.4. MOTION PLANNING 29 unknown (more details are provided in the next chapter). However, the linearization may reduce the accuracy and the force sensors required for each tube are cumbersome.

• A model which considers the deflection of the robot tip in interaction with the environment has been developed in [START_REF] Mahvash | Stiffness Control of Surgical Continuum Manipulators[END_REF] and used for stiffness control. This model is also based on Cosserat Rod model.

Motion planning

In order to orient the robot to a specific target with the presence of obstacles (sensitive structures), a motion planner based on constrained environment is required. This planner calculates the optimal actuation sequence which allows the robot tip to move to the desired position. The calculation of the optimal actuation sequence is highly dependent on the robot kinematics. Motion planning algorithms for CTR have been developed and evolved with the kinematics evolution of CTRs and the algorithms have been dedicated to different medical applications. In this section, existing motion planning algorithms for CTRs are presented, and the specification of each work is discussed.

The first work in motion planning of CTR has been introduced by Lyons et al. [START_REF] Lyons | Motion planning for active cannulas[END_REF]. This algorithm calculates the optimal translations and rotations of the tubes which guide the robot to reach its target with avoidance of the anatomical obstacles (spherical obstacles). The desired translations and rotations are calculated based on the target position, the initial tip position and orientation, and the geometric representation of the obstacles.

The motion planning problem is formulated as a constrained nonlinear optimization problem which is then transformed into a sequence of unconstrained optimization problems using penalty method. The cost function formulates the obstacles avoidance and the target reaching with constraints imposed on the orientation angles and section length. This work utilizes a basic kinematic model which neglects the beam mechanics [START_REF] Furusho | Curved multi-tube systems for fetal blood sampling and treatments of organs like brain and breast[END_REF].

In [START_REF] Lyons | Planning active cannula configurations through tubular anatomy[END_REF], another optimization-based planing algorithm has been proposed.

It computes the configuration of the robot which navigates through tubular structures to reach a target. This work is an extension of the previous one by providing two main improvements which are the calculation of the robot shape by minimizing the interaction energy between the tubes, and the consideration of the motion through tubular structure rather than motion with

CHAPTER 2. FROM DESIGN TO CONTROL OF CONCENTRIC TUBE ROBOTS: A LITERATURE SURVEY

avoiding spherical obstacles which is a more realistic scenario. The inputs of the algorithm are similar to the previous work. However, the kinematic model used in the formulation of the optimization problem is a beam mechanics model [R. J. Webster et al., 2009] which is more accurate compared to the previous one. This algorithm has been simulated using extracted human CT scan of the lung bronchi for the lung biopsy navigation scenario.

In [START_REF] Torres | Motion planning for concentric tube robots using mechanics-based models[END_REF], Torres et al. have used the Rapidly-Exploring Roadmap (RRM) [START_REF] Alterovitz | Rapidly-exploring roadmaps: Weighing exploration vs. refinement in optimal motion planning[END_REF] to calculate the control inputs that guide the robot to the target with minimal collision probability. Unlike previous works, this motion planner uses TCM model which considers the bending and the torsion of the tubes in both straight and curved parts. The clinical application addressed in this work is neurosurgery.

Another work which solves the combined problem of design optimization with motion planing has already been introduced in the previous section [START_REF] Luis G Torres | Task-oriented design of concentric tube robots using mechanics-based models[END_REF].

In order to avoid the elastic instability which may cause problems of control, a specific motion planner has been developed by [START_REF] Bergeles | Planning stable paths for concentric tube robots[END_REF]. This algorithm aims to find an elastically stable path respecting the anatomical constraints for robots which exhibit only local stability using RRT* algorithm [START_REF] Karaman | Incremental sampling-based algorithms for optimal motion planning[END_REF].

The latest work so far [START_REF] Luis G Torres | Interactive-rate Motion Planning for Concentric Tube Robots[END_REF] in path planning for CTR uses TCM model and enables a collision-free navigation of the robot tip and shaft. The proposed approach combines off-line pre-computation of a collision-free road-map with on-line position control.

The road-map is built using modified Rapidly Exploring Random Graph (RRG) [START_REF] Karaman | Sampling-based algorithms for optimal motion planning[END_REF]. The medical scenario of this work is neuro-surgery via a nasal access.

All these works are summarized in the following Table .2.3.

Control of CTR

The complexity and hence high computation time of the CTR forward and inverse kinematics makes its real-time control difficult. It is obvious that the control performance depends on the accuracy of the kinematic model used. For CTR, the precision of the direct and inverse kinematics modelling is inversely proportional to the total length of the robot [START_REF] Rucker | THE MECHANICS OF CONTINUUM ROBOTS[END_REF].

The existing works in control of CTR are presented in this section and are summarized in is considered as the robot length); the implementation in real-time and the cycle frequency.

The performances of each architecture in Table 2.4 are evaluated using the Root Mean Square Error (RMS).

In 2009, Webster et al. [Webster et al., 2009] have proposed to control the position of

CHAPTER 2. FROM DESIGN TO CONTROL OF CONCENTRIC TUBE ROBOTS: A LITERATURE SURVEY

CTR through visual servoing. In this work, one rotation joint has been fixed, which reduces the dimension of the Jacobian matrix and makes it square as only three actuators have been used: two translations and one rotation. The inverse of the Jacobian matrix has been used to transform the Cartesian errors to the joint errors.

Dupont et al. in [Dupont et al., 2010b] and [Dupont et al., 2010a] have proposed a teleoperated control architecture of the position and orientation angle of the CTR tip. The Cartesian target given by the manipulator from the master side is interpreted to joint target through its inverse kinematic model. Because of the kinematics complexity, and consequently its heavy calculation requirement, an offline approximation has been carried out. A lookup table containing the joint variables and their corresponding Cartesian poses has been used to identify a Fourier series function which models the forward kinematic model. The inverse kinematic model is then a root-finding problem of the approximation function. It is solved using Gauss-Newton algorithm. The control at the slave level is ensured by a PID controller in the joint space, which makes the control performances sensitive and highly depending on the kinematics modelling errors. The RMS error of the forward kinematic model used that works is about mm. Moreover, the approximation function affects the forward kinematic modelling accuracy and consequently the inverse model.

The stiffness control has been applied to control the contact with a soft environment by [START_REF] Mahvash | Stiffness Control of Surgical Continuum Manipulators[END_REF]. The system is teleoperated and the master controller calculates the desired forces to be exerted by the tip (3 DOFs). This calculation is based on the error between the feedback tip position and the interpretation of the teleoperator hand motion. Through a deflection model, the desired forces are used to calculate the desired joint positions. The joint control is then done at the slave level using a PD controller.

The reference position in [START_REF] Burgner | A bimanual teleoperated system for endonasal skull base surgery[END_REF] and [START_REF] Burgner | [END_REF] is obtained by interpreting the hand motion of the teleoperator at the master side. The desired joint positions corresponding to the desired Cartesian position are obtained through the inverse kinematics calculation. The approach used for the inverse kinematics calculation is the damped least square [START_REF] Wampler | Manipulator Inverse Kinematic Solutions Based on Vector Formulations and Damped Least-Squares Methods[END_REF] of the Jacobian matrix. Based on this method, the joint velocities obtained ensure the avoidance of the singular configurations.

Sliding mode for CTR position control has been used in [START_REF] Comber | Sliding Mode Control of an MRI-Compatible Pneumatically Actuated Robot[END_REF]. This control architecture provides good control performances, but it has been studied at the level of joint tracking (rotations and translations). The control performances of the robot tip Cartesian position has not been evaluated in this work.

The fast kinematic model proposed in [START_REF] Xu | A fast torsionally compliant kinematic model of concentric-tube robots[END_REF] has been derived to calculate the differential kinematic model. The obtained Jacobian matrix has been used for an inverse kinematics control of CTR in [START_REF] Xu | Position control of concentrictube continuum robots using a modified Jacobian-based approach[END_REF]. Later, this model has been improved in order to control the position of CTR under external loads [START_REF] Xu | Real-time Trajectory Tracking for Externally Loaded Concentric-tube Robots[END_REF].

In order to avoid elastically instable configurations and the joint limits (tube translation limits), Azimian et al. in [START_REF] Azimian | Closed-loop inverse kinematics under inequality constraints : Application to concentric-tube manipulators[END_REF] have proposed an inverse kinematic controller under inequality constraints. All the aforementioned works are summarized in Table 2.4.

Introduction

The kinematic modelling is required for both tube design optimization and control. For CTR, the forward kinematic model is a mapping from the tubes translations and rotations to the tip pose passing by a set of sampling points from the robot shaft. The inverse kinematic model is to find the rotations and translations of the tubes given a pose in task space. The differential kinematic model defines the relation between the spatial tip velocities and the rotations and translations velocities of the tubes thanks to the Jacobian matrix.

The forward kinematic model is used for the reachability and elastic stability evaluation in the Chapter 4 which discusses the tubes design optimization. The differential kinematic model is also employed in the task-space control scheme proposed in Chapter 5. The kinematic modelling developed in this chapter is therefore to be used in the next two chapters. The CTR considered in Chapter 4 includes three tubes which are required to ensure the navigation and the surgical workspace coverage. In Chapter 5, two-tube CTR is used to validate the proposed controller through simulation. The principal of kinematic modelling is similar for CTRs with two tubes and three tubes. The difference lies in the fact that the mathematical formulation complexity is higher in the case of three tubes. For clarity of presentation and ease of understanding, the case of two tubes with a straight and a curved part of each tube is considered in all the kinematics modelling presented in this chapter.

Forward kinematic model of concentric tube robot

In this section, the existing forward kinematic models for CTR with their mechanics principles used are introduced starting from the model neglecting interactions between tubes to the model including bending and torsion interactions between tubes. In this order, the modelling accuracy, formulation complexity and calculation time increase accordingly.

Infinite stiffness model

The simplest and least accurate forward kinematic model of CTR was proposed in [START_REF] Furusho | Development of a curved multi-tube (CMT) catheter for percutaneous umbilical blood sampling and control methods of CMT catheters for solid organs[END_REF] and [START_REF] Smith | Computed tomography-based tissue-engineered scaffolds in craniomaxillofacial surgery[END_REF] and assumes that the bending and the torsion due to the interactions between tubes are neglected. The robot arm is split into a set of constant

FORWARD KINEMATIC MODEL OF CONCENTRIC TUBE ROBOT

37 curvature sections. The geometrical shape of a section is assumed to be identical to the preshape of the outer-most tube surrounding it regardless the number and the pre-shapes of the inner tubes. As it has been explained in [R. J. Webster and Jones, 2010] and illustrated in

f  f  (α, ι) = f (θ ) (ρ, l, ϕ) = f  (α, ι) X = f  (ρ, l, ϕ)
Joint space (actuators)

Configuration space (sections)

Task space (end-effector) Figure 3.1 -Tow mappings defining the forward kinematic model: f  from the joint space (α i , ι i ) to the configuration space (ρ j , l j , ϕ j ) and f  from the configuration space to the task space (tip pose X). This figure is adapted from [R. J. Webster and Jones, 2010]. • The first mapping f  is to split the robot into a set of constant curvature sections, define the dominant tube in each section, and then to calculate the configuration parameters of each section. The configuration parameters of a section j are: the curvature ρ j ; the length l j ; and the orientation angle ϕ j ( j = , , ...m with m being the number of sections) with respect to the previous section j -. They are calculated in terms of the joint parameters (rotations α i and translations ι i , i = , , ...n with n is the number of tubes) which depend directly on the actuator angles θ k (k = , , ... n). The number of sections m and their lengths l j are defined in terms of the tubes translations ι i , the lengths of straight parts L S i and curved parts L C i of the tubes. The orientation angle ϕ j of a section j is obtained by subtraction of the dominant tube rotation angle from the dominant tube rotation angle of the previous section j -.

Consider the case of two tubes inserted concentrically with the configuration showed in figure 3.2.The limits of the sections are defined by the start and the end points of the straight and the curved parts (L S i , L C i ) of the dominant tube at each section. The number of sections CHAPTER 3. KINEMATIC MODELLING OF CONCENTRIC TUBE ROBOT

l  l  l  H b H 1 z  x  H 2 z  x  H 3 z  x  Ht z t x t

FORWARD KINEMATIC MODEL OF CONCENTRIC TUBE ROBOT

39

The lengths of the sections are calculated as follows:

       l  = L S l  = L C l  = (L S + L C ) -(L S + L C ) + (ι  -ι  ) (3.1)
The tubes are numerated from the outer-most to the inner-most as shown in Figure .3.2.

The orientation angles of the three sections are given by:

       ϕ  = α  ϕ  =  ϕ  = α  -α  (3.2)
where α i is the Z axial rotation angle of the tube i which is assumed to be constant along the whole tube length since the torsion effect is neglected.

The curvature of a section is equal to the pre-cuvature of the outer-most tube present in it:

       ρ  = O ρ  = /r  ρ  = /r  (3.3)
where ρ j is the curvature of the section j and r i is the curvature radius of the tube i.

• The second mapping f  is to use the configuration parameters (l j , ϕ j , ρ j ) to calculate the homogeneous transformation matrices between the section frames. The position and the orientation of the end effector are then obtained in the base frame by a successive product of these matrices.

The first mapping expressed by the homogeneous matrix b H 1 ∈ SE() (SE() is the special Euclidean group) is then a pure rotation (Rot z ) around the Z axis of the base frame given by:

b H 1 =        [Rot z (ϕ  )]             (3.4) CHAPTER 3.

KINEMATIC MODELLING OF CONCENTRIC TUBE ROBOT

For the first straight section, the mapping from the frame H  to the frame H  attached to the base of the second section is given by a linear translation along the Z axis with a quantity equal to the section length. The  H  ∈ SE() homogeneous matrix expressing the mapping from (H  ) the base of the first section to (H  ) the base of the second section is given by:

1 H 2 =                  l            (3.5)
In the case of a curved section, the mapping from its base to its distal extremity is defined by a rotation ϕ around the Z axis of the base frame followed by a second transformation as shown in Figure 3.3. The second transformation defines the coordinates of a point from the planar section of radius r in the frame (X ′ ,Y ′ , Z ′ ), then the translation between the former and the section base frame is considered. When the rotation ϕ = , the coordinates of the distal extremity point as can be seen from the Figure .3.3 are: (, r cos(lρ)r, r sin(lρ)) and the orientation of the frame attached to distal extremity point of the section is expressed by a rotation lρ around X where the curvature ρ =  r and l is the section length. The mapping  H  therefore is expressed as follows:

2 H 3 =        [Rot z (ϕ  )]                    [Rot x (l  ρ  )] (cos(l  ρ  )-) ρ  sin(l  ρ  ) ρ            (3.6)
Similarly for the last section we obtain: The forward kinematic model is then given by the product of these transformation matrices in this order:

3 H t =        [Rot z (ϕ  )]                     [Rot x (l  ρ  )] (cos(l  ρ  )-) ρ  sin(l  ρ  ) ρ             (3.7) ϕ X Y Z X ′ Y ′ Z ′ l r rsin(lρ) r - rcos(lρ) (ϕ = ) (ϕ>)
b H t = b H 1 1 H 2 2 H 3 3 H t = j j-1 H j (3.8)

Bending model

The model proposed in [R. J. Webster et al., 2009] and in [START_REF] Sears | A Steerable Needle Technology Using Curved Concentric Tubes[END_REF] takes into consideration the bending interaction between tubes. This bending depends on the tubes pre-curvatures, their stiffnesses, and their rotation angles. The outer-most tube is not considered as the dominant tube in terms of stiffness and all the existing tubes in one section are used to calculate its overall curvature and planar angle. Consequently, the number of sections increases compared to the previous model. In our case of study, the section  in figure 3.2 is split into two sections as the pre-curvature of the inner tube changes from straight to curved part, and the bending interaction between tubes is considered rather than considering the outer tube as the dominant stiffness tube in the previous model. The sections are then delimited by the straight and curved parts of all the tubes, as shown in figure 3.4.

l  l  l  l  H b H 1 z  x  H 2 z  x  H 3 z  x  H 4 z  x  Ht z t x t
Figure 3.4 -Two-tube CTR with four constant curvature sections, from the base to the tip: the first section is composed of two straight tubes; the second section is composed of an outer curved tube and an inner straight tube; the third section is composed by two curved tubes and the fourth section is a single curved tube.

The total number of sections is four, and for each section j three parameters are required to define the homogeneous transformation matrix between its distal and proximal extremity.

FORWARD KINEMATIC MODEL OF CONCENTRIC TUBE ROBOT

43

These parameters are the length l j , the curvature ρ j , and the orientation angle ϕ j with respect to the previous section j -. The lengths of the sections are expressed by:

               l  = L S l  = (L S + ι  ) -(L S + ι  ) l  = (L S + L C + ι  ) -(L S + ι  ) l  = (L S + L C + ι  ) -(L S + L C + ι  ) (3.9)
For the orientation angle and curvature of a section j, they are obtained after decomposition of the section curvature into two coordinates in the (X,Y ) plan [START_REF] Webster | Design and Mechanics of Continuum Robot for Surgery[END_REF]. These (X,Y ) coordinates, denoted by (χ j , ν j ), are expressed in the base frame of each section. They are obtained from the Bernoulli-Euler beam mechanics and are given by:

χ j = i E i I i ρi, j cos α i, j i E i I i , ν j = i E i I i ρi, j sin α i, j i E i I i (3.10)
where E i is the modulus of elasticity (Young's modulus), I i is the second moment of area and the ρi, j is the pre-curvature of the tube i in the section j. The limits of the sums are defined by the number of the tubes present in the section. Then, the orientation angle and curvature are calculated by using these (χ j , ν j ) coordinates :

ϕ j = tan - ν j χ j ρ j = χ  j + ν  j (3.11)
For each section j, the curvature does not vary along the whole section. The transformation matrix of a section defining the mapping from its distal to its proximal extremity expressed in the proximal frame is given by the following relationship:

j-1 H j =        [Rot z (ϕ j -ϕ j- )]                     [Rot x (l j ρ j )] (cos(l j ρ j )-) ρ j sin(l j ρ j ) ρ j            (3.12)
and the pose of the tip expressed in the base frame is given by the product of these transformations:

b H t = j j-1 H j (3.13)
where j variates from the base to the tip to include all the sections.

Straight torsion model

The flexibility of the tubes' material leads to torsion problem of the tubes which affects their rotation angles and makes them change from a section to another. In the model of [R. J. Webster et al., 2009], the tube torsions are included in the first straight part of the tubes before the first curvature. This model allows to calculate the rotation angles β i of all the tubes (i = ...n) after the torsion deformation.

In the considered case study, the angles β i represent the rotations of the tubes at the distal extremity of the first section (or the base of the second section). The quantity of the elastic energy developed during the tubes interactions is due to their torsions and flexures. The bending energy ϒ bend stored in a given section j is the sum of the bending energy along the X direction (ϒ x, j ) and the bending energy along the Y direction (ϒ y, j ):

ϒ bend j = ϒ x, j +ϒ y, j , (3.14)
where:

ϒ x, j = i E i I i l j  (χ j -ρi, j cos β i )  , ϒ y, j = i E i I i l j  (ν j -ρi, j sin β i )  (3.15)
As has been previously explained, the only torsion energy included in the elastic energy is the one before the distal extremity of the first straight section given by:

ϒ tor = n G n J n L i (α i -β i )  (3.16)
where L i is the length from the tube base to the end of the first section, G i is the shear modulus and J i is the polar moment of inertia. The total elastic energy is the sum of the torsion energy and the bending energy of all sections:

ϒ = ϒ bend +ϒ tor .
(3.17)

For the torsion energy of the present case study, it is given by:

ϒ tor = G  J  L  (α  -β  )  + G  J  L  (α  -β  )  (3.18)
By applying the energy calculation (3.15) for the second and third section yields:

ϒ bend = l  ρ  E  I  E  I  (E  I  + E  I  ) + E  I  E  I  l  (E  I  + E  I  ) ( ρ  + ρ  - ρ ρ cos(β  -β  )) (3.19)
which expresses the total bending energy stored in the robot.

Once the tubes interact with each other, the resulting shape is found such that the elastic energy ϒ is locally minimal. The angles β i are calculated at the equilibrium by setting the gradient ∇ϒ with respect to β i to zero which leads to the following system:

       - G  J  L  (α  -β  ) + ρ ρ E  I  E  I  E  I  + E  I  l  sin(β  -β  ) =  - G  J  L  (α  -β  ) - ρ ρ E  I  E  I  E  I  + E  I  l  sin(β  -β  ) =  (3.20)
This system can be solved using a variety of numerical solvers. The one that has been used in this work is the Newton-Raphson algorithm [START_REF] Doron | Introduction to Numerical Analysis[END_REF]. After obtaining the angles β i , they are used in (3.10) and (3.11) to calculate the resulting curvatures and orientation angles of each section after including the torsion as has been explained earlier. Later, the obtained curvatures and orientation angles are used in the formula (3.12) for each section, and the pose of the tip is calculated using the product (3.13).

Torsionally compliant model

In order to improve the accuracy of the kinematics, the torsion effect in both straight and curved parts of the tubes has been included in the models introduced in [START_REF] Dupont | Torsional kinematic model for concentric tube robots[END_REF] [ Rucker et al., 2010a] [Dupont et al., 2010a] [Rucker et al., 2010b] instead of considering it only in the straight part in the previous model. The model proposed in [Dupont et al., 2010a] is derived from the special Cosserat rod model [START_REF] Cosserat | Théorie des corps déformables[END_REF]. The curvature ρ i ∈ R  of a tube i has three components (ρ i,x , ρ i,y , ρ i,z ) along (X,Y, Z) axes where the first two are due to the bending, and the last one is due to the torsion. For previous models, the curvatures are defined section by section and assumed to be constant along the whole section. However, for this model, the curvatures vary continuously according to the robot

length s ∈ [, L].
As torsion is considered, the rotation angles α i (s) of the tubes depend on the distance s along the robot shaft from the robot base. The relative twist angle between the two tubes is defined by:

γ(s) = α  (s) -α  (s) (3.21)
The bending curvature components (ρ ,x , ρ ,y ) and (ρ ,x , ρ ,y ) of the first and second tube respectively are calculated after defining the equilibrium moment conditions and ensuring the same bending shape of tubes once assembled. After the mathematical formulation of these conditions presented in [Dupont et al., 2010a], it leads to:

ρ  | x,y = (Rot z (γ(s))ρ  )| x,y ρ  | x,y = (K  + K  ) - (Rot T z (γ(s)K  ρ (s) + K  ρ (s))| x,y (3.22)
where ρ i | x,y are the (X,Y ) curvature components of the tube i and

K i = diag(K i,x , K i,y , K i,z ) =     E i I i    E i I i    G i J i     (3.23)
is the stiffness tensor of the tube i.

The torsion curvature components are obtained after including the torsional twist, satisfying the equilibrium of the torsional moments, and using the equilibrium equation of the special Cosserat road model [START_REF] Antman | Nonlinear Problems of Elasticity[END_REF].

FORWARD KINEMATIC MODEL OF CONCENTRIC TUBE ROBOT

47

The development of these conditions yields:

ρ ,z = - K ,z K ,z ρ ,z dρ ,z ds = K ,xy K ,z ( ρ,y ρ ,x -ρ,x ρ ,y ) (3.24) with K ,xy = K ,x = K ,y .
The twist rate is expressed by:

γ = ρ ,z -ρ ,z (3.25)
The resulting equations (3.22), (3.24), and (3.25) form a nonlinear differential equation system with (γ, ρ ,z ) as integration variables.

At the robot base s = , γ is known as it is calculated from (3.21), however, ρ ,z is unknown. At the distal extremity of the robot s = L, the torsion curvature is zero, which means ρ ,z is zero. The system is then a nonlinear differential equation with boundary conditions. The bvpc matlab function allows to solve this kind of systems, which integrates from the initial limit to the final one and then integrates back and repeats this go-and-back integration until both initial and final conditions are satisfied. Once this system is integrated, the obtained curvatures are used to define the robot shape. The curvatures are then integrated analogously to the twist velocity integration by preserving the group structure on SE(). A variety of methods are presented in [START_REF] Park | Geometric Integration on Euclidean Group With Application to Articulated Multibody Systems[END_REF].

A less complicated alternative method is to approximate the robot shape by a set of small constant curvature sections of length ∆ s, then it is integrated using the Rodrigues' formula.

Consider H as an homogeneous matrix expressing the mapping between two frames and V is the twist between them. Recall the following well-known relationship in kinematics of rigid body [START_REF] Murray | A Mathematical Introduction to Robotic Manipulation[END_REF]:

⌈V⌉ = H -1 Ḣ (3.26)
where the operator ⌈.⌉ denotes the mapping from R  to se() (se() is the Lie Algebra of SE() [START_REF] Park | Geometric Integration on Euclidean Group With Application to Articulated Multibody Systems[END_REF]).

Neglecting the elongation of the tubes, the reference twist is given by:

ξ (s) = [  ] ρ(s) T T (3.27)

CHAPTER 3. KINEMATIC MODELLING OF CONCENTRIC TUBE ROBOT

Analogous to (3.26), the following is obtained:

⌈ξ (s)⌉ = H -1 (s) Ḣ(s) (3.28)
As the curvatures ρ(s) have already been calculated after integrating the previous system, and considering these curvatures constant at each sampling section of length ∆ s, the robot shape can be obtained recursively by the Rodrigues' formula:

H(s) = H(s -∆ s)expm(⌈ξ (s)⌉∆ s) (3.29)
where expm is a matlab function calculating the exponential of a matrix.

Differential kinematic model of concentric tube robot

In this section, the calculations of the differential kinematic model ending model, straight torsion model, and torsionally compliant model are presented in details. For a constant curvature section j of the robot, the spacial velocity vector V j ∈ R  of its distal extremity with respect to its base expressed in its frame base can be obtained by integrating the following equation:

⌈V j ⌉ = ḢH -1 (3.30)

DIFFERENTIAL KINEMATIC MODEL OF CONCENTRIC TUBE ROBOT

49

Substituting the homogeneous matrix H by the general expression (3.12) in equation (3.30) yields:

                                                       V j =             ρ ′ sin(ϕ j -ϕ j- )(cos(ρl) -) / ρ  ρ ′ cos(ϕ j -ϕ j- )(cos(ρl) -) / ρ  l ′ -ρ ′ (sin(ρl) -ρl) / ρ  cos(ϕ j -ϕ j- )(ρ ′ l + l ′ ρ) sin(ϕ j -ϕ j- )(ρ ′ l + l ′ ρ) ϕ ′ j -ϕ ′ j-             q = J j q if ρ = 0 V j =               l ′   ϕ ′ j -ϕ ′ j-             if ρ =  (3.31)
with: q = [α 1 ...α i ι 1 ...ι i ] T is the joints vector, ρ ′ = ∂ ρ ∂ q ,l ′ = ∂ l ∂ q , and ϕ ′ = ∂ ϕ ∂ q . (l ′ , ρ ′ , ϕ ′ ) are vectors of derivatives with respect to the joint variables q.

Once the velocity of each section is obtained, the general spacial velocity of the robot tip expressed in the base frame is given by:

V = J q (3.32) with J = J 1 + Ad( 0 H 1 )J 2 + ... + J j Ad( 0 H j ) (3.33)
being the spacial Jacobian. Ad(H) is the adjoint transformation of the matrix H [START_REF] Murray | A Mathematical Introduction to Robotic Manipulation[END_REF].

• In the case of bending model, the vectors (ρ ′ , l ′ , ϕ ′ ) can be obtained by deriving the expressions of (l, ρ, ϕ).

• In the case of straight torsion model, the derivatives of (ρ, ϕ) can not be obtained by a straightforward derivation.

Consider β * as the solution of the energy gradient set to zero (see (3.19)), it has:

       ∂ ρ j ∂ q = ∂ ρ j ∂ β * ∂ β * ∂ q ∂ ϕ j ∂ q = ∂ ϕ j ∂ β * ∂ β * ∂ q (3.34)
The derivatives (

∂ ρ j ∂ β * , ∂ ϕ j ∂ β *
) can be obtained by a straightforward derivation of the expressions of (ρ, ϕ) given by (3.10).

At β * , the gradient ∇ϒ = Φ(q, β * )=  and D q Φ = D  Φ + (D  Φ)(D q β * ) = , (D is the differential operator), which yields:

D q β * = -(D  Φ) - (D  Φ).
(3.35) Substituting (3.35) into (3.34) gives:

       ∂ ρ j ∂ q = ∂ ρ j ∂ β (- ∂ Φ ∂ β * ) - ∂ Φ ∂ q ∂ ϕ j ∂ q = ∂ ϕ j ∂ β (- ∂ Φ ∂ β * ) - ∂ Φ ∂ q (3.36)
• The Jacobian matrix for the case of torsionally compliant model has no closed form. As introduced in [START_REF] Rucker | Computing jacobians and compliance matrices for externally loaded continuum robots[END_REF]], a column J s i of the Jacobian corresponding to an element q i of the joint vector and at a distance s from the robot shape can be approximated by the following relationship:

J s i ≈ ( H i (s) -H(s) ∆ q i H(s) -1 ) ∨ (3.37)
with H i (s) being the pose of the end-effector at the joint position (q i + ∆ q i ) and the operator allows the mapping from se() to R  (see [START_REF] Antman | Nonlinear Problems of Elasticity[END_REF]). This approximation is valid only under the conditions of small variations ∆ q i .

Inverse kinematic model

The inverse kinematic model is required in inverse kinematic control strategies which are widely used in CTR, especially in teleoperated architectures in which the master target

INVERSE KINEMATIC MODEL

51

pose is interpreted to desired joints at the slave side (CTR) [Dupont et al., 2010b]. The inverse kinematic model is also used in design optimization in order to avoid obstacles [START_REF] Bedell | Design optimization of concentric tube robots based on task and anatomical constraints[END_REF]. For the bending model [START_REF] Sears | A Steerable Needle Technology Using Curved Concentric Tubes[END_REF], the inverse kinematic model introduced in [START_REF] Sears | Inverse Kinematics of Concentric Tube Steerable Needles[END_REF] is split into two mappings: the first mapping calculates the section parameters represented by the bending components (χ j , ν j ) using the generalized inverse of the Jacobian matrix linking the tip velocity V to the bending parameter velocities (χ j , ν j ). The second mapping is a closed-form expression of the joints (α i , ι i ) in terms of the bending components (χ j , ν j ).

For other models, the most oftenly used techniques for inverse kinematics calculation are not in a closed form. Root finding function has been used in [START_REF] Su | A MRI-guided concentric tube continuum robot with piezoelectric actuation: A feasibility study[END_REF] which allows to solve the inverse kinematic problem after implementation of the forward kinematic model.

In [START_REF] Su | A MRI-guided concentric tube continuum robot with piezoelectric actuation: A feasibility study[END_REF], a lookup table of joint and position samples is used to initialize the rootfinding algorithm. The Fourier series approximation function of the forward kinematic model mentioned in Section 2.5 is then inverted using a root-finding algorithm.

The pseudo-inverse of the Jacobian J + in [START_REF] Torres | Motion planning for concentric tube robots using mechanics-based models[END_REF] has been used to calculate the joint variations ∆ q approximatively corresponding to the tip position changes ∆ X in task space:

∆ q = J + ∆ X (3.38)
where (3.39) this pseudo inverse corresponds to the least square solution minimizing the joint velocities norm q and the errors Ẋ -J q .

J + = J T (JJ T ) -1 ,
The damped least square (see [START_REF] Wampler | Manipulator Inverse Kinematic Solutions Based on Vector Formulations and Damped Least-Squares Methods[END_REF]) of the Jacobian matrix J † has been largely used [START_REF] Hendrick | A multi-arm hand-held robotic system for transurethral laser prostate surgery[END_REF]] [Burgner et al., 2014[START_REF] Xu | Position control of concentrictube continuum robots using a modified Jacobian-based approach[END_REF][START_REF] Burgner | A bimanual teleoperated system for endonasal skull base surgery[END_REF] as in the case of serial robots where the inverse kinematic model has no closed form and for the same objective as mentioned above: (3.40) where:

∆ q = J † ∆ X,
J † = (J T J + λ I) -1 J T , (3.41)
where λ is the damping factor. This solution ensures a minimization of the errors Ẋ -J q 2 + λ q 2 . Inequality constraints can be added to this approximation problem in order to avoid the singularities and the joint limits.

Bifurcation phenomenon

α  () -α  () ( • ) α  (L) - α  (L) ( • )
Figure 3.5 -Effect of tube curvatures on the elastic stability: from the blue (***) to the red (***) dot-curve, the curvatures of the two tubes are increased as /, /, /, / mm - , and the lengths of both tubes are set to  mm. When a CTR is in motion and maneuvering, the interactions between tubes develop elastic energy dues to bending and torsion. The resulting shape of a set of interacting tubes at the equilibrium is conformed such that the elastic energy is minimal. In some configurations, the elastic energy allows two minimum solutions, which causes an instable state of the robot in the form of 'Jumps' or also named 'Bifurcations' [Dupont et al., 2010a], [START_REF] Webster | Design and Mechanics of Continuum Robot for Surgery[END_REF].

For the case of two tubes, these instable states appear once the relative tube rotations is close to π. The elastic stability can be evaluated by studying the evolution of the tubes angle difference at the distal extremity (α  (L) -α  (L)) according to the difference tube angles at the base (α  () -α  ()). The elastic stability depends on the tube lengths and curvatures

α  () -α  () ( • ) α  (L) - α  (L) ( • )
Figure 3.6 -Effect of tube lengths on the elastic stability: from the blue (***) to the red (***) dot-curve, the lengths of the two tubes are increased as (, , , ) mm, and the curvatures of both tubes are set to / mm - . as well. A simulation study has been carried out in order to show the effect of the lengths and curvatures of the tubes on the elastic stability. In Figure 3.5, the pre-curvature for both two tubes forming the CTR has been increased. It can be seen from the figure that the green dot-curve is at the limit of stability and the distal angle difference admits two solutions, which creates the jumps at the tip position for the red dot-curve. The amplitude of jumps increase once the curvature increases. The length has the same impact on the elastic stability as it can be seen from Figure 3.6. The characteristics of the two tubes used in this simulation case study are given in Table 3.1.

Conclusion

In this chapter, the forward, inverse and differential kinematic models for CTR have been introduced in detail. From this chapter, it can be concluded that the kinematic modelling of CTR represents a balance between complexity and accuracy which should be considered in the CTR motion control. The bifurcations presents a serious problem which should be

Introduction

In order to better benefit new less or non-invasive surgical procedures (Minimally Invasive Surgery, Single-Incision Laparoscopic Surgery, and Natural Orifice Transluminal Endoscopic Surgery), flexible medical devices are favorable [START_REF] Jessica | Continuum Robots for Medical Applications: A Survey[END_REF]. They are more convenient and safer in maneuver and navigation through body cavities than rigid devices.

Moreover, flexible devices offer more accessibility to the surgical site in case of complex path from the insertion body point to the surgical site. CTR is flexible robot with the main concept of inserting pre-shaped flexible Nitinol tubes such that each of them can rotate and translate [START_REF] Sears | A Steerable Needle Technology Using Curved Concentric Tubes[END_REF].

CTR has been used in many medical applications, such as neurosurgery [START_REF] Anor | Algorithms for design of continuum robots using the concentric tubes approach: A neurosurgical example[END_REF], intracardiac surgery [START_REF] Bedell | Design optimization of concentric tube robots based on task and anatomical constraints[END_REF], endonasal skull base surgery [START_REF] Burgner | A bimanual teleoperated system for endonasal skull base surgery[END_REF] and lung biopsy [START_REF] Luis G Torres | Task-oriented design of concentric tube robots using mechanics-based models[END_REF]. It is obvious that the required design of the robot tubes differs from one medical application to another. The tubes dimensions required for a given application depend mainly on the specific task and the involved anatomical structure. To achieve good operation performances, the robot tubes design should be optimized accordingly.

The evolution of CTR design optimization presented previously in a general context can be traced from three aspects. Firstly, the design optimization follows the advances in kinematics modelling. In the early works on tubes design optimization, only basic kinematic model was used [START_REF] Bedell | Design optimization of concentric tube robots based on task and anatomical constraints[END_REF][START_REF] Burgner | A bimanual teleoperated system for endonasal skull base surgery[END_REF]. Later on, more accurate kinematic model has been involved [Burgner et al., 2013a]. With more accurate kinematics, better robot design optimization could be obtained since the kinematic model is involved in the optimization calculation. Secondly, the design optimization evolution could also be analyzed from the view

INTRODUCTION

57

of clinical applications addressed seeing that many specific design optimization algorithms have been proposed for different clinical procedures [START_REF] Anor | Algorithms for design of continuum robots using the concentric tubes approach: A neurosurgical example[END_REF][START_REF] Bedell | Design optimization of concentric tube robots based on task and anatomical constraints[END_REF] [ [START_REF] Burgner | A bimanual teleoperated system for endonasal skull base surgery[END_REF][START_REF] Luis G Torres | Task-oriented design of concentric tube robots using mechanics-based models[END_REF]. Another axis of evolution is the transition from considering a discrete set of points in the cost function evaluating the coverage of the desired workspace to a volume-based cost function [Burgner et al., 2013a]. This evolution can be illustrated in Figure 4.1. The problem of elastic instability was not taken into account in CTR tube design optimization in literature until the work of [START_REF] Ha | Achieving elastic stability of concentric tube robots through optimization of tube precurvature[END_REF], but it is limited to a particular case of two planar tubes with variable pre-curvature of the tubes. [START_REF] Bergeles | Concentric tube robot design and optimization based on task and anatomical constraints[END_REF] have shown a way to analyze the stability when the number of tubes exceeds two, but it is not applicable for arbitrary number of tubes. The elastic instability may occur when the difference of tube angles at the base of the robot is at the neighborhood of π. The rotation angle difference of the tubes at the distal part has no unique solution in the case of instability, in which the solutions are analyzed in a direct way from graphics illustrations of the solutions distribution. The problem of pre-weighting the objective functions is then avoided.

Scalarization technique

Among all methods to solve Pareto optimal solutions for multi-criteria optimization, the scalarization method is the most commonly used one. In scalarisation method, the objective functions to be minimized are transformed to a single-objective scalar function by weighting and summing them [START_REF] Kasprzak | An approach to facilitate decision tradeoffs in Pareto solution sets[END_REF]. The minimization of the new objective function is defined as:

m x in n i= ω i • Φ i (x), n i= ω i = , (ω i >) subject to x ∈ ϖ (4.3)
The minimum of this new optimization objective belongs to the Pareto curve. This method is easy to understand and to implement. However, as the weights are pre-selected with no exact knowledge of the importance of each objective function, they may not correspond to the relative importance of the objective functions. It is also difficult for the decision maker

Grid-searching technique

As an alternative to the scalarization technique, grid-searching method generates a sampling of the Pareto curve and allows to avoid the drawbacks of scalarization previously introduced due to the pre-selection of the weights. Three steps are to be followed in order to define the best trade-off between different objective functions:

• First, the Pareto set is generated by taking a sampling grid of the decision space x and then to calculate their corresponding objective functions. The obtained points are plotted in the objective space as illustrated in Figure 4.4 for a two-function case, where the Pareto curve can be approximated by a polynomial. If the grid resolution is fine enough, this method guarantees to find both convex and non-convex portions of the Pareto curve.

• The next issue after tracing the Pareto set is to find the suitable solution on it. Many ap-TUMOR SURGERY proaches have been proposed in literature to solve this problem. The approach used in this work is the no-preference approach assuming "equal importance" of each objective function where the mostly used way to find the best solution is to calculate the L  norm for all the points as follows:

D L  (x) = n i= (Φ i (x) -Φ min i )  , (4.4)
where Φ min i is the minimum of Φ i , and then to select the point corresponding to the minimum of this distance [START_REF] Kasprzak | An approach to facilitate decision tradeoffs in Pareto solution sets[END_REF]. In case of maximization of the objective functions, the maximum of the L  norm is used. Note that if the equal importance of the objective functions is not considered, then one of the plot axes should be scaled before calculating the L  norm, and this scale corresponds to the weight. This is different from scalarization methods where the weights are pre-selected without analysis of the solutions distribution in the objective space and this grid-searching technique where the scaling of the objective functions (weighting) is performed after analysis of the solutions distribution.

• The last step is to map from the Pareto set back to the decision space. There are many methods in literature which ensure this mapping such as intersection of constant objective contours, etc. In this work we select the decision variables that correspond to the optimal objective solution with highest L  norm obtained in the previous step. In the previously presented technique, the correspondence between the weights and the solutions is not known i.e. the weights se-

TARGETED SURGICAL APPLICATION

63

lected by the decision maker do not necessarily correspond to the relative importance of the objective functions. However, in this method the scaling in the objective space is done in a straightforward way with aid of graphics and corresponds to the relative importance of the objective functions. Moreover, in scalarization method, one weights vector corresponds to only one point of the Pareto curve which requires several optimization times to obtain a portion from the Pareto curve imposing computational burden, while in grid-searching method this problem is avoided and both convex and non-convex Pareto curve can be found (see [START_REF] Caramia | Multi-objective management in freight logistics[END_REF] for more details).

Targeted surgical application

Frontal lobe to accomplish the resection procedure cf. Figure 4.5.b-d. Following this procedure, the brain will swell up after the operation, which creates a brain pressure and causes a large trauma for the patient. To alleviate this problem, an endonasal approach is proposed by using CTR which is flexible and adaptable with the anatomical complexity. Note that, similar as the case in [START_REF] Burgner | A bimanual teleoperated system for endonasal skull base surgery[END_REF] which addresses pituitary gland tumor removal, this operation requires bone drilling in order to access to the surgical site and reach the tumor. But the insertion path in our case from the nostril to the anterior wall of the frontal sinus (location of the drilling point) is surrounded with structures of high sensitivity cf. to ≈ mm for bigger ellipse diameter) and a length of ≈ mm for the patient sample studied.

These allowed workspace borders are illustrated by six green ellipses in Figure 4.6. These anatomical constraints (orbits, optic canal, nerves etc) should be considered in an eligibility criterion. In order to reduce the friction during the robot insertion, a large passive tube is pre-inserted before introducing the robot. This facilitates the insertion of the robot in the first part of its path which is characterized by high friction according to the surgeons. Once the navigation path is traversed, the distal section of the robot accesses to the anterior skull base to manipulate at the tumor location.

A scenario of a sphere tumor defined by its center position coordinates c and its radius a =  mm is shown by purple color in Figure 4.6. The frequent scenario in reality is that the tumor is located at the left of at the right of the frontal lobe. To perform the task of tumor removal, various tools can be attached to the distal extremity of the robot for purposes like ablation, retraction, blood stopping, etc.

CTR configuration selection

By considering the aforementioned anatomical description and operation requirements, a robot configuration is proposed in this section.

• Kinematics:

To fulfill accessibility and dexterity requirement, the candidate robot to be optimized is a three-tubes CTR. Only one outer tube can ensure the navigation from the nostril access to the entry of the desired surgical workspace because this path is almost a straight path after preparation from the surgeon; remove some nonsensitive structures (mycosis, etc) to facilitate the robot insertion. A variable curvature section is required to manipulate the robot tip as it offers more manipulability and coverage of the desired workspace. Two pre-curved inners tubes are required to form the variable curvature section. The outer-most tube T  is straight and the other two tubes T  and T  have two parts (straight and curved) as shown in Figure 4.7. The tubes deployment decides the number of tubes existing in a section, the number of sections, their order, lengths and curvatures (constant or variable curvature).

In terms of reachability, the proposed configuration allows to benefit from both advantages of constant curvature and variable curvature of the robot distal section. The distal section can The translation limits, the curved and straight tube lengths define this transition constant/variable curvature of the distal section. The distal section manipulates the tool to accomplish the requested task (grasping, aspiration, electrocoagulation etc). For the proposed robot structure, the constant curvature section allows to reach the regions far from the centerline of the robot shaft and the regions close to it could be reached by a variable curvature distal section as its curvature can decrease and then comes closer to the centerline.

• Actuation:

Concerning the actuation block it is composed of three blocks with one for each tube.

Note that for the proposed actuation unit and for each tube, the rotation and translation are decoupled. Thanks to this decoupling, each rotation is controlled by one dedicated motor, which improves the control performance, facilitates the calibration procedure and reduces the modelling complexity. The kinematic model used in this work for forward kinematics calculation is the torsionally compliant model [START_REF] Dupont | Torsional kinematic model for concentric tube robots[END_REF] presented previously. In this model, the pose of a point G(s) from the robot shaft and at a distance s from the robot base is obtained by two mappings h  and h  . The first mapping h  is to calculate the curvature u(s) in terms of the joint variables q as noted by u(s) = h 1 (q). The second mapping is to use the curvatures u(s) in order to define the pose G(s) = h 2 (u(s)).

The forward kinematic modelling is then split into two steps: the first one is to calculate the curvature u(s) in term of the joint variables q as noted by u(s) = h 1 (q); the second one is to use the curvatures u(s) in order to define the robot shaft G(s) = h 2 (u(s)).

Formulation of grid-searching optimization algorithm for CTR

Recall that previous CTR design multi-criteria optimization algorithms use scalarization methods, which lumps the optimization objective functions in one single-objective function by assigning weights empirically to each objective function. The challenge of these methods is to properly choose the criteria weights which defines well their priorities and consequently

k  =  : ∆ k : mm - ,k  =  : ∆ k : mm -
influences the optimization results [START_REF] Unal | Multi-criteria Design Optimization of Parallel Robots[END_REF]. However, there is no systematic and rigorous way reported to choose the weight values. Comparatively, the proposed grid-searching method [START_REF] Caramia | Multi-objective management in freight logistics[END_REF] allows an interpretation of the optimization results in transparency, and the decision is taken by studying a set of possible solutions. The objective functions are evaluated along a large range of each optimization parameter.

Decision space of the CTR optimization problem

For our specific surgical application, the manipulation part of the robot extends from the drilling point to the posterior part of the tumor. For such a distance, the lengths of the curved parts of the two inner tubes do not have significant effect on the optimization results. Moreover, reducing the number of decision parameters allows to decrease the computation time significantly. Therefor the lengths of tubes are pre-selected suitable for the task respecting the workspace and anatomical constraints. The curvatures k  and k  of the tubes T  and T  respectively are considered as the decision variables. The localization of their ranges is firstly estimated from the geometry of the anatomical constraints, then obtained by executing the Algorithm.1 as introduce later.

FORMULATION OF GRID-SEARCHING OPTIMIZATION ALGORITHM FOR CTR 71

Objective space of the CTR optimization problem

During the operation, the displacement of the robot tip inside the whole desired workspace is important. This displacement should be as long as possible and at the same time should be safe without sudden jumps. The objective space is then defined by two criteria (objective functions):

• Reachability criterion: In order to evaluate the reachability of the surgical workspace (tumor), this desired workspace is sampled homogeneously to a set of points W(i)

(i = , ..., n,
where n is the total number of sample points). If the Euclidean distance between the workspace center and the tip position is less than or equal to the work-space radius a, we check which point from W(i) is closest to this tip position and then increase the score of reachability of this point W(i). For a given configuration (k  , k  ), a point from W(i), i = , ..., n can be reached more than once. Considering m is the number of the reached points from W(i),

and recalling n is the number of the work-space sample points, the reachability is evaluated as:

C reach (k  , k  ) = m n ×  (4.5)
• Stability criterion: For the elastic stability analysis, previous works rely on the evolution of the difference of distal tubes angles in terms of the base difference tubes angles [Dupont et al., 2010a]. The case of three tubes was treated by [START_REF] Bergeles | Concentric tube robot design and optimization based on task and anatomical constraints[END_REF]. There is no systematic and general approach so far which allows to evaluate the elastic stability once the tube number exceed three.

In this work, a new way of stability evaluation is proposed. It relies on numerical evaluation of the tip trajectory smoothness. In case of elastic instability, jumps appear in the robot tip trajectory even if the tubes rotations variations are continuous and smooth. The elastic stability could be quantified by evaluating the Euclidean distance √ χ   + χ   covered by the robot tip, since the robot shaft rotation is around the Z axis of the robot frame. In order to find out the jumps (case of elastic instability) on the distal angle difference in ascending and descending (see Figure 4.9), the bvpinit Matlab function is used to detect the symmetric solution (where the system admits two solutions) which shows a distal angle jump leading to a tip position jump (see blue and red arrows if Figure 4.9). The newly proposed approach is to quantify these jumps or discontinuities of the robot tip trajectory. Firstly, the numerical derivative of the Euclidean distance in term of the joint variation (rotation sampling step) is

∆ √ χ   + χ   △q ,
then the standard deviations of these derivatives are calculated. This mean allows to quantify the distance between two successive points of the trajectory along the whole set of points forming the trajectory and thus smoothness. The stability of a given configuration

(k  , k  )
is inversely proportional to the derivative standard deviation as the former increases with the occurrence of the jumps.

C stab (k  , k  ) = normalization(/std( ∆ √ χ   + χ   △q )) (4.6)

Constraints

The first constraint is to consider only the decision variables satisfying (k  <k  ) during the grid sampling. This heuristic condition is recommended in [Dupont et al., 2010a] ellipse. This condition is to be checked before evaluating the reachability and stability. The configurations outside of the borders are rejected and will not be evaluated in the reachability and stability criteria.

Grid-searching algorithm

Based on the selection criteria and the constraints aforementioned, a grid-searching algorithm can be developed for the CTR design optimization problem.The algorithm is explained as in following. As illustrated in Algorithm 1, for all the configurations of the grid (decision space), the forward kinematic model is calculated for all possible joint configurations (resolution sampling in rotation of π/ rad, and in translation of  mm is chosen to allow the coverage evaluation of the desired workspace regarding its dimension. The work of [START_REF] Burgner | A bimanual teleoperated system for endonasal skull base surgery[END_REF] helps for this resolution selection). At each joint configuration, the limits defined by the ellipses are checked, then if the robot shaft satisfies these borders, the tip position is stored, else, this position is rejected. At the end of the algorithm, the stored positions of each configuration from the grid are evaluated in reachability and stability.

After the first execution of the Algorithm 1 with large sampling interval (mm - ) which allows to locate the variation ranges of the decision variables (k  , k  ); where the optimal decision variables intervals corresponding to the highest reachability and stability may lie in. A second execution is conducted with smaller sampling interval ( mm - ) to ensure a finer mapping from the decision variables (k  , k  ) to the objective space (C stab ,C reach ). The resolution of the grid in both large and small sampling intervals is estimated from this range of variation of the decision variables (k  , k  ) as given in The grid-searching optimization method as explained in last section was executed twice: the first time with larger parameters interval (k  , k  ) and curvatures sampling resolution; then the second time to refine the optimization procedure with smaller parameters interval and curvature sampling. In order to accelerate the computation, the algorithm has been run on a server with the characteristics: Intel(R) Xeon(R), CPU E5-2695-v3 2.30GHz, (27 cores / 396 Go ram).

The results of the new stability evaluation approach obtained through the second execution are interpolated and shown in Figure 4.10. It can be seen from this illustration that the stability criterion index increases along both X and Y axis, i.e. when the curvatures (k  , k  ) decrease (curvature radii increase) which is consistent with the evolution stability in terms of curvatures presented in [Dupont et al., 2010a]. The variation amplitude of the stability criterion is around

..
Theoretically, the reachability can be improved by increasing the curvatures. However, in our study this rule is not well satisfied as the anatomical obstacles are considered along the whole robot shaft, i.e. many points from the workspace are eliminated by the eligibility condition. From Figure 4.11 it can be seen that simple linear relationship does not exist between reachability and curvature. In some cases the reachability increases by increasing the curvature and it decreases in other cases by increasing the curvature, which is explained by the anatomical constraints, i.e. the elimination of the tip positions which do not satisfy the eligibility condition and excluded from the reachable set of points. It can also be seen that the variation amplitude of the reachability is around % which justifies the importance of tubes design optimization.

From the mapping to the objective space shown in Figure 4.12, the solution with highest TUMOR SURGERY Remark 3.1: The used optimization method addresses two-objective problem, however it can be extended to more than two objectives. If the number of objective functions exceeds three, the description of the Pareto set is hyperplane and more complicated compared to the case of two objectives.

/k  (mm) /k  (mm) C stab
P opt P min x (k  =   , k  =   ) DL  =  ,     C stab C reach
Remark 3.2: Implementing of all optimization algorithms used for CTR existing in literature to carry out a comparison study is not a trivial work, especially for multi-objective scalarization methods in which one vector of weights provides only one point of the Pareto curve. This requires tremendous execution times for each algorithm with different combinations of the weights vector which is a huge task. The kinematics calculation is also time consuming, and to obtain the inverse kinematic model is an optimization problem itself and hence computationally heavy. For this reason, we focus on analyzing the advantages and the drawbacks of the technique used in literature and our proposed technique in this work from the working mechanism point of view and only numeric analyses of the latter are presented.

Conclusion

The main contribution of the work presented in this chapter is to provide an alternative method which avoids using scalarization techniques for multi-objective CTR tubes design optimization. These scalarization techniques as has been explained require weights pre-selection without any a priori information about the objective functions distribution. On contrary, gridsearching method has been used with a direct interpretation of the optimization results. Moreover, a new method of elastic stability evaluation has been used as an objective function. In this new approach, the stability is evaluated explicitly based on the smoothness of the robot tip trajectory instead of evaluating the tubes angles difference as in literature.

The obtained optimized tubes parameters will be used to curve the tubes to be integrated to our CTR platform, and the addressed deep tumor removal surgical procedure will be tested through ex vivo lab test and in vivo cadaver test. Task-space Position Control of CTR

Using Approximated Jacobian Preamble

The kinematics modelling of CTR has been presented in Chapter 3. It should be noted that so far the developed kinematic models are either simplified and hence inaccurate, or complex and hence computationally heavy. Existing works in control of CTR are based on inverse kinematics calculation and thus the control performance largely relies on the accuracy of the kinematic model used. In this chapter, a new control method designed at the actuator level is proposed which shows that the control design of actuator input in task-space with approximate Jacobian matrix provides more flexibility and robustness in handling inaccuracy in kinematic model and maintains good control performance at the same time [START_REF] Boushaki | Task-space position control of concentric-tube robot with inaccurate kinematics using approximate jacobian[END_REF].

Introduction

Mechanics-based models with different assumptions have been developed to realize real time applications of CTR control. The first model proposed assumes that the outer most tube has infinity stiffness compared to inner tubes [START_REF] Furusho | Curved multi-tube systems for fetal blood sampling and treatments of organs like brain and breast[END_REF]. Later on, more realistic physical phenomena previously neglected, such as bending and torsion, have been taken

into consideration in order to improve the kinematics [D. Caleb and R. J, 2009b]. Based on the developed kinematic models, several motion planning algorithms have been developed [START_REF] Lyons | Motion planning for active cannulas[END_REF][START_REF] Lyons | Planning active cannula configurations through tubular anatomy[END_REF]. These motion planning algorithms resort to inverse kinematics calculation to obtain the actuator joint angles to achieve desired task space positions. However, the inverse kinematics calculation of CTR is not straightforward due to the nonlinear mapping between relative tube displacements and tip configuration as well as due to the multiplicity of solutions. A geometric approach is given for single and multiple sections in [START_REF] Neppalli | Closed-form inverse kinematics for continuum manipulators[END_REF] by applying an analytical process to solve inverse kinematics based on modelling of each section with a spherical joint and a straight rigid link. Jacobianbased methods represent another approach to robot-independent inverse kinematics. In the Jacobian-based inverse kinematics strategy, it is possible to build actuator limits into the control law so that the robot trajectory is always physically realizable [R. J. Webster and Jones, 2010], e.g. [START_REF] Sears | Inverse Kinematics of Concentric Tube Steerable Needles[END_REF] for active cannula results. Root finding function are also used to solve the inverse kinematics calculation [START_REF] Su | A MRI-guided concentric tube continuum robot with piezoelectric actuation: A feasibility study[END_REF].

It is noticed that existing works in literature on CTR position control are based on inverse kinematics calculation to get the necessary actuator joint angle first in order to achieve control tasks that are usually defined in operational space (image space, Cartesian space, etc).

Consequently, the control performances are sensitive to the kinematic modelling inaccuracy.

However, so far there exists no closed form inverse kinematic of CTR when the number of tubes exceeds two, and the modelling challenge lies in the physical phenomena due to elasticity interaction between tubes such as bending, friction and torsion. With the existence of kinematics modelling errors, the inverse kinematics based control methods may suffer position errors and lead to safety problems. Although in teleoperation scenario this kind of positioning deviation could be partially compensated by the surgeon, it is not always guaranteed to work as the surgeon intuitive manual correction command could be mis-interpreted by the inaccurate kinematics used in the inverse kinematics calculation.

FORMULATION OF CTR POSITION CONTROL TASK AT ACTUATOR LEVEL 81

To alleviate this problem, the CTR position control from the actuator input level is investigated in this work. It is shown that even if only approximate Jacobian matrix is available, the position control task can still be accomplished with a simple PID form control torque input based on task space feedback position errors. The proposed control design is rigorously analysed through theoretical proof. Simulation studies based on a realistic experimental platform (see Figure 5.2) are shown to justify the advantage of proposed task-space control method over inverse kinematics based design method in handling kinematics modelling inaccuracy.

Formulation of CTR position control task at actuator level

In literature, the position control of concentric tube robot is mostly formulated as a motion planning problem based on inverse kinematics calculation. Given the desired position or trajectory x d in operational space, which can be obtained by interpreting the tele-operator hands motions in case of tele-operation scenario (refer to Figure 5.1), the corresponding desired joint positions or trajectories q d are calculated by solving the inverse kinematics:

q d = f InvKin (x d ), (5.1)
and it is assumed that the actuator for each joint will faithfully generate the desired joint motion q d . However, if the robot kinematics is not accurately modeled, the joint position obtained through calculation (5.1) may not guarantee the desired operational space position.

x d q d

Figure 5

.1 -The task-space desired trajectory at the master side is converted through the inverse kinematic model f InvKin to a desired joint trajectory at the slave side.

CHAPTER 5. TASK-SPACE POSITION CONTROL OF CTR USING APPROXIMATED JACOBIAN

In this work, the kinematics modelling error from the actuator control input level are compensated and the inverse kinematics calculation is bypassed, inspired by the works on uncertain kinematics task-space control for serial robot in literature [Liu and Cheah, 2006] [ Liu et al., 2006].

For existing CTR design, the joints are normally actuated by current driven motors which possess the following dynamics:

M i θ +C i θ = K i u i - τ ei r i (5.2)
where θ denotes the motor rotor shaft angle, M i denotes the rotor inertia moment, C i denotes the friction coefficient and K i the torque constant. u i is the control current input, τ ei is the external load torque. r i = θ i /q i is the transmission gear ratio from the motor angle θ to the CTR joint variables q, which is usually quite big and thus the effect of external load τ ei could be neglected. Considering the linear relationship between actuator and joint variables θ i = r i q i and, for simplicity of expression and without loss of generality, setting K i r i = , the dynamic relationship between CTR joint variables and actuator input can be written from (5.2) as:

M q + C q = u (5.3) where M = diag([M  , • • • , M n ]), C = diag([C  , • • • ,C n ]) and u= diag([u  , • • • , u n ]).
Now the position control task could be specified at the actuator level: design the actuator current input u such that the tip position x of CTR can reach the desired operational space position (∆ x = xx d = ) even with inaccurate kinematics model. Compared to joint level control which is directly influenced by the kinematics errors, the lower level control design provides another solution which possesses more flexibility to handle such kind of problem.

Task-space position control of CTR with inaccurate kinematics

The operational space velocity ẋ and the CTR joint velocity q could be determined by Jacobian matrix J(q) as [R. J. Webster and Jones, 2010] ẋ = J(q) q.

(5.4)

TASK-SPACE POSITION CONTROL OF CTR WITH INACCURATE KINEMATICS 83

With the presence of kinematics errors, accurate Jacobian matrix is unknown and only estimated (best-guess) Jacobian Ĵ(q) is available.

In this study, it is assumed that the Jacobian estimation error Ĵ(q) -J(q) is bounded with an upper limit β as Ĵ(q) -J(q) ≤ β .

(5.5)

Based on the approximate Jacobian matrix Ĵ(q), the control current input can be designed with the use of task-space position error as in the following form:

u = -ĴT (q)K p ∆ x -K d q -k i t  ( q + α ĴT (q)∆ x)d t (5.6)
where K p = k p I, K d = k d I and k i are positive definite control gains, α is a positive scalar.

Since the CTR joint velocities are not directly measurable and noting the fact that q = r -1 θ , the control input is actually implemented using easily available and more accurate actuator joint velocity θ as

u = -ĴT (q)K p ∆ x -K ′ d θ -K ′ i t  ( θ + α ′ ĴT (q)∆ x)d t (5.7)
which is equivalent to the control design in (5.6) with K

′ d = K d r -1 , K ′ i = k i r -1 , α ′ = rα.
It is seen from the proposed control input (5.6) or (5.7) that no inverse Jacobian calculation is required in this control method which is in general much more complicated than transpose Jacobian calculation. The task-space position error can be obtained online through different imaging modalities (camera, ultrasound, MRI, magnetic tracker, etc) according to different clinical applications. Next, the system positioning performance through Lyapunov analysis are investigated. Substituting control input (5.6) into the dynamics equation (5.3), the closed-loop system dynamics has

M q = -ĴT (q)K p ∆ x -(C + K d ) q -k i t  ( q + α ĴT (q)∆ x)d t .
(5.8)

Propose a Lyapunov function candidate as

V =   qT M q + α qT M ĴT (q)∆ x +   ∆ x T K p ∆ x +   z T k i z (5.9)
where z = t  ( q + α ĴT (q)∆ x)d t . control design. On the other hand, it should be noted that the Lyapunov analysis is very conservative and so are the conditions deduced meaning that in many cases even if conditions

(5.18) and (5.19) are not satisfied the control task could be still achieved.

Remark 5.2: Although the proposed actuator current control input also possesses a PID form, it is totally different from the PID control of the actuator which is used to generate the calculated joint position as in inverse kinematics based methods. The difference lies in that the proposed PID actuator controller takes no desired actuator joint position through inverse kinematics as reference but the direct operational space positioning error as feedback signal, and this makes the difference between joint-space controller and task-space controller. In this simulation, two kinematic models ares used: the Transmissional Torsion Model (TTM) as the true kinematic model for forward kinematics calculation to feedback the position, and the simpler and less accurate Torsion-Free Model (TFM) which is assumed to be the model that is known to the designer and used to carry out the inverse kinematics based control. These two models with different accuracy are used to simulate the case of inaccurate modelling. Based on TFM the Jacobian matrix is calculated which represents the inaccurate Jacobian used in the proposed actuator current input design and in the inverse kinematic control. Here the more accurate but more complicated torsionally compliant model is not used as the true kinematics model since the purpose of this simulation study is to investigate the capability of the proposed control method in handling kinematics inaccuracy, in this sense the choice of true model does not make real difference as both are simplified approximations of the reality.

Simulation studies

CHAPTER 5. TASK-SPACE POSITION CONTROL OF CTR USING APPROXIMATED JACOBIAN

Simulation studies are carried out to compare the performances of the traditional inverse kinematics based control method with the proposed task-space controller in face of different kinematics uncertainties. The Tubes parameters used in this comparison simulation between task-space control and inverse kinematics control are given in Table 5.1. For the task-space (TS) control proposed, the torques are calculated using the transpose of the Jacobian matrix, the Cartesian errors, and the joint velocities. The Jacobian is calculated based on the TFM model and the Cartesian position feed-back are calculated by the TTM model after integration of the actuators dynamics. The actuator dynamics parameters used are the same as calibrated from the real prototype.

- + X d Ĵ+ (q) TFM ∆ X q i = q i-1 + ∆ q X = TTM(q) ∆ q X a) - + X d τ(∆ X, ĴT (q)) TFM ∆ X M q + C q = τ X = TTM(q) τ X b)
For the position control task, the desired trajectory is defined by 40 points in cannula base frame (as seen in Figure 5.2) along a circular trajectory of 10 mm radius which is within the workspace and defined by:

x d = cosφ , y d =  √ sinφ -, z d = sinφ +  √ .
(5.20)

A position sensor is assumed to be available to provide online update of task-space positioning errors (in simulation, this information is provided by the TTM model) and it is also assumed that this sensor information is used in the inverse kinematic calculation at the same sampling frequency.

Two types of simulations have been carried out to compare the position control performance of the two control strategies used (IK control and TS control):

-In the first type of simulation study, the same kinematic model (TFM) is used for both control methods and no kinematics parameter inaccuracy is considered.

-In the second type of simulation study, except using the same simplified model (TFM), kinematics parameter inaccuracy is also introduced to show the control performance of both methods against both structural and parametric inaccuracies in kinematics.

By using the control algorithm proposed as (5.6), the control performances of the taskspace approximate Jacobian position control method with the same kinematic uncertainties as for inverse kinematics based control method are illustrated in Figure 5. analysis presented before: the proposed control method can only tolerate uncertainty (β ) to certain level as specified by conditions (5.18) and (5.19) which is natural as no control method can tolerate infinitely large uncertainties. As the theoretical based on Lyapunov method which is conservative, the exact value (threshold) of kinematic inaccuracy can not be explicitly calculated. Remark 5.3: Although in this simulation study the inverse kinematics based control performance using TFM model is shown for comparison, it does not mean that proposed method is compared to this specific model. In fact, through the simulation study we want to show the fact that with the existence of kinematics modelling error, which is the case for all existing kinematics models, the proposed method can provide automatic compensation at the actuation level to achieve the control task. As mentioned in the introduction, although with manual compensation of the surgeon existing models may be enough for current teleoperated control task, it is still desirable if this manual compensation work could be left to the lower level motion controller such that the surgeon can operate in more intuitive way. To achieve that, the control design should provide the robotic instrument capabilities to tolerate uncertainties and disturbances within certain extent. Remark 5.4: The proposed control method relies on Jacobian matrix of the concentrictube robot. It is noticed that for most existing kinematics models the Jacobian matrix, i.e. the mapping from the operational space velocity to joint velocity, can be obtained. For very complicated models (existing or to come), an optimized approximation could be used for the Jacobian matrix. In any case, as shown in this work, the proposed control method may work with only an approximated Jacobian without need of exact kinematic information.

Conclusion

Existing motion control methods for concentric-tube robot are mainly based on inverse kinematics calculation which is prone to inaccuracy in kinematics modelling. Moreover, the use of the most accurate kinematic model is complex and time consuming. Through this work, 5.5. CONCLUSION 95 we show that this problem can be tackled by using the approximate Jacobian matrix and direct measurement of operational space error in the actuator control input design. The positioning errors are guaranteed to converge asymptotically through rigorous Lyapunov analysis given that the control parameters are properly chosen to satisfy certain conditions. Therefore, accurate CTR motion control can be realized without resorting to highly complicated kinematic model but a simpler alternative, which represents the biggest advantage of this method. This result presents a new control method that has not been explored for concentric-tube robot motion control. It is verified through simulation studies that the control performance of the proposed method with the presence of kinematics inaccuracy is better than inverse kinematics based method.

CHAPTER

6

Conclusions and Perspectives

Conclusions

In this thesis, the principle of CTR has been explained and existing research CTR works in literature have been presented. Mainly the CTR research works can be categorized into:

the design of actuation unit and tubes; the kinematic modelling; the path planning and the control. The works on tubes design and path planing are usually associated to a given surgical application. The tube design works optimize the dimensions of the tubes according to the surgical work space and the anatomical constraints. Feasibility of using CTR in various surgical applications (lung biopsy, neuro-surgery, etc) has been studied and validated through in vivo tests for some prototypes.

The forward, inverse, and differential kinematic models have been presented with details.

The evolution of kinematic modelling from assuming infinite stiffness of the outer tube to introducing the bending and torsion tubes interaction in the model and the cause and effect of 'bifurcations' have been presented. The second contribution is about compensation of kinematic modelling error from the actuator control input level and avoiding the inverse kinematics calculation for traditional CTR motion control. This contribution is motivated by the fact that the kinematic modelling of CTR is challenging because of complicated physical phenomena caused by the elasticity interaction between tubes. A new control method from the actuator level has been proposed

and has shown that the control design of actuator input in task-space with approximate Jacobian matrix provides more flexibility and robustness in handling inaccuracy in kinematic model since no accurate Jacobian matrix is obligatory. It has been shown through simulation study that the proposed control method presents better performance compared with traditional inverse kinematics based control method used in literature in face of kinematics inaccuracy.

Perspectives

As a new area of research in medical robotics, many issues in design, modelling, sensing and control are still open. For all existing prototypes, the actuation unit is not optimized to meet the requirements of the considered surgical application, e.g. limited space inside MRI tunnel. Compact robot design in weight, length, form, and assembly/disassembly in the operating room is required which can be one of the future directions in CTR design.

Regarding optimization tf the tubes design, other parameters such as the outer and the inner diameter can be added to the decision space. Moreover, the dexterity and force generation of the end-effector are also important criteria to evaluate the performance capability. They can be introduced in the objective function of design optimization algorithms. As for material of the robot tube, rather than Nitinol, other materials may be used such as thermoplastic materials.

The work of [START_REF] Amanov | Additive manufacturing of patient-specific tubular continuum manipulators[END_REF] is a start of this new research direction in CTR design.

Other materials could offer full elastically stable workspace. Sensing system, visual systems are usually used to feedback the robot shape (eg. in [Croom et al., 2010]), however, during the navigation the robot is not visible. The interventional MRI feedback or even 3D Ultrasound 6.2. PERSPECTIVES 99 feedback [START_REF] Ren | Tubular enhanced geodesic active contours for continuum robot detection using 3d ultrasound[END_REF]] are known to be slow for real-time applications. For this reason, other shape detection sensors embedded in the robot could be used while considering the complexity of integration inside the tubes of small diameter. For CTR motion control, the work presented in [START_REF] Bibliography | Real-time adaptive kinematic model estimation of concentric tube robots[END_REF] shows the promising use of adaptive laws in future either at the kinematic stage or at the control level in real-time for enhancing position control performances of CTR.

In future, the following works are considered to be further pursued :

-Lab experimental tests are scheduled as a direct continuity of this work. The tubes are to be curved accordingly to the obtained optimal curvatures in order to validate the concept of the addressed surgical procedures. For the second contribution of this thesis, the kinematic modelling has already been implemented in the real platform. The controller proposed is to be implemented and compared to the inverse kinematic control.

-In the preliminary tests, it is noticed that the effect of frictions existing in the system is quite significant. Future works will investigate how to deal with system frictions and improve operation accuracy and safety.

-An extension of the optimization algorithm used in this work will be proposed by including the outer and inner diameters of the tubes in the decision space and by considering dexterity in the objective space.

-Adaptive control for CTR is to be studied which is another strategy to compensate the kinematic inaccuracies. The unknown dynamics of the system can be estimated and adapted in real-time. The obtained results will be compared to the second contribution presented in this thesis.

  contrôle assure une robustesse en présence des incertitudes cinématiques au niveau de la matrice Jacobienne et permet d'obtenir des bonnes performances de contrôle en terme d'erreur de poursuite. assisted surgery and flexible robots

1. 1

 1 .a and Figure 1.1.b). The tubes are actuated in axial rotation and linear displacement. Actuating the tubes creates an elastic interaction between them and creates a deformation of the robot shape. The shaft of the robot is composed of a set of sections with different curvatures. The frontiers of the sections are defined by the starting and the end points of each straight 1.3. MOTIVATIONS OF THE THESIS 9 or curved part of all the tubes (cf. Figure 1.1.c). A more detailed description of the CTR kinematics is provided in Chapter 3. L s  L c  (a) Outer tube composed of straight part L s  and curved part L c  L s  L c  (b) Inner tube composed of straight part L s  and curved part L c  Sec  S e c  S e c  S e c  (c) Combination of the inner and the outer tube forms a two-tube Robot: the robot is composed of four sections(Sec  , Sec  , Sec  , Sec  ). The start and the end points of these sections are defined by the start and the end points of the straits and the curved parts of the tubes.

Figure 1

 1 Figure 1.1 -Conception of a concentric tube robot.

CHAPTER 2 .Figure 2

 22 Figure 2.1 -Actuation unit with serial translation blocks: a) structure, b) kinematic scheme.

CHAPTER 2 .

 2 FROM DESIGN TO CONTROL OF CONCENTRIC TUBE ROBOTS: A LITERATURE SURVEY in Figure 2.4.a and Figure 2.4.b respectively.

Figure 2

 2 Figure 2.4 -Two-tubes CTR with a) pneumatic actuation unit[START_REF] Cardona | A mri compatible concentric tube continuum robot with pneumatic actuation[END_REF] and b) piezoelectric actuation unit[START_REF] Su | A MRI-guided concentric tube continuum robot with piezoelectric actuation: A feasibility study[END_REF] non coupled.

CHAPTER 2 .Figure 2

 22 Figure 2.6 -CTR actuation unit with parallel translation blocks: a) structure b) kinematic graph.

Figure 2 . 7 -

 27 Figure 2.7 -RTC a) [Dupont et al., 2010a] b) [Gosline et al., 2012a] with three tubes decoupled and actuated in parallel.

Figure 2 . 8 -

 28 Figure 2.8 -Two-tubes CTR with parallel actuation unit[START_REF] Xu | Position control of concentrictube continuum robots using a modified Jacobian-based approach[END_REF]

Figure 3

 3 Figure 3.1, two mappings are required to calculate the pose of the end-effector of the CTR arm.

Figure 3 . 2 -

 32 Figure 3.2 -Two-tube CTR with three constant curvature sections, from the base to the tip: the first section of length l  delimited by the frame H 1 (attached to the base of the outer tube) and by the frame H 2 (attached to the distal extremity of the straight part of the outer tube) is composed of two straight tubes; the second section of length l  delimited by the frame H 2 and the frame H 3 (attached to the distal extremity of the curved part of the outer tube) is composed of a curved dominant outer tube and an inner tube with a straight and a curved part; the third section of length l  is composed of only one curved tube and is delimited by the frames H 3 and H t which is attached to the tip.

Figure 3

 3 Figure 3.3 -Section geometry, adapted from [R. J. Webster and Jones, 2010].

  works on the tubes design optimization of CTR do not consider the elastic stability in the objective function. The only work which formulates the elastic stability in the objective function is based on scalarization method which is the main tool used in existing multi-objective CTR design optimization works. The optimization objective in scalarization method is to maximize/minimize one single objective function composed of a set of weighted objective functions. In this method, the selection of function weights is crucial as the optimization results are greatly affected by them and could be misleading if these weights are improperly chosen. As an alternative optimization technique, Pareto grid-searching method can be used to avoid this problem and allow a straightforward interpretation of the optimization results. This chapter shows a three-tube CTR design based on Pareto grid-searching method in order to maximize the reachability and the elastic stability of the CTR within a specific curvature range dedicated to the deep anterior brain tumor removal surgery.

Contents

Figure 4

 4 Figure 4.1 -Evolution of CTR tubes design. (ISM: Infinite Stiffness Model, TFM: Torsion Free Model, TCM: Transmissional Compliant Model)

Figure 4 . 2 -

 42 Figure 4.2 -Example of a Pareto front in case of two objective functions.

  Figure 4.3 -Convex and non-convex Pareto curve.

  Figure 4.4 -Mapping from decision space to objective space.

Figure 4 . 5 -

 45 Figure 4.5 -Transcranial intervention for frontal lobe tumor resection [YOU]

  Figure 4.6.a-c. Another challenge of the targeted surgical application is the no-coaxiality of the insertion path with the tumor location increases the complexity of accessibility to the surgical site.4.5 Workspace characterization and anatomical constraintsGiven the three planes CT scan images (axial, coronal, and sagittal) of a real patient case study, a base frame F b is attached to the patient head with its origin O close to the nose and centered at the median between the two nostrils as shown in Figure4.6. The Z axis is oriented from the nose to the head back, the Y axis is oriented from the nose to the head top. The insertion point of the robot is at point P in Figure4.6, with P defined as the origin of the robot frame F r . The orientation of F r is decided by two rotations: φ which orients the robot from the nostril to the frontal sinus in the Y -Z plan of the frame F b and allows to avoid the nasal turbinates; and the ψ which orients the robot to the left or to the right of the cerebral falx depending on the tumor location.Based on CT scan images and with more clarification and orientation of the surgeon, the space of the passage used to introduce the CTR from the nostril to the frontal sinus is identified roughly as a straight tunnel with variable diameter (from ≈ mm for smaller ellipse diameter 4

Figure 4

 4 Figure 4.6 -Work-space characterization (a-c) and inner tubes insertion (d).

Figure 4 . 7 -

 47 Figure 4.7 -Robot shaft sections

Figure 4 . 8 -

 48 Figure 4.8 -CTR Robot: CAD and real platform.

F b :

 b base frame with origin O, F r : robot frame with origin P φ : X-rotation angle of Rob, ψ: Y -rotation angle of Rob after its X-rotation c: tumor center coordinates h  : solves the nonlinear differential equation system of the forward kinematic model to obtain the curvatures u(s) h  : numerical integration of the curvatures u(s) to define the robot shape G(s): homogeneous matrix of the pose of a point s from the robot shaft k  ,k  : curvatures of tubes T  and T  expressed in mm - ∆ q: sampling step of the rotation joints, ∆ s: sampling length χ: Cartesian tip position in F b , χ = [χ  , χ  , χ  ] a: radius of the sphere tumor =  mm l s : length of the straight part of T  =  mm, l s : length of the straight part of T  =  mm l c : length of the straight part of T  =  mm, l s : length of the curved part of T  =  mm l s : length of the curved part of T  =  mm ∆ k: curvature sampling step in the refining = mm -

Figure 4

 4 Figure 4.9 -Joint jumps

Figure

  Figure 4.10 -Stability evaluation

Figure 4 .

 4 Figure 4.12 -Stability/Reachability selection

  . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 80 5.2 Formulation of CTR position control task at actuator level . . . . . . . . . 81 5.3 Task-space position control of CTR with inaccurate kinematics . . . . . . 82 5.4 Simulation studies . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 86 5.5 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 94 80 CHAPTER 5. TASK-SPACE POSITION CONTROL OF CTR USING APPROXIMATED JACOBIAN

Figure 5 . 2 -

 52 Figure 5.2 -CTR Prototype: (top) Zoomed in picture of the concentric tubes (bottom) Whole system with actuators.

  Figure 5.3 -Workspace

Figure 5 . 4 -

 54 Figure 5.4 -Control schemes: a) Inverse kinematic control, b) Task-space control

  .4.a), the Cartesian errors are converted to joint errors trough the pseudo inverse of the Jacobian matrix which is based on the TFM. The corresponding actual task-space positions are calculated through the TTM 5.4. SIMULATION STUDIES 89 forward kinematics which is considered as the true model, and then task-space position errors are fed back to update the joint variables through TFM inverse kinematics in order to reach next desired position.

  5 and Figure 5.6. It is seen that the first type of simulation results are illustrated in Figure 5.5.a and 5.5.b. In Figure 5.5.c-f and Figure 5.6, kinematic parameter inaccuracies of deferent levels (on tube curvatures k  and k  ) are also included. CHAPTER 5. TASK-SPACE POSITION CONTROL OF CTR USING APPROXIMATED JACOBIAN IK Control with only modelling diference innacuracy.

  Figure5.5 -The performances of the two control strategies with modelling inaccuracy and uncertainty in the parameter k  .

Figure 5 . 7 -

 57 Figure 5.7 -TS Control with (k′  = % • k  , k ′  = % • k  ).The upper limit of the kinematic uncertainties bound is exceeded, and thus the tip position diverges.

Figure 5 . 8 -

 58 Figure 5.8 -CTR Shape During Control

  After literature survey and fundamental knowledge introduction, two contributions have been presented in this thesis. The first one is a concept study of using CTR for frontal lobe tumor removal. An optimization algorithm based on grid-searching Pareto technique which allows to avoid the limitations of scaliarization techniques used in literature has been used se-97 98 CHAPTER 6. CONCLUSIONS AND PERSPECTIVES lect tube dimension parameters in order to maximize the reachability of the desired workspace and to maximize the elastic stability. A new prototype has been designed to satisfy this specific surgical application requirements (way of tubes insertion, translations dimensions, number of tubes etc).

Table 2

 2

			.1 -Existing CTR prototypes.
	Prototype	Number Number	Number Actuation	MRI compatibility
		of arms of actuated of DOF coupling	and actuators
			tubes			nature
	[Webster III, 2007]	1	3	6	Rot/Trans	No
						DC Motors
	[Dupont et al., 2010a]	1	3	6	No	No
	[Gosline et al., 2012a]					DC Motors
	[Burgner et al., 2011]	2	3x2	6x2	No	No
	[Burgner et al., 2014]					DC Motors
	[Cardona, 2012]	1	2	3	No	Yes
						Pneumatic
	[Su et al., 2012]	1	2	3	No	Yes
						Piezoelectric
	[Comber et al., 2012]	1	3	6	Trans/Trans Yes
						Pneumatic
	[Burgner et al., 2012]	1	0	-	No	No
						Manual
	[Butler et al., 2012]	1	2	3	No	No
						DC Motors
	[Swaney et al., 2012]	4	4x3	4x6	No	No
						DC Motors
	[Xu et al., 2013]	1	2	4	No	No
	[Xu et al., 2014]					DC Motors
	[Hendrick et al., 2014] 2	3+2	6+3	No	No
						DC Motors
	[Burgner et al., 2013b] 1	2	3	No	No
						DC Motors
	[Azimian et al., 2014]	1	2	4	No	No
						DC Motors
	[Swaney et al., 2015]	1	2 + 1 wire 4+2	Non	No
						DC Motors
	between them (bending and torsion) occur. All these make the shape and the length of the
	robot arm change.
	It is obvious that a specific clinical task requires specific robot shape and length. The shape

Table 2 .

 2

		2 -Algorithms of tubes design optimization for RTC
	Optimization	Method	Objective	Optimization Model	Medical
	Algorithm	used	function	variables	used	application
	[Bedell et al., 2011]	GPS	• Min(lengths and curvatures of sections)	lengths and curvatures -free Torsion-	Intra-cardiac
			• obstacles avoidance	of sections	model
	[Anor et al., 2011]	PSA	• Min(lengths and curvatures of sections)	lengths and curvatures stiffness Infinite	Neurology
			• obstacles avoidance	of sections	model
		Nonlinear Min(unreachable points) Lengths	Torsion-	Skull
	[Burgner et al., 2011]	Simplex		and curvatures -free	base
				of the tubes	model
			Design and path	Length	Torsionally
	[Torres et al., 2012]	RRT	ensuring the access	and curvatures compliant Lung
			to the desired workspace of the tubes	model
		Nelder-	Max(coverage of the	Lengths	Torsionally Skull
	[Burgner et al., 2013a] -Mead	desired workspace)	and curvatures compliant	base
		Simplex		of the tubes	model
		Study		Curvatures	Torsionally
	[Ha et al., 2014]	of the	Max(elastic stability)	in terms of	compliant	-
		monotony		the lengths	model
	[Bergeles et al., 2015] Nelder-	• Min(Length and curvature of the robot)	Curvatures and lengths	Torsionally compliant	-
		-Mead	• elastic stability	of the tubes	model
	2.3 Kinematic models of CTR

Since the concept of CTR was proposed in 2006

[Sears and Dupont, 2006] [Webster III et al., 2006]

, many research efforts have been made to study its kinematics behavior [D. Caleb CHAPTER 2. FROM DESIGN TO CONTROL OF CONCENTRIC TUBE ROBOTS: A LITERATURE SURVEY and R.

Table 2 .

 2 4 specifying: the control strategy used; the robot length (the inner-most tube length

	2.5. CONTROL OF CTR

Table 2 .

 2 3 -Path planing algorithms for CTRs

	Planner	Method	Inputs	Outputs	Model	Medical
		used	of the algorithm	of the algorithm used	application
	[Lyons et al., 2009]	Constrained nonlinear programing	• Target position • Initial pose • Obstacles	Translations and rotations of the tubes	Infinite stiffness model	-
	[Lyons et al., 2010]	Constrained nonlinear programming • Obstacles • Desired position • Initial pose	Translations and rotations of the tubes	Bending only model	lungs
	[Torres and Alterovitz, 2011] RRM	• Desired position • Initial pose • Obstacles	Translations and rotations of the tubes	Bending and torsion -base Skull-model surgery
	[Torres et al., 2012]	RRT	Vertices • Initial configuration and edges • Design of tubes of the graph • Number of joint samples constituting	Bending and torsion Lung model
			• Obstacles	the trajectory
	[Bergeles and Dupont, 2013] RRT*	• Initial configuration A sequence • Number of samples of stable • Obstacles configurations	Bending and torsion Lung Model
				Translations
	[Torres et al., 2014]	RRG	• New actual and desired positions	and rotations allowing	Bending and torsion Lung
			at each iteration	obstacles	model
			• Obstacles	avoidance of the whole
				robot

Table 2 .

 2 4 -Control architectures of RTC. Differential kinematic model of concentric tube robot . . . . . . . . . . . . 48 3.4 Inverse kinematic model . . . . . . . . . . . . . . . . . . . . . . . . . . . 50 3.5 Bifurcation phenomenon . . . . . . . . . . . . . . . . . . . . . . . . . . 52 3.6 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 54

	Controller	Control law	Teleoperated L(mm) R-T F cycle (Hz) RMS (mm)
	[Webster et al., 2009]	IKC V. servoing	No	128,6	Yes 15	0,674
	[Dupont et al., 2010b]	PID	Yes	200	Yes 1000	4,2
	[Mahvash and Dupont, 2011] Force control, PD	Yes	150	Yes 1000	-
	[Burgner et al., 2014]	Inverse kinematic control Yes	275	Yes 400	1.5
	[Comber et al., 2012]	Sliding mode	No	-	-	-	-
	[Xu et al., 2013]	Inverse kinematic control No	400	Yes -	3,86
	[Xu et al., 2014]	Under external loads	No	303	Yes -	1,4
	[Azimian et al., 2014]	IKC under constraints	No	312	Yes 1000	-
	2.6 Conclusion
	More discussions about existing works in tubes design and control of CTR are also necessary
	before introducing the two main contributions of this thesis and therefore are provided in the
	introductions of Chapter 4 and 5.

This chapter has summarized the main works in literature and shown the directions of CTR evolution. Deep background on kinematics modelling is detailed in the next chapter.

Preamble

The evolution of CTR kinematic modelling has been introduced briefly in the previous chapter. Since kinematics serves as the basis for CTR design and control, in this chapter the mathematical frameworks of the most oftenly used forward, differential, and inverse kinematic models are presented in detail for the case of two tubes. The pre-shape of the tubes considered is the most oftenly used one which has a straight part followed by a planar constant curved part. Bifurcation, as a special phenomenon for CTR, presents big challenge for motion planning and implementation safety and is presented in the end of this chapter. Contents 3.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36 3.2 Forward kinematic model of concentric tube robot . . . . . . . . . . . . . 36 3.3

Table 3 .

 3 1 -Tubes parameters used in simulation of length and curvature effects on the elastic stability.

		Outer tube Inner tube
	Outer diameter (mm)	, 	, 
	Inner diameter (mm)	, 	.
	Young modulus (G Pa)		
	Shear modulus (G Pa)	, 	, 

considered at the level of tube design or at the level of path planing. The next two chapters of contributions are based on the kinematics introduced in this chapter and also take bifurcation into consideration for practical implementations.

  Results and discussions . . . . . . . . . . . . . . . . . . . . . . . . . . . 74 4.9 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 78

	CHAPTER 4. OPTIMIZATION OF CTR DESIGN FOR DEEP ANTERIOR BRAIN
	56	TUMOR SURGERY
	4.3 Techniques used to solve multi-objective optimization problem . . . . . . 59
	4.4 Targeted surgical application . . . . . . . . . . . . . . . . . . . . . . . . . 63
	4.5 Workspace characterization and anatomical constraints . . . . . . . . . . . 64
	4.6 CTR configuration selection . . . . . . . . . . . . . . . . . . . . . . . . . 66
	4.7 Formulation of grid-searching optimization algorithm for CTR . . . . . . . 69
	4.8
	55

4.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56 4.2 Multi-objective optimization using Pareto methods . . . . . . . . . . . . . 59

Table 4 .

 4 1 -Notations and Dimensions

  be-4.7. FORMULATION OF GRID-SEARCHING OPTIMIZATION ALGORITHM FOR CTR 73 cause the proximal sections are constrained by the navigation path and then should be of small curvature. Moreover, the distal sections which are composed of inners tubes dedicated to manipulating the tip should be of high curvatures. The second constraint is defined by an eligibility condition under the notation C eli . The confliction of the robot shaft with the anatomical structure is checked using the ellipses borders. Once the Y coordinate χ  of the robot is equal to the Y coordinate as an ellipse, this shaft robot point is checked to see if it is inside the ellipse or not by using the two radii of the ellipse to calculate the distance from the borders. The robot shaft is considered satisfying the eligibility condition if and only if each point of the robot which has the same Y coordinate of an ellipse is inside its corresponding

  Table 4.1. For each configuration, (k  ,k  ) Table 4.2 -Tubes parameters used in simulation of length and curvature effects on the elastic stability.

		T 	T 	T 
	Outer diameter (mm)	, 	, 	, 
	Inner diameter (mm)	, 	, 	, 
	Young modulus (G Pa)			
	Shear modulus (G Pa) ,  ,  , 

Table 5

 5

	.1 -Tubes parameters.
		Outer tube Inner tube
	Outer diameter (mm)	, 	, 
	Inner diameter (mm)	, 	.
	Young modulus (G Pa)		
	Shear modulus (G Pa)	, 	, 
	Length of the straight part (mm)		
	Length of the curved part (mm)		
	For inverse kinematics (IK) based control method (see Figure 5

ACKNOWLEDGEMENTS

I would like to thank my advisors, Prof. Philippe Poignet and Dr. Chao Liu, for giving me the opportunity to be a Doctor, fort guiding me, and for all the scientific discussions and

CHAPTER 4. OPTIMIZATION OF CTR DESIGN FOR DEEP ANTERIOR BRAIN TUMOR SURGERY

which generates a jump at this distal proximity [Dupont et al., 2010a]. As discussed in the previous chapter and shown in Figures 3.5 and 3.6, the elastic stability decreases once the curvatures and/or the lengths increase, and when it exceeds the stability limit the snapping phenomenon appears. This leads to a sudden variation of the robot end effector pose even if the rotation angle velocities of the tubes are low and smooth. This could lead to serious safety problem during the robotized intervention if it is not considered in the robot design or in the path planning.

The targeted surgical application considered in this work is deep frontal lobe brain tumors where elastic instability is of high risk to the operation safety. Therefore, in this work the elastic stability must be explicitly considered as an objective function in the optimization problem formulation. Instead of using scalarization methods in which the weights are pre-selected without any a priori informations, a Pareto grid-searching method has been used which incorporates all performance criteria within the optimization process and addresses them simultaneously to find a set of optimal designs in the objective space (the space constituted by the objective functions) [START_REF] Unal | Multi-criteria Design Optimization of Parallel Robots[END_REF], and therefore avoids the weight choice problem.

The objective of optimization in this work is to find the tube curvatures which suit best the specifications of the addressed surgical application. The robot structure and other parameters such as the lengths of the tubes and variation intervals of the curvatures to be optimized are pre-selected based on CT scan images in order to respect the anatomical constraints during the insertion and to be as close as possible to the optimal parameters with the help of surgeons.

Another contribution of this work is the use of forward kinematic model in the objective function instead of the inverse kinematics which has been used in [START_REF] Bedell | Design optimization of concentric tube robots based on task and anatomical constraints[END_REF] and [START_REF] Bergeles | Concentric tube robot design and optimization based on task and anatomical constraints[END_REF] to avoid the anatomical obstacles. The inverse kinematic model is complex and time consuming to calculate, and it requires a numerical optimization to be solved. The forward kinematic model used in this work is the torsionally compliant model developed in [Dupont et al., 2010a].

The optimal design presented in this chapter belongs to multi-objective (reachability and stability) optimization problem under constraint (anatomical constraints). The formulation of multi-objective optimization problem and the techniques used to solve it are briefly introduced in the next two sections.

MULTI-OBJECTIVE OPTIMIZATION USING PARETO METHODS

59

Multi-objective optimization using Pareto methods

Multi-objective optimization considers a set of objective functions to be optimized simultaneously. These objective functions are usually in conflict with each other which causes a trade-off in selection of the optimal solution. A perfect solution that optimizes all the objective functions simultaneously is in general not available. The solution making a balance between all the objective functions is selected from a set of solutions in the objective space (formed by the set of objective functions) called Pareto front or Pareto curve.

Consider the following mathematical formulation of a multi-objective optimization problem:

where n> is the number of objective functions, x ∈ R m is the decision variables vector of dimension m, Φ i=,...,n are the objective functions constituting the objective space, and ϖ is the set of constraints [START_REF] Caramia | Multi-objective management in freight logistics[END_REF]. The constraints are expressed in terms of the decision variables and can have an equality ϖ eq (x) and/or an inequality ϖ inq (x) form:

A decision variable vector x * ∈ ϖ is the Pareto optimal vector for a multi-objective problem if all other vectors x ∈ ϖ cause higher values for at least one objective function Φ i or have the same values for all Φ i . The Pareto front is formed by the collection of all efficient solutions (see the red curve in Figure 4.2). The trade-off between objective functions is concluded from this Pareto curve. and their corresponding (k  , k  ) are defined as the optimal decision variables.

Initialization:

Reachable(con f ig, ind) + + save(χ) to be used in C stab evaluation save(Reachable(con f ig)) to be used in C reach evaluation end if end if end if end for end for Calculate C stab (con f ig) for all configurations Calculate C reach (con f ig) for all configurations Algorithm 1: Grid-searching algorithm execution

Results and discussions

The tubes used in this study are made of Nitinol (% Nickel-% Titanum, see [EUR]) with the outer tube straight and the two inner tubes composed of two parts (a straight part followed by curved part). Their parameters are given in the Table.4.2.

CHAPTER 5. TASK-SPACE POSITION CONTROL OF CTR USING APPROXIMATED JACOBIAN

Rewrite V into the form

Since matrix Ĵ(q)M ĴT (q) is positive semi-definite, if α is chosen small such that the following matrix is positive definite

which is positive definite in q + α ĴT (q)∆ x and z.

Differentiate V with respect to time, it has V = ( q + α ĴT (q)∆ x) T M q + α qT M JT (q)∆ x + ∆ x T K p J(q) q +α qT M ĴT (q)J(q) q + ( q + α ĴT (q)∆ x) T k i z (5.13)

Substituting the closed-loop system dynamics (5.8) into (5. [START_REF] Burgner | A bimanual teleoperated system for endonasal skull base surgery[END_REF], it has

TASK-SPACE POSITION CONTROL OF CTR WITH INACCURATE KINEMATICS 85

By using the Jacobian estimation bound (5.5), it has

where λ min [•] denotes the minimum eigenvalue of a matrix.

Since

Inequality (5.15) can be further written as

(5.17) So if α is chosen small enough, K p and K d are chosen large enough, the following conditions can be satisfied for a certain range of kinematics estimation error β :

Considering that V is positive definite, according to Invariant Set Theorem [START_REF] Slotine | Applied Nonlinear Control[END_REF], it can be concluded from (5.17) that q →  and ∆ x →  asymptotically.

Remark 5.1: From conditions (5.18) and (5.19), it can be seen that the proposed method can handle kinematics error but with certain range. When the kinematics estimation error is too big, meaning β is big, the two conditions may no longer be satisfied and the system performance is not guaranteed, which is logical that automatic controller cannot handle infinite large modelling errors without online updation or adaptation that lead to more complicated

CHAPTER 5. TASK-SPACE POSITION CONTROL OF CTR USING APPROXIMATED JACOBIAN

To evaluate the performances of each control strategy, the Cartesian Root Mean Square Errors (RMSE) along X,Y, Z axis of task space axis are calculated in all cases of study presented above and shown in the

3.9103 2.9382 3.6250 6.0880 TS Control with:

4.1804 3.5221 3.9358 6.7358 TS Control with:

It can be seen from the Table 5.2 that the Cartesian RMSE for the case of task-space control is in range of . mm and . m when the kinematic uncertainties increase. This RMSE is in the rage of  mm and  m for the case of inverse kinematics control. This means that the impact of kinematics inaccuracy is much higher on the inverse kinematic control (from  to  mm) than on the proposed task-space controller (from . to . mm). From Figure 5.5, Figure 5.6, and Table 5.2, it can be concluded that with the presence of kinematics modelling error, the task-space approximate Jacobian control method provides better positioning performance over inverse kinematics based control method.

It is also noted that once the kinematic inaccuracies exceed certain limit, the position diverges as can bee seen from Figure 5.7. This observation is consistent with the theoretical

SIMULATION STUDIES

APPENDIX

A LIRMM CTR Prototype

In LIRMM, a latest CTR prototype has been designed and assembled for future experimental test. The complete system setup is shown in List of Publications

List of Figures

3.1 Tow mappings defining the forward kinematic model: f  from the joint space (α i , ι i ) to the configuration space (ρ j , l j , ϕ j ) and f  from the configuration space to the task space (tip pose X). This figure is adapted from [R. J.

List of Tables