Faisal Islam

Sébastien Joannès
email: sebastien.joannes@mines-paristech.fr

Anthony Bunsell

Lucien Laiarinandrasana

TWENTY-SECOND INTERNATIONAL CONFERENCE ON COMPOSITE MATERIALS (ICCM22) ADAPTATION OF WEIBULL ANALYSIS TO REPRESENT STRENGTH BEHAVIOUR OF BRITTLE FIBRES

Keywords: Fibre strength, statistics, preselection, truncation

Fibre properties are used for estimating the damage and strength behaviour of composite materials and structures, as the mechanical properties of fibre reinforced composites are primarily dependent on those of fibres. Fibre reinforced composites owe their outstanding properties to their constituents, mainly the fibres which reinforce matrices in composites. The accuracy of the predictions made by composite strength models would be highly dependent on the accuracy of the fibre strength distribution which is used as input. Any error in input fibre strength properties would lead to inaccurate model predictions. It has however been highlighted by many studies that there is discrepancy in experimental data even for the same type of fibres. The fibre strength data used to determine fibre strength distributions may be inaccurate. Due to the fibre preselection effect, experimental data is usually incomplete and does not appropriately represent the true fibre strength variation. An adaptation to the Weibull analysis has been presented to predict a more accurate statistical distribution for the representation of the strength behaviour of fibres inside composites. This analysis considers the preselection effect, which is used in modelling the fibre strength behaviour. The described analysis has the potential to significantly improve the accuracy of fibre strength distribution determined using experimentally generated fibre strength data.

INTRODUCTION

When fibres are used to reinforce polymeric matrices, it results in significant improvements of their mechanical behaviour. This provides benefits such as light weight, high strength to weight ratio, excellent weathering stabilities and enhanced dimensional stabilities to the resulting composite material. These properties make them ideal to be used for very critical structural applications such as pressure vessels and components for aerospace applications. Due to the widespread use of fibre reinforced composites, several computational models have been developed for predicting the behaviour of composite materials and structures. These models use constituent properties as inputs and predict the properties of composite materials and structures. The reliability of the model predictions depends on the quality of the input data used, mainly that of the properties of fibres, since they are the principal load bearing constituent of the composites. Thus, for failure modelling and to simulate the effective properties of composite products, detailed information about fibre tensile strength distribution is required. Due to the widespread use of fibre reinforced composites, it is very important to determine the characteristics of fibres and matrices and their interaction with each other in order to have a clear understanding of their effect on the composite material and the end-products. Correctly capturing fibre strength statistics and failure kinetics is also crucial to the success of composite stochastic strength models. The accuracy of the predictions made by composite strength models would be highly dependent on the accuracy of the fibre strength distribution which is used as input. Any error in input fibre strength properties would lead to inaccurate model predictions.

In practice, the characterization of fibre properties is known to be very challenging. Their very small cross-sections, morphological variations, or simply the sampling randomness and sample size used for analysis, directly affect determination of fibre properties especially tensile strength. They make the characterization uncertain and may lead to significant variations in results. Apart from the limited precision of measurement tools used, uncertainties could also arise from many other sources such as the inconsistency of the individual performing the measurement, choice of sample size used, sampling error, etc. The representativeness of the experimental fibre tensile strength data set that is used for determining a statistical distribution to characterize the tensile strength behaviour can be a major source of uncertainty in distribution parameters. Implementation of advanced data processing techniques on fibre strength results obtained from experiments may help in drawing attention to unknown limitations in the present system of analysing test results. It may also significantly enhance the reliability of the obtained results and help in extracting more information from the limited resources available.

Technical fibres which are commonly used as reinforcements in composite materials are usually brittle in nature. The stochastic distribution of defects inside them controls their strength. Representing the strength of such fibres thus requires a statistical approach. The most popularly used statistical function for representing strength of brittle fibres is the Weibull distribution. The parameters of this distribution represent the strength variation of fibres and are determined by appropriately fitting experimentally generated fibre strength data with this statistical function. Depending upon the type of fibre strength data, authors usually choose between the generic 3-parameter Weibull distribution, and a specific 2-parameter distribution obtained by fixing the value for one of the parameters. Using either of the two functions is associated with corresponding limitations which have been discussed in this article. It has often been observed that single fibre strength data that is experimentally generated does not conform to the variations of the standard Weibull functions, i.e. there is a misfit between the experimental data points and the Weibull model, on the graphical plot. To overcome this misrepresentation, many authors have proposed using different modified versions of Weibull distributions which have been shown to statistically fit the fibre strength experimental data well. Based on the goodness of fit, conclusions on the physical nature of the fibre strength behaviour have been made. This may have serious consequences as correlation may not necessarily mean causation.

The problem does not lie entirely with the statistical distribution but mainly with the data used for determining the distribution. The fibre strength data which is experimentally generated using the single fibre testing process is usually incomplete and does not accurately represent the strength of all the fibres present in composite materials and structures. This is due to the effect of fibre preselection, the causes of which will also be discussed in this article. To address this issue, a generalized Weibull analysis will be presented. This analysis takes into consideration the fibre preselection effect which would be analytically represented in the Weibull distribution. Analysing experimentally generated fibre strength data using the presented methodology can significantly enhance the accuracy of the results obtained and would help in extracting more useful information from the limited resources available. This would therefore also improve the reliability of the predictions made by composite strength models on the behaviour of composite materials and structures.

This work aims at finding an appropriate distribution that could represent the actual fibre strength behaviour accurately. This advanced knowledge may enable the design of new and improved composite materials and products with more accurate information about their strength behaviour.

The strength of a single fibre cannot be captured in one single average value. Single fibres normally exhibit a distribution of strengths that could be represented by using a probabilistic approach. The statistical function that best describes strength behaviour of brittle fibres such as carbon is the Weibull distribution.

The Weibull distribution is considered to be a good representation of the weakest link behaviour of brittle single fibres. The classical Weibull distribution function for fibre strengths is given by Eq. 1.

P(σ) = 1-exp{-( L L 0 ) ( σ-σ u σ 0 ) m } (1)
Where, P(σ f ) is probability of fibre failure, L being the characteristic gauge length, L 0 the reference gauge length, σ the fibre strength, σ u the threshold strength, σ 0 the scale parameter and m the shape parameter or Weibull modulus. Eq. 1 is also commonly known as the 3-parameter Weibull distribution. The value for σ u is usually set to be zero and the resultant is known as a 2-parameter distribution (Eq. 2).

P(σ) = 1-exp{-( L L 0 ) ( σ σ 0 ) m } (2)
Certain authors have speculated that multiple strength determining flaw populations may exist inside fibres and have recommended the use of a multimodal Weibull distributions to represent the strength of a fibre population. A power law accelerated Weibull distribution has also been proposed by some authors. It tries to represent the gauge length dependence of fibre strengths by including an additional parameter in the distribution. However, these models are mainly based on curve fitting and do not have a strong physical explanations.

EXPERIMENTATION AND RESULTS

T700 carbon fibres are very popularly used in fibre reinforced composites for structural applications and are brittle in nature. Single fibre tensile tests using T700 fibres were conducted to determine the tensile strength of the fibres. The procedure described in ASTM C1557-14 was followed for all tests. The fibre diameters were carefully measured using a laser diffraction system. The failure force was divided by the fibre cross-sectional area to determine the failure stress or strength of the fibre. The test was performed at the recommended crosshead speed of 1 mm/min. A tensile strength data set for 50 fibres was determined at a gauge length of 30 mm. It is a very common practice to use the 2-parameter distribution to represent single fibre strengths and for deriving conclusions on the physical behaviour of fibres. However, when a sufficiently large single fibre strength data set is represented on a Weibull plot using a 2-parameter distribution, there is a clear misfit between the experimental points and the Weibull model, as shown in Fig. 1. This deviation can partially be explained using the fibre preselection effect, which means that the weak fibres break during the process of fibre extraction and testing. This has also been observed by other researchers. The weak data points are thus not represented in the data set used for analysis. This elimination of weak fibre strength values from the data set causes a deviation in linearity on the Weibull plot, as shown. The 2-parameter Weibull distribution is unable to capture this non-linear behaviour accurately. It may be tempting to use the 3-parameter Weibull distribution (Eq. 1) as it would be able to present a distribution that would fit the experimental data points very well. However, this distribution would assume that the minimum strength value in the experimental data set represents the weakest fibre in the population which is an incorrect assumption. It is not desirable to only obtain the best fit distribution but the objective is to find a representation of the actual strength behaviour of the fibres.

DISCUSSIONS

Using the theory of conditional probability, a generalized Weibull distribution can be derived. This distribution is capable of representing the preselected fibre strength data obtained from experiments. The generalized equation is given by Eq. 3.

P(σ) = 1-exp{-( L L 0 ) ( (σ-(1-ω)σ r ) 𝑚 -(σ r -(1-ω)σ r ) 𝑚 σ 0 𝑚 )} (3)
Where, σ 𝑟 is a truncation limit and ω is a constant. The quantities,σ 𝑟 , ω and the threshold strength σ 𝑢 are related by the empirical relationship given by Eq. 4. The value of ω can be determined by a user based on the observed behaviour of fibre strengths.

σ 𝑢 = (1 -ω)σ 𝑟 (4)
For certain specific conditions, Eq. 3 simplifies to return the classical 2 and 3-parameter Weibull distributions. (Eqs. 1 and 2). The generalized Weibull distribution contemplates the incompleteness of the experimental data set and can further be used to determine an accurate representation of the actual fibre strength behaviour, otherwise not possible using the classical distributions.

The inadequacy of the popular statistical distributions for representing fibre strength has been highlighted and the experimental data set is analyzed statistically. This is due to the problems in experimental data set and not necessarily due to the inadequacy of the statistical distribution used.. Experimentally generated fibre tensile strength data sets are not always reliable representatives for the strength of the fibres inside a fibre bundle, composite material or a composite structure. Weaker fibres are unable to be tested due to the effect of fibre preselection. They are thus not represented in the experimentally determined data set. This results in a truncation in the obtained fibre strength data which also propagated into the resulting strength distribution. Analyzing this incomplete data with the standard Weibull functions may result in incorrect conclusions about the physical nature of the strength behavior of fibres. This error would also propagate into the predictions made by composite strength models which use the fibre strength information as input.

It is shown that inclusion of a truncation parameter in the Weibull distribution could help in estimating the strength behavior of a given fibre population. Inclusion of a user defined constant ω allows the experimental results for different types of fibres to be addressed. By varying ω, the effect of the fibre strength distribution on the composite strength models can also be simulated. Ideally, this distribution may also have been obtained if the entire fibre strength spectrum had been experimentally accessible. However, since that is practically impossible, the presented truncated Weibull analysis can be used to represent the actual strength behavior of fibres inside composite materials.

Figure. 1 .

 1 Figure. 1. Weibull plot with 50 fibre strength data points for T700 carbon fibres (Determined at Dia-Stron Ltd.).

ACKNOWLEDGEMENTS

The research leading to these results has been done within the framework of the FiBreMoD project and has received funding from the European Union's Horizon 2020 research and innovation programme under the Marie Sklodowska-Curie grant agreement No 722626.


Thibault Sénac

Arnaud Lelevé

Richard Moreau

Laurent Krahenbuhl

Florent Sigwalt

Christian Bauer

Quentin Rouby

Designing an accurate and customizable epidural anesthesia haptic simulator

Epidural anesthesia, despite being a relatively common medical procedure, remains quite demanding in terms of skills as it is mostly blind and thus heavily reliant on the haptic sensations. Although some training support solutions exist, anesthetists consider them mostly inefficient or impractical. A few attempts at creating a simulator for this particular procedure exist but each one lacks one of the important requirements of the procedure. This article introduces a haptic simulator featuring a more complete and realistic simulation of the procedure than we could observe in existing simulators. The simulator is composed of a generic electrical haptic interface coupled with a pneumatic cylinder.

I. Introduction

A. Medical context

Epidural anesthesia is a commonly used medical gesture performed mainly during childbirth delivery, but it is very complex and hard to learn. This is mainly due to the fact that the procedure is mostly performed blind, relying solely on haptic cues and their interpretation by the anesthetist. To perform such a gesture, the practitioner has to insert a Tuohy needle between two vertebrae while injecting a fluid using an epidural syringe. The Tuohy needle has to cross several physiological layers, namely derma, supraspinous ligamentum, intraspinous ligamentum, and ligamentum flavum, to finally reach the epidural space. Throughout the insertion, the anesthetist experiences an increasing resistance felt both from the needle insertion and from the fluid injection. This resistance reaches its maximum when the needle tip crosses the ligamentum flavum, and suddenly plummets, creating what is commonly called the loss of resistance (LOR) principle. It is this phenomenon that allows physicists to know they reached the epidural space and this principle is specific to the epidural anesthesia procedure. Due to the relatively high forces exerted simultaneously on the needle and on the plunger of the syringe, the procedure is really demanding in term of precision as the epidural space is usually only about 4 mm wide, and about 6 mm deep.

In consequence, this particular medical gesture has a quite steep learning curve and may require up to 90 attempts to be performed with a success rate of only 80% [START_REF] Vaughan | A review of epidural simulators: Where are we today?[END_REF], which does not suffice for healthcare applications. The need for an efficient teaching tool has been strengthened by public institutes such as the French H.A.S. [START_REF] Granry | État de l'art (national et international) enmatière de pratiques de simulation dans le domaine de la santé[END_REF] (Haute Autorité de la Santé), which prescribes to "never do the first time on a patient".

B. State of the art

In the past years, some training solutions to support this teaching have been developed such as manikins or complete robotic simulators. Most of the proposed solutions have been listed and compared in [START_REF] Vaughan | A review of epidural simulators: Where are we today?[END_REF], which considers 17 manikin based solutions and 14 computer based ones. According to this study, despite promising solutions, none of them responds to every training requirement. On the one hand, manikins generally offer an appropriate level of realism concerning the general look and feel of the back of the patient. However, steep prices, the lack of easy customization, the non-repeatability of the attempts, the lack of objective feedback for the trainers to perform a supervisory control, usually limit the practicability of these solutions. On the other hand, the computer based solutions lack realism in the patient back reproduction. Nevertheless, computer based simulators allow for a deeper and easier customisation of the type of patient and objective assessment. Some computer based solutions have been developed to help trainees learn the epidural anesthesia. We can highlight, for instance, the works of Manoharan et al. [START_REF] Manoharan | Design and validation of an epidural needle insertion simulator with haptic feedback for training resident anaesthesiologists[END_REF], Magill et al. [START_REF] Magill | A Novel Actuator for Simulation of Epidural Anesthesia and Other Needle Insertion Procedures[END_REF], Dubey et al. [START_REF] Dubey | Biomedical Engineering in Epidural Anaesthesia Research[END_REF] or Thao et al. [START_REF] Thao | Development and Evaluation of an Epidural Injection Simulator with Force Feedback for Medical Training[END_REF], who developed haptic simulators for this particular procedure. These simulators however, provide an incomplete representation of epidural anesthesia. Concerning needle insertion, they often feature a custom haptic interface to reproduce a simplified version of the necessary haptic cues of the procedure. The procedure may, for example, be limited in terms of degrees of freedom, or only consider crossing the different physical layers without taking into account a potential bone contact. Regarding the LOR, the simulations are very simplified (if present), only using on/off solenoid valves.

In the literature, needle insertion modeling has been thoroughly studied. The insertion forces are often described as a sum of various contributions [START_REF] Okamura | Force Modeling for Needle Insertion Into Soft Tissue[END_REF] which are: a constant cutting force, a friction force and an elastic force. Moreover, as the epidural space is located in the spine, one may hit a vertebrae during the insertion of the needle. In that case, simulators must reproduce the feel of an interaction with a hard body. Thereon, the most common approach is to use a proxy system (also called god object [START_REF] Zilles | A constraint-based godobject method for haptic display[END_REF]). This solution is usually quite efficient and simple for surface exploration but has shown some limitations when it comes to representing collisions, as Kuchenbecker et al. [START_REF] Kuchenbecker | Improving contact realism through event-based haptic feedback[END_REF] pointed out in their study. They therefore proposed using event based haptic feedback to generate a vibration cue when collision occurs. Then, once the contact is established and the transient vibration force has damped, it goes back to a classical proxy based approach to generate the resting force.

C. Required features

Regarding epidural anesthesia simulator, Vaughan et al. [START_REF] Vaughan | A review of epidural simulators: Where are we today?[END_REF] has exposed the requirements for a perfect epidural anesthesia simulator. According to them, the simulator has to be customisable and should contain a haptic interface paired with a 3D interface. Finally, the study insists that the simulator should emulate not only the needle insertion forces, but also the LOR felt using a syringe. However, most of the current simulation options only reproduce what is called 100% LOR, which is encountered only in the "average" and "easy" cases.

D. Contributions of this study

In this paper, we introduce an hapting training simulator for epidural needle insertion covering the requirements proposed in [START_REF] Vaughan | A review of epidural simulators: Where are we today?[END_REF]. It features a customizable needle insertion haptic simulator coupled with a LOR syringe haptic simulation. The simulator uses a off-the-shelf haptic interface handling the reproduction of the needle insertion, and a pneumatic cylinder which reproduces the classical LOR syringe behaviour during this kind of procedure. The novelty lies in the use of a pneumatic cylinder fitted with a servovalve to precisely reproduce the LOR forces, instead of a simple solenoid valve in existing solutions. Also, an enhanced model of needle insertion is introduced in order to take into account the potential contact with bones. This simulator reproduces realistic epidural anesthesia procedures, records data to assess the trainee level, and allows the repeatability of the attempts. This paper is organized as follows. Section II introduces the way this simulator has been designed. Section III depicts the experimental setup for validation. Section IV reports the results and provides a discussion. Finally, some concluding remarks and future directions are provided in Section V.

II. The simulator and its implementation

In this simulator, the haptic electric interface reproduces the haptic cues generated by the needle insertion only. The pneumatic cylinder is employed to emulate the syringe which is used to detect the end of the procedure through the LOR approach.

A. Prototype setup

A latex based human back phantom is used, to reproduce the back of the patient. Some landforms represent the vertebrae. As the choice of the insertion point is out of the pedagogical objectives of this simulator, the needle insertion location is static. The prototype features a Virtuose™ 6D [START_REF]Virtuose™ 6d -HAPTION SA[END_REF] from Haption ® as the haptic interface for the Tuohy needle driving. The needle is reproduced by a plain metallic rod which is attached to the tip of the haptic interface and crosses the fake patient back. An Airpel ® anti stiction double acting pneumatic cylinder (M9D50.0UT1 [START_REF]Double-Acting Archives[END_REF]) with a shape close to the real LOR syringe, represents the latter and is attached to the "external" end of the rod. In order to keep the instrumentation of the cylinder as light as possible and as small as possible, we opted for a magnetic position sensor (MPS-064TSNU0 [START_REF]MPS-064tsnu0 | Magnetic cylinder sensors SICK[END_REF]) from SICK ® and two miniature pressure sensors (EPRB-1 [START_REF]Miniature Ratiometric Pressure Transducer | EPRB Series | TE Connectivity[END_REF]) from TE ® . A dSpace ® DS1103 controller commands a Festo MPYE ® -5-M5-010 B servo-valve to control the air mass flow rates entering this cylinder. An Arduino Uno ensures the communication between the haptic interface and the pneumatic cylinder controllers. The prototype, manipulated by a trainee, is visible in Fig. 1. We added some mechanical parts to bring the prototype closer to a real case in term of look and feel. A set of 3D printed parts dresses the pneumatic cylinder to resemble a real syringe. Moreover, a curved plate represents the back of the patient.

The control laws concerning the pneumatic cylinder were written with Simulink. A real-time user interaction with them is performed with ControlDesk. The haptic interface application was developed in C/C++ using the Haption API coupled with a physical simulation engine based on the ODE library [START_REF]Open Dynamics Engine -Browse /ODE/0.13 at SourceForge[END_REF], and with a virtual 3D environment.

B. Pneumatic cylinder control 1) Objectives:

The aim of the pneumatic cylinder is to reproduce the forces resisting to the injection of a fluid throughout the insertion of the needle. The main objectives are:

• clearly reproducing the LOR phenomenon, creating different resistant forces according to the physiological layer simultaneously crossed by the virtual needle,

• allowing the user to empty (and potentially refill) the syringe.

The pneumatic cylinder is fitted with a servo-valve instead of a simple solenoid valve, to reproduce a larger variety of haptic cues during the LOR experience. Indeed, the variability of the patient morphology implies various force levels as well as various speeds of transition when the force plummets at the end of the procedure. Also, it is necessary to emulate the emptying/re-filling of the virtual syringe through the pneumatic cylinder and its control law.

2) Cylinder control:

To implement the LOR behavior, two control laws were needed, as well as an algorithm selecting the right law and switching them in a stable fashion. The first control law needed to produce a resistant force (as the real syringe which hardly empties itself as long as the needle has not reached the epidural space) and the second control law had to keep the plunger at a set position (to set the plunger at the right positions during the procedure). We apply the control laws defined in [START_REF] Herzig | Nonlinear position and stiffness Backstepping controller for a two Degrees of Freedom pneumatic robot[END_REF], [START_REF] Senac | Control laws for pneumatic cylinder in order to emulate the Loss Of Resistance principle[END_REF] based on the A-T transform [START_REF] Abry | Nonlinear position control of a pneumatic actuator with closed-loop stiffness and damping tuning[END_REF] and a backstepping approach to control either the position or the pneumatic force generated by the system. We then added a switching algorithm that selects the control law to activate in real-time. To ensure that the system remains stable during the whole procedure, we use an average dwell time approach to set a minimum switching time between the two laws in our switching algorithm [START_REF] Hespanha | Stability of switched systems with average dwell-time[END_REF]. It guarantees that the laws stay active long enough for the whole system to remain stable.

C. Haptic interface control 1) Objectives:

The main objective of the electric haptic interface is to reproduce the forces encountered when inserting a needle in tissues. According to Okamura et al. [START_REF] Okamura | Force Modeling for Needle Insertion Into Soft Tissue[END_REF], one can describe the insertion forces in soft tissues (muscles, ligaments and skin in particular) as a sum of the following contributions:

F needle = F stif f ness + F f riction + F cutting ( 1
)
where F stif f ness represents the force corresponding to the puncture of a tissue, F f riction the force generated by the friction between the tissues and the needle and F cutting the force at the tip of the needle once the tissue has been punctured. In our case, the stiffness part is impactful mostly before the needle has punctured the derma. The cutting force depends on the tissue the needle tip crosses, and the friction force grows with the distance traveled inside the tissues. This behaviour is suitable for a standard insertion, where the needle only crosses ligaments, muscles or skin for example but it has to be enhanced to represent some bone contact.

To reproduce the haptic cue corresponding to a contact between the needle and a bone, we referred to [START_REF] Kuchenbecker | Improving contact realism through event-based haptic feedback[END_REF] which provides various desired force patterns to reproduce whenever a contact occurs. Out of these patterns, we chose to re-use the one based on a decaying sine as it is presented as being relatively simple, identifiable as a solid-solid contact for the user, and quite realistic. The trajectory is the following :

F contact = -Ae ln(0.01)t/d sin(2πf t) v in , for 0 < t ≤ d
(2) with A being a nominal magnitude, d the duration of the pattern, f a set frequency, and v in the velocity vector before hitting the surface. Moreover, once the needle collides with a bone, the needle tip cannot slip from the contact point. It has to be pulled back. Lastly, once the needle has started to cut through tissues, we consider that it remains on a set straight trajectory. This trajectory may be considered as curved in some cases depending on the type of needle used. However in our case, with the pierced tissues being relatively dense and because of the form of the needle tip, the trajectory remains relatively straight. Each of the previously forces are further detailed in the following part.

2) Implementation:

The Virtuose 6D is designed to operate mostly using a proxy based/virtual mechanism system. In this type of simulation, a virtual object is attached to the haptic interface. This virtual object receives the position of the interface handle and sends back forces to the interface through a virtual mechanism (usually a mass, spring, damper mechanism, as show in Fig. 2). These forces may be reaction forces when the virtual object touches another another object in a virtual scene, or forces rendering the gravity or some viscosity. This kind of approach is widely used to reproduce virtual object contacts as it ensures the stability of the interface as long as the stiffness and the damping are appropriately chosen [START_REF] Diolaiti | A criterion for the passivity of haptic devices[END_REF], it limits vibrations and it has low CPU usage [START_REF] Zilles | A constraint-based godobject method for haptic display[END_REF].

In this application, the Virtuose 6D control soft-

Interface Position

Virtual mechanism

Virtual object position

Object surface

F d
Figure 2: Example of a proxy and virtual mechanism on a contact situation ware features two ways to generate the desired force applied by the haptic interface on the fake needle. On the one hand, we make use of the aforementioned proxy virtual object approach to model the virtual needle. The latter interacts with several other 3D virtual objects in order to generate each part of the needle insertion forces: the various physiological layers, two vertebrae and a plane to generate the stiffness (F stif f ness ) part of the reaction force. The friction force has not been implemented so far, as it did not seem crucial on initial examination. On the other hand, some haptic cues require fast force patterns to be fully efficient (typically F stif f ness or F contact ) as they are based on a sensation of vibration. To reproduce these cues, we superimpose the corresponding force patterns directly on the haptic interface. To implement the cutting force we used nested boxes, as shown in Fig. 3, to generate a constant force opposed to the current motion of the virtual object. According to Okamura et al. [START_REF] Okamura | Force Modeling for Needle Insertion Into Soft Tissue[END_REF], the force is created as follows:

F cutting = F layer • v || v|| , if ||v|| > , ( 3
)
with F layer being the sum of each force from each layer crossed and a velocity threshold to trigger this force behavior. We used a nested environment instead of stacked spaces to avoid any issues at the interfaces between each layer, where one could cross either an empty zone or a heterogeneous area which would create undesirable high force variations. Concerning the stiffness force, according to Okamura et al. [START_REF] Okamura | Force Modeling for Needle Insertion Into Soft Tissue[END_REF], it has the following formulation:

F stif f ness = (a 1 • z + a 2 • z 2 ) • n, ( 4
)
n being the normal to the plane the needle is starting to cross, a 1 and a 2 two positive constants, and z the distance between the object center and its orthogonal projection. This force is computed using the virtual object position but is applied directly to the haptic interface without going through the proxy because we need the feel of a rupture of the capsule when the needle has entered enough.

The contact force is computed according to (2) and also applied directly to the haptic interface, using the position and velocity of the coupled virtual object. At the beginning of the procedure, we attach the virtual object to a virtual line which is aligned with the orientation of the fake needle (which coincides with the virtual needle), in order to keep the needle trajectory straight. We let the virtual mechanism generate the forces needed for the users to realize that they are locked on this straight trajectory. Finally, we parameterized the virtual hard object so that it can roll without slipping. When coupled with the set straight insertion trajectory, this helps users realize they are stuck on a bone when the virtual object contacts one.

III. Experimental setup and first tests

A. Objectives

First, we tried the prototype with skilled anesthetists to make sure the simulation was accurate enough. It allowed us to apply some necessary adjustments to some simulation parameters, such as the values of F layer or the resistant force created by the pneumatic cylinder. Once the parameters were set to obtain realistic simulations, we proceeded to recording of complete procedures. The aim of the first experimental study was to know whether it would be possible to differentiate an experienced user from a novice through the use of this simulator.

B. Test settings

We set up three kinds of patients for the first test sessions. The first patient type was the average patient, which represents an "easy" case and serves as a basis from which we derived the others. The parameters changed for the other patient types were: the length of the derma layer, the length of the epidural space layer, the value of F dura , and the LOR speed (i.e. the speed of transition between the highest resistant force in the ligamentum flavum to the lowest one in the epidural space). The patient characteristics are provided in table I. Their dimensions derive from measures on MRI images, the forces from an article by Tran et al. [START_REF] Tran | Instrumentation of the Loss-of-Resistance Technique for Epidural Needle Insertion[END_REF] and the rest from the knowledge and feedback from the skilled anesthetists working on the project.

The procedure of the test was:

1) quickly present the aim of the epidural anesthesia and the risks it bears, 2) show an example of the correct way to perform the procedure (the main objective was to introduce the right way to handle the syringe),

3) let the subject perform the procedure himself without recording so they can get an idea of the forces they will have to face, 4) record at least one attempt on each patient case we implemented.

The tests were done with a relatively low number of subjects, namely 2 experts and 6 novices.

IV. Experimental results

A. Data analysis

Once a few tests were recorded, we analyzed the data featuring a total of 24 variables: including in particular the mean and variance of both velocity and acceleration of the haptic interface handle along three axes, the mean and variance of both velocity and acceleration of the cylinder plunger. We also recorded the number of bone contacts, needle insertions, the insertion angles, the proportion of syringe volume used, as well as the overshoot which represents the distance traveled inside the epidural space. We also developed an indicator (syr rate ) based on the amount of syringe emptied relative to the distance traveled by the needle in the back. The formula for this indicator is:

syr rate = Proportion of syringe emptied d (5)
with d being the distance travelled by the haptic interface from the puncture of the derma to the end of the procedure along the elastic plane normal vector. Finally, we recorded whether the attempt ended up puncturing the dura or not as a way to describe whether it was a success or a failure.

As a first result we decided to analyse the success rate per patient type, the velocity of the needle tip along the insertion path, and finally the aforementioned syr rate to try differentiating a novice user from an experienced one.

B. Results

The experienced users who tried the simulator were quite satisfied with the overall realism of the proposed simulator. Some room for improvement was highlighted in some areas such as bone contact, ability to steer the needle while in the derma and some variability in the space between each vertebrae. The results shown in Fig. 4 highlight some patterns that separate the unskilled users from the trained experts. The success rate shows that the only users who failed were the unskilled ones. The overall results show around 60% of success rate for the novices. Their success rates are better on the calcified/average cases and worse on the overweight ones. This indicates that the LOR speed played a great role in how the unskilled users detected the epidural space. This strengthens the need of a customisable LOR speed instead of the existing binary 100% LOR.

Concerning the speed, the data show that the experts perform the gesture with a relatively low but constant speed, regardless of the patient type. The unskilled users showed a much more variable speed that greatly depends on the patient type. Hence, the consistency in speed may indicate a stronger mastery of the procedure.

Finally, the experts tend to have a much higher syr rate (indicator defined in (5)) than the novices, showing that they usually use a greater portion of the cylinder length (d being constant for a given patient type). This might be a way for the skilled users to gather as much haptic information as possible. This may help them guess the current position of the needle tip with much more precision, besides detecting the LOR. Typically, they are able to differentiate each crossed layer, from the force feedback provided by the syringe. Furthermore, concerning the novices, the standard deviation of the syr rate is also quite high. This shows a wide variety of use.

It is noteworthy that the low number of experts attempts may have caused the very high standard deviation for the calcified case.

V. Conclusion, perspective and future work

In this paper, we introduced a realistic and customisable haptic simulator dedicated to the hands-on training on the epidural anesthesia procedure, based on the ideal simulator depicted by Vaughan et al. [START_REF] Vaughan | A review of epidural simulators: Where are we today?[END_REF]. The overall results of this first study are quite encouraging and indicate that such a simulator might be realistic enough to provide an efficient training tool in the future. However some areas need improvements, such as introducing a very distinct feel of cutting through the ligamentum flavum which, according to our experts, feels like friction and cracklings. Moreover, according to their feedback, when the needle is in the ligamentum flavum, the syringe plunger feels like it is locked. Finally, we would like to perform an automatic classification of the user once the tests are completed. This requires many more trials to generate a training set for such an algorithm.

VI. Acknowledgment

Authors would like to thank the Neuro-immersion tech facility of the CRNL (Centre de Recherche en Neurosciences de Lyon) for their material support and the ANR (French National Research Agency) for financing SAMSEI project (ANR-11-IDFI-0034) under the supervision of Pr. X. Martin.

Figure 1 :

 1 Figure 1: (a) Prototype viewed from the side (b) Complete prototype and pneumatic cylinder setup

Figure 3 :

 3 Figure 3: Diagram of the 3D model used to replicate the intervertebral space with force values of each layer from a standard case

Figure 4 :

 4 Figure 4: Needle speed, syr rate and success rate achieved depending on the patient type on separated by experience of the user (med: average, cal: calcified, ove: overweight) with two populations : novices and experts

Table I :

 I Patient types

	Patient	Derma	Epidural space	F dura	LOR
		length (cm) length (mm)	(N)	speed
	Average	1.0	6.0	8.0	average
	Overweight 3.0	4.0	6.0	slow
	Calcified	1.0	4.0	14.0	quick