Kejun Hu
email: kejun.hu@edu.univ-fcomte.fr

Emmanuelle Jacquet
email: emmanuelle.jacquet@univ-fcomte.fr

Najib Kacem
email: najib.kacem@univ-fcomte.fr

Noureddine Bouhaddi
email: noureddine.bouhaddi@univ-fcomte.fr

Jérôme Chambert
email: jerome.chambert@univ-fcomte.fr

Mechanical characterization and identification of a dielectric elastomer

à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

device to control the environmental conditions (temperature and humidity). Tensile test have been performed with an imposed loading-unloading strain rate for various temperatures. The measured forcedisplacement data have been filtered by moving mean. The force-displacement response exhibits an energy storage within the material and a dependence on temperature. Assuming incompressibility, three hyperelastic models have been considered : the Neo-Hookean model with one parameter, the Mooney-Rivlin model with two parameters and the Ogden model with two parameters. Each set of material parameters has been identified by using a non-linear inverse method based on non-linear least squares. It has been observed that, not only the Neo-Hookean model does not fit well with experimental data, , but also the Mooney-Rivlin model is more precise than the Ogden model. Thus, it has been found that the Mooney-Rivlin model is the best compromise between complexity and suitability. The relationship between initial shear modulus and temperature has been approximated by a linear relationship with a sufficiently small residual error.

Keywords : Elastomer -Mechanical characterization -Thermo-mechanical effect -Parameters identification

Introduction

Mechanical structures are subjected to ambient vibrations or acoustic pressures which can be used to produce energy using adequate transducers that convert mechanical stresses into electrical polarization. The energy produced in these cases can be stored and used in non-energy-intensive applications for which the autonomy is a priority, irrespective of the external conditions under which these structures are located. This approach to the production of clean and secure energy (endurable, wireless) now occupies a prominent place in many applications. It is very well adapted to intelligent systems of small sizes such as sensor networks or embedded systems, or mobile phones for various applications for example, surveillance, diagnosis and control in a transportation facility, aeronautics, biomedical, the environment, etc. In addition, the use of intelligent materials in the recovery of energy from the mechanical vibrations of the surrounding environment is a recent and promising concept. For this, electro-active polymers (PEA) are materials that have the ability to change shape and size when a certain voltage or currant intensity is applied to them. The use of these polymers has recently emerged to take advantage of their mechanical behavior and especially the use of dielectric elastomer (DGEs) as an energy harvesting device can dramatically improve transduction performance in storage energy density and non-linear properties (large displacements and large deformations). In this multi-physics and multidisciplinary project, we are interested in the characterizing and the modeling of a material marketed by 3M company under the name of "VHB 4910" (see Ref. [1]) for different temperatures.

Experiments

Specimen preparation

The 3M VHB 4910 is in the form of an adhesive tape with a thickness of 1 mm and a width of 5 cm. The gelatinous aspect of this elastomer leads to some difficulties during the realization of the tensile tests, especially in the cutting of the samples. It is very easy to impose a pre-tension within the materials if we cut it with a cutter which could trigger the non-repeatability of test results. The chosen solution is to cut the sample along its width direction instead of its length direction. The sample dimensions are finally fixed at a length of 40 mm, a width of 50 cm and a thickness of 2 mm (consisting of the superposition of two adhesive tape layers). Figure 1 shows some samples which have been tested.

Experimental set-up

Experimental processing

The tensile tests have been done by the following steps : (1) Prepare the sample carefully. (2) Prepare the input conditions (applied displacement) within the computer software. (3) Put the sample on the load cell, and tighten the crosshead properly with a wrench while keeping the crosshead parallel to avoid twisting. (4) Close the thermo-conditioning device and configure the tested temperature, wait until the environmental conditions stabilize. (5) Initialize the measured force to zero, launch the test, stop when the measured force returns to zero. Then Save the test result. (6) Open the thermo-conditioning device and remove the sample (7) Reinitialize the position of crosshead for the next test.

Hyperelastic modelling

A hyperelastic material is defined by an non-linear elastic behaviour for which the stress-strain relationship derives from a strain-energy density function W . For isotropic material, the strain-energy function depends upon, either the principal invariants (I 1 , I 2 , I 3 ) of the Cauchy-Green tensors, or the principal stretches (λ 1 , λ 2 , λ 3 ) :

I 1 = λ 2 1 + λ 2 2 + λ 2 3 ; I 2 = λ 2 1 λ 2 2 + λ 2 2 λ 2 3 + λ 2 3 λ 2 1 ; I 3 = λ 2 1 λ 2 2 λ 2 3
If the hyperelastic material is also incompressible, then I 3 = 1.

Selected models for incompressible materials

Three strain-energy density functions, namely Neo-Hookean, Mooney-Rivlin and Ogden have been previously used in Ref. [START_REF] Abbad | Développement d'un traitement acoustique basses-fréquences à base de résonateurs d'Helmholtz intégrés à membrane électroactive[END_REF] for describing the mechanical behaviour of VHB 4910 elastomer at room temperature. For a complete review on hyperelastic materials, see Ref. [START_REF] Holzapfel | Nonlinear Solid Mechanics -A Continuum Approach for Engineering[END_REF].

The Neo-Hookean strain-energy potential (hereinafter called 'NH model') is given by :

W NH (I 1 ) = µ 2 (I 1 -3) (1)
where the material constant µ is the initial shear modulus.

The strain-energy function for 2-parameter Mooney-Rivlin model (further referred to as 'MR model') is defined by :

W MR (I 1 , I 2 ) = C 10 (I 1 -3) + C 01 (I 2 -3) (2)
where C 10 and C 01 are material constants. The initial shear modulus is given by µ = 2(C 10 + C 01 ).

The strain-energy function proposed by Ogden has the general form :

W Ogden (λ 1 , λ 2 , λ 3 ) = N i=1 µ i α i (λ α i 1 + λ α i 2 + λ α i 3 -3) (3)
where the material constants N , µ i , α i (i = 1, 2, • • • , N ) are related to the initial shear modulus by :

µ = 1 2 N i=0 µ i α i .
In this paper, the N -parameter is fixed to 1, so that we consider the Ogden model with 2 parameters (µ 1 and α 1 ).

Uniaxial tension tests for incompressible materials

It is assumed that the 3M VHB 4910 polymer has a mechanical behaviour similar to isotopic incompressible hyperelastic materials, which can be described by the NH, MR and Ogden models. The material parameters of these models are identified by uniaxial tension tests.

In case of uniaxial tension, the principal stretch along the loading direction is λ 1 = λ. The incompressibility condition implies that λ 2 = λ 3 = λ -1/2 . The strain principal invariants are defined by :

I 1 = λ 2 + 2 λ -1 ; I 2 = 2 λ + λ -2 ; I 3 = 1
The first Piola-Kirchoff stress tensor P (also denoted engineering or nominal stress tensor) is then given by : P = P e 1 ⊗ e 1 , where P is the unique component of P, the symbol ⊗ stands for the tensor product and e 1 is the unit basis vector along the loading direction.

In the framework of a uniaxial tension mode for isotropic incompressible hyperelastic materials, the stress component P for the successive NH, MR and Ogden models is given by :

P NH = µ λ -λ -2 ; P MR = 2 1 -λ -3 (λC 10 + C 01 ); P Ogden = µ 1 λ α 1 -1 -λ -α 1 2 -1 (4)
4 Results and discussion

Experimental results

Experimental force-displacement curves do not follow the same trajectory during loading and unloading, it is an indicator that shows the existence of energy storage within the material or the existence of structural damping. Figure 2 illustrates experimental curves in terms of the first Piola-Kirchoff stress P (engineering stress) versus stretch ratio λ after filtering by moving mean for various temperatures. First, it is observed that temperature has an influence on the mechanical behaviour of 3M VHB 4910 material. Second, this figure shows that the repeatability of test results is quite good, except for the tests with temperature of 55 • C, for which discrepancies between the two tests occur at low and high stretches.

In the following of this paper, we will consider only the first test for each temperature (see Figure 2) to perform the parameter identification of material behaviour. F 2 -Experimental data in terms of engineering stress versus stretch ratio after filter process.

Parameter identification of hyperelastic models

Three hyper-elastic models defined by Equations (4) have been used for simulating the uniaxial tensile tests. Non-linear least squares method based on Levenberg-Marquardt algorithm has been performed in MATLAB R (function 'lsqnonlin') to minimize the discrepancy between experimental data and model results. This optimization process extracts the material parameters for the best fit between the test data and the considered model.

The identified parameters and residual sums of squares (RSS) are presented on Table 1. The residual sum of squares (RSS) is defined by :

RSS = n i=1 P data i -P model i 2 (5)
where P data i are the experimental data stresses, P model i the modelled data stresses, and n is the number of data points. A perfect fit would yield a small RSS.

Temperature influence on initial shear modulus

The initial shear moduli on Table 1 using MR and Ogden models have been calculated according to equations presented in Section 3.1. Figure 4 shows that the initial shear modulus decreases with increasing temperature whatever the used model. In addition, a relatively good fitting is obtained by using a linear relationship between initial shear modulus and temperature.

Conclusion

The main objective of this paper was to perform a parameters identification of material properties for a dielectric elastomer named VHB 4930 used as a vibration-based acoustic energy harvester. According to a previous research study [START_REF] Abbad | Développement d'un traitement acoustique basses-fréquences à base de résonateurs d'Helmholtz intégrés à membrane électroactive[END_REF], three hyperelastic material models have been considered to describe the mechanical behaviour using strain-energy functions for rubber-like materials.

Mechanical tests have been carried out using tensile testing machine with stable thermo-environment conditions. The test manipulation was also presented. A series of force-displacement responses have been measured and filtered to identify the model parameters according to each model. A comparison of model performance between those models has been realized. The comparison analysis has shown that : (1) The Neo-Hookean model does not fit well the experimental data. [START_REF] Abbad | Développement d'un traitement acoustique basses-fréquences à base de résonateurs d'Helmholtz intégrés à membrane électroactive[END_REF] The Mooney-Rivlin model is more precise than the two-parameter Ogden model except for high temperature. (3) Material damping has been observed which means the energy dispersion should not be neglected in dynamic analysis. By comparing experimental data with theoretical models, it has been noticed that the initial shear modulus is linearly dependent on temperature in the range of studied temperature. In this paper, we have studied three hyperelastic models and the Mooney-Rivlin model has been chosen as the material model to describe the material behavior for energy harvesting applications.

F 1 -

 1 Examples of tested samples.

A

  universal testing machine (MTS Criterion R : Model 45 with 'MTSS test suite' as the output software) is used for characterizing the mechanical behavior of 3M VHB 4910 at different temperatures. All the tensile tests has been performed at a imposed strain rate of ε = 0.125 s -1 . The machine is equipped with a thermo-conditioning device to control the environmental condition which provides a stable temperature and humidity for the test. A temperature series of [25, 35, 45, 55, 65] • C with a humidity rate of 50 % has been chosen as the test environment conditions. Two tensile tests have been carried out for each temperature. It should be pointed out that tensile tests of 3M VHB 4910 have been previously studied at room temperature in Ref.[START_REF] Abbad | Développement d'un traitement acoustique basses-fréquences à base de résonateurs d'Helmholtz intégrés à membrane électroactive[END_REF].

F 3 -

 3 Stress-stretch curves for various hyperelastic models and temperatures after optimization. (a) : T = 25 • C, (b) : T = 35 • C, (c) : T = 45 • C, (d) : T = 55 • C and (e) : T = 65 • C.

Figure 3

 3 Figure 3 shows the parameter identification results in terms of stress-stretch curves with temperatures of 25 • C to 65 • C. The results from Table 1 and Figure 3 show that : (1) Whatever the test temperature, the Neo-Hookean model is not capable to represent the non-linear mechanical behavior of VHB 4910 elastomer. (2) Regardless of test temperature, the Mooney-Rivlin and two-parameter Ogden models can capture the general behaviour of the stress-stretch curves. (3) By considering RSS values, the Mooney-

  of shear modulus versus temperature. Dotted lines are linear regressions with a coefficient of determination R 2 equal to 0.9324 for NH, 0.9802 for MR and 0.9769 for Ogden.

  Rivlin model is more precise than the two-parameter Ogden model except for the temperature of 65 • C.

			NH			MR				Ogden
	T	µ	RSS	C 10	C 01	µ	RSS	µ 1	α 1	µ	RSS
	[ • C] [kPa]	[kPa 2 ] [kPa]	[kPa]	[kPa]	[kPa 2 ] [kPa]	[-]	[kPa]	[kPa 2 ]
	25	47.7653 12300. 9.8883 33.0891 85.9548 29.75	141.9725 1.0208 72.4628 221.6
	35	41.3455 9010.7 9.0705 27.8029 73.7468 27.5	115.5615 1.0742 62.0681 152.5
	45	37.8530 6809.3 8.3316 24.8727 66.4086 11.5	105.1759 1.0683 56.1797 144.5
	55	33.0365 6572.7 5.5927 25.1361 61.4576 15.8	122.4686 0.8500 52.0492 83.5
	65	32.8303 2745.2 9.1798 16.3530 51.0656 110.5	70.5793	1.2598 44.4579 53.0
		T	1 -Results of parameters identification for various models and temperatures.


Aflah Elouneg

Danas Sutula

Marco Sensale

Franz Chouly

Jérôme Chambert

Arnaud Lejeune

Davide Baroli

Paul Hauseux

Stéphane Bordas

Emmanuelle Jacquet

Mechanical parameters identification of keloid and surrounding healthy skin using Digital Image Correlation measurements in vivo

Keywords: keloid, parameter identification, Digital Image Correlation, Finite Element Model Updating

à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Mechanical parameters identification of keloid and surrounding healthy skin using Digital Image Correlation measurements in vivo

A. ELOUNEG a , D. SUTULA a, b , M. SENSALE Corresponding author, emmanuelle.jacquet@univ-fcomte.fr

Résumé :

La peau humaine se comporte comme une membrane élastique précontrainte. La présence de sites anatomiques favorables à l'apparition de certaines tumeurs, une chéloïde dans notre cas, alors que d'autres en sont systématiquement dépourvus atteste de l'importance de l'environnement mécanique du tissu. Par conséquent, une caractérisation de la peau avec tumeur est nécessaire pour comprendre l'expansion de la chéloïde d'un point de

Abstract :

The human skin behaves as an elastic membrane initially prestressed but not uniformly. The presence of anatomical sites favorable to the appearance of some tumors, a keloid in our case, while other sites never develop them attests to the importance of the mechanical environment of the tissue. Thus, a mechanical characterization of the tumored skin is necessary to understand the keloid expansion from a mechanical point of view. Our case study consists in modeling a bi-material structure composed of a keloid skin surrounded by healthy skin located on upper left arm of a young female. From the experimental measurements in vivo, by combining force sensor, displacement sensor and Digital Image Correlation techniques, we perform a mechanical analysis to characterize the mechanical stress fields over the entire area and on the interface 'healthy skin/keloid skin'. Since the mechanical behavior of the tumorous skin is unknown, many physical models can be implemented and assessed very easily inside the specific digital software to fit with the real data. Once a set of mechanical parameters for both the healthy skin and the keloid skin are identified, the stress fields around the keloid are calculated. Next steps consist in determining matching preferential directions in order to define as precisely as possible the specifications of a device for preventing the growth of keloids.

Introduction

Our problem consist in identifying realistic parameters of a bi-material structure by using both force sensor and Digital Image Correlation measurements in vivo. There are many applications where this problem is important. for example, keloid growth beyond original defect margin have been linked to several supposed mechanisms which include skin mechanical properties [START_REF] Ramakrishnan | Study of 1,000 patients with keloids in South India[END_REF], [START_REF] Wang | Increased transcriptional response to mechanical strain in keloid fibroblasts due to increased focal adhesion complex formation[END_REF]. Moreover, given the fact that specific anatomic sites submitted to high forces and usual body movements are known to be pro-keloid sites [START_REF] Ogawa | Keloid and hypertrophic scarring may result from a mechanoreceptor or mechanosensitive nociceptor disorder[END_REF], researchers have been encouraged to focus on mechanical behaviour of such structures. The main target of the paper is to develop the generic methodology used to identify the mechanical properties of the tumored skin. The whole procedure from in vivo experimental test to stress field computation would be described. The corresponding displacement field is measured by Digital Image Correlation on both healthy skin and keloid skin separetely, and the reaction forces by sensors.

Then, an inverse method, known as Finite Element Model Updating (FEMU) based on the inverse least square method, is performed on these data. The completed model with successfully identified parameter will be used finally to calculate the stress field. Hereafter, the experimental device and the uniaxial measurements are described. Then, a numerical modeling is presented and resulting parameters are provided. The originality of the method consists in an objective function built from displacement fields and taking into account the reaction force value necessary to extension test as a constraint. In conclusion, several improvements are explained.

Materials

Uniaxial tensile test

The keloid specimen investigated is a butterlfly-shaped keloid scar located on the upper left arm of a young Caucasian female (figure 1a). Healthy skin reference on the contralateral site is investigated as an homogeneous reference material. Mean values of the thickness of each structure composed of healthy skin or keloid skin have been obtained by high frequency ultrasound. An uniaxial extension test has been conducted with an ultra-light homemade portable extensometer which minimizes the disturbance of measuring area induced by peripheral zones. Two double pads are sticked on the in vivo skin surface. The extension consists in moving one pad in opposite direction to the other one, which remains static, with a controlled servomotor. The reaction force on the moving pad is measured by a strain gauge sensor. A high resolution camera is synchronized to the mechanical test to record the loaded area during the test. For image correlation clarity reasons, a speckle pattern is applied onto the skin within an area of interest (figure 1b).

Experimental data

The mean value of healthy skin thickness is 1.46 + /-0.1mm measured by analyzing ultrasound images. Whereas, the mean value of keloid skin thickness is 2.93 + /-1.26mm. The extension test provides the reactional force values exerted on the measurement pads and the displacement between the pads during all the mechanical tests (Figure 2). The displacement fields of the zone of interest have been identified by digital image correlation. X-displacement field and Y-displacement field for healthy skin are shown on figure [START_REF] Ogawa | Keloid and hypertrophic scarring may result from a mechanoreceptor or mechanosensitive nociceptor disorder[END_REF]. Note that the zone of interest covers all the area located between the two double-pads 3× larger than the one submitted to a quasi-traction test [START_REF] Jacquet | Ultra-light extensometer for the assessment of the mechanical properties of the human skin in vivo[END_REF]. Hence, these two kinds of measurement data, reactional force and displacement fields, are the exprimental inputs of our numerical model that will be used to identify the physical model.

Process scheme

We recall that the main objective of this work is producing somehow a 'black box' which identify the physical parameters of the keloid, by taking in account obviously the behavior of the peripheral healthy skin. As shown in figure (4), the numerical tool adapted of the uniaxial tensile test can deal with many multi-material models since the geometry and the behavior law are part of the inputs, which makes it more flexible so far. The whole procedure of the inverse identification problem would be described. As the corresponding displacement field is measured by Digital Image Correlation and the reactional force by sensors, experimental data need to be filtered to reduce measurement noise. Then, an inverse method, known as Finite Element Model Updating (FEMU) containing, is performed on these data. In addition, the objective cost function (ref the equation) takes account of the difference between numerical and experimental reaction force as a constraint. The most interesting sub-blocks of the FEMU algorithm are the non-linear FEM solver and the inverse solver. By setting an initial guess of the parameters, the displacements over the nodes are computed for each load step k = 1, 2, ..., N E by solving the variational problem below

F(u (k) ; v) = 0 ∀v ∈ V (1)
where V is a suitable function space of admissible virtual displacement that satisfies boundary conditions on u, and F is the variational form of the quasi-static equilibrium

F(u; v) = dΠ(u + v) d =0 (2)
Π is total potential energy, it's given by

Π = ∫ Ω ψ(u)dx - ∫ Ω B • dx - ∫ ∂Ω T • ds (3)
B and T are respectively the body and the traction forces. In our uniaxial tensile test configuration, they are null. While ψ is the elastic strain energy density which characterizes the behavior law of the materials. We need to set it as an input data, because the main issue of fitting with the experimental force-stretch curve is the choice of the behavior equation. On the other hand, by linking the simulated displacement fields to the experimental ones through the optimization formulation, the inverse solver

Constitutive model

Many recent studies showed that the biological soft tissues, human skin for instance, behave as a hyperelastic material. According to Limbert's review ( [START_REF] Limbert | Mathematical and computational modelling of skin biophysics: a review[END_REF]), the mathematical and computational models of the human skin can be classified into three main categories : phenomenological, structural and structurally based phenomenological models. The first one is purely mathematical, the microstructure and multiple phases as well their associated mechanical properties are ignored. If one considers mechanical behaviour only, a phenomenological model is a set of mathematical relations that describe the evolution of stress as a function of strain, not linearly most of time. It is generally always possible to fit such a constitutive law to a set of experimental data. This approach has two major drawbacks : the resulting constitutive parameters often do not have a direct physical interpretation, and the displacement solution obtained by minimizing the strain energy can be not unique. Although, our objective is to determine the stress field whatever the nature of the behavior law that fit to the experimental curve. We choose to implement then a phenomenological model in our program [START_REF] Holzapfel | Nonlinear Solid Mechanics: A Continuum Approach for Engineering[END_REF]. The other classes of constitutive models can be theoretically input into the numerical model, but they contains too many parameters, which may increase horrifically the computation cost. The underlying assumptions for both materials, healthy skin and keloid, should be the same. In other words, they present the same behavior law but different parameters. We assume that theses soft tissues are hyperelastic, homogenoeus and isotropic. Their large deformations are assumed to be quasi-static.

Inverse identification

The objective function used to identify the material parameters θ = {θ healthy skin ; θ keloid } is the sum over all loading time steps of the absolute difference between the predicted and the experimental displacements, by taking account of the difference between the predicted forces, calculated from the simulated displacements, and the ones gathered from the test. An additional parameter Lagrange multiplier links the two mismatches. The N θ + 1 parameters are calculated by minimizing the following term

J(θ, λ) = N E k=1 ∫ Γ u msr u (k) (θ) -u (k) msr 2 dx + λ ∫ Γ f msr f (u (k) (θ), θ) -f (k) msr dx (4)
with Γ u msr a part of the domain Γ E where the displacement field is measured (i.e. using DIC) and Γ f msr the boundary where the force is measured using deformation gauge. Among many non-linear least squares methods detailed in Harb's Thesis [START_REF] Harb | Identification inverse de paramètres biomécaniques en hyperélasticité anisotrope[END_REF], we have selected Newton-Gauss algorithm to be implemented into the inverse solver. The latter method is enough accurate and not so expensive providing that we set a valid initial guess of the material parameters.

Validation

The entire inverse problem numerical tool is validated in two steps. First of all, we validate separately the FEM non-linear solver by comparing the nodal solutions vis-à-vis commercial FEM codes, COMSOL MULTIPYSICS and ANSYS, and the inverse solver by using the simulated data as the experimental inputs (The material parameters values used for making data are computed through the inverse solver factory). By fixing randomly values too much far from the targeted parameters (up to ×500), the inverse solver was able to converge to the right set of parameters with reasonable number of iterations. Secondly, the computations are performed on real data. We have chosen the alternative compressible hyperelastic Neo-Hookean model ( 5) as an input behavior law to make it simpler for validation. Moreover, the deformed structure shape deriving from the response of the Neo-Hookean model is so realistic as shown in figure [START_REF] Limbert | Mathematical and computational modelling of skin biophysics: a review[END_REF]. Consequently, we have identified 4 material parameters in addition to Lagrange multiplier θ = {µ healthy skin , λ healthy skin , µ keloid , λ keloid }. Further behavior laws should be tested very easily such that the model is in good agreement with available observations.

ψ = µ 2 (I 1 -3 -2 ln J) + λ 2 (ln J) 2 (5)
Figure 5: Predicted displacement fields of the Neo-Hookean alternative model

Postprocessing

At the end of the identification process, the parameters are again introduced into the FEM solver to simulate the model responses as close as possible to the real data. Two final operations take place then: quantifying the model parameters sensitivities with respect to the measurement noise and calculating the stress fields. The output stresses, in the form of Cauchy or Piola-Kirchhoff I (3), will help us to understand better at the first level the behavior in the interface medium healthy-skin/keloid.

Conclusions

Our methodology, concertized in a sophisticated open-source computational framework, is able to provide the mechanical properties, true stress field for instance, of a keloid scar. The latter is growing at the expense of the surrounding healthy skin. Thus, we developed a bi-materials model inside the program to observe the behavior on the interface zone between the keloid and the healthy skin using the FEM Updating algorithms. The inputs to the program are the experimental data, gathered from a uniaxial test and DIC measurement, and the constitutive model as a strain energy density. We have validated and analyzed the reliability of the outputs using faked and real data. In perspectives, several phenomenological material model will be tested to provide results fitting very well with the experimental inputs. And in the long term, our valid results will help to define specifications of a device for preventing the growth of keloids.

Figure 1 :

 1 Figure 1: Uniaxial tensile testing on keloid-skin

Figure 2 :

 2 Figure 2: Force-Displacement curve during extension test

Figure 3 :

 3 Figure 3: Experimental data of keloid-skin undergoing large deformation

Figure 4 :

 4 Figure 4: Process scheme of the inverse problem

  . Univ. Bourgogne Franche Comté, FEMTO-ST institute, UFC/CNRS/ENSMM/UTBM, Departement of Applied Mechanics, Besançon, France. b. Technical University of Liberec, Faculty of Tissue Engineering, Liberec, Czech Republic. c. University of Luxembourg, Institute of Computational Engineering, Luxembourg. d. Univ. Bourgogne Franche Comté, Institut de Mathématiques de Bourgogne, Dijon, France.

	CHAMBERT a , A. LEJEUNE a , D. BAROLI c , P. HAUSEUX c , S. BORDAS c ,
	E. JACQUET a, ,

a, c , F. CHOULY d , J. a


Danas Sutula

Aflah Elouneg

Marco Sensale

Franz Chouly

Jérôme Chambert

Arnaud Lejeune

Davide Baroli

Paul Hauseux

Stéphane Bordas

Emmanuelle Jacquet

Parameter identification problem in bimaterial human skin and sensitivity analysis : Uncertainties in biomechanics of skin

niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Corresponding author, emmanuelle.jacquet@univ-fcomte.fr

Résumé :

Le travail proposé concerne la prédiction de la réponse numérique d'une structure biomécanique soumise à un état de chargement externe inconnu. La méthodologie s'appuie sur la modélisation de structures homogènes puis hétérogènes, telles que des tissus cutanés sains ou pathologiques accessibles expérimentalement et soumis à des conditions aux limites partiellement connues. Les

Abstract :

The proposed paper concerns the prediction of the numerical response of a biomechanical structure submitted to an unknown external loading state.

Key words : Mechanical behaviour ; Human skin, Inverse method 1 Introduction

Biomechanics nowdays largely covers many fields of mechanics and combines the various complexities related to.

Our focus is on the patient-specific biomechanical characterisation of human skin with a keloid scar (Fig. 1). Keloids are non-cancerous tumours [START_REF] Ogawa | Mechanobiology of cutaneous wound healing and scarring[END_REF] that grow continuously on the skin. The evolution of keloids is known to be related to genetic, biological, biophysical, and biomechanical factors [START_REF] Chambert | Numerical analysis of keloid scar in the presternal area[END_REF]. In our work, we focus on the latter type of influence ; for example, the state of stress inside a keloid and in the surrounding skin is known to play a significant role in the growth of the keloid [START_REF] Chambert | Numerical analysis of keloid scar in the presternal area[END_REF].

Materials and Methods

Material problem

To characterize the mechanical behaviour of a bi-material structure composed of healthy skin and keloid skin, a suitable material model has been chosen. In vivo experiments are carried out using a custommade extensometer [START_REF] Jacquet | Ultra-light extensometer for the assessment of the mechanical properties of the human skin in vivo[END_REF] to obtain N E data points : k = 1, 2, . . . , N E , of the applied force f k and the resulting surface displacement field ūk . The surface displacements are determined using a digital image correlation (DIC) technique (Fig. 2). As shown in figure [START_REF] Jacquet | Ultra-light extensometer for the assessment of the mechanical properties of the human skin in vivo[END_REF], the DIC measurement domain limits are far beyond the limited zone where the inverse identification is performed. This would help to reduce noises. The numerical model consists of a finite element discretisation of the part of the skin with the keloid. A keloid is observed on a patient using ultrasound imaging and optical microscopy. The 3D image of the surface is used to develop the keloid-skin geometrical model (Fig. 3). The model domain Ω is divided into two sub-domains

Ω 1 ∪ Ω 2 = Ω.

Optimization problem

For the chosen hyperelastic material model(s) with a total of N m model parameters we look for the parameters m that minimise the mismatch between the experimental u (k) msr and the numerical u (k) observations of the displacements on the surface of the skin Γ E subject to quasi-static equilibrium F (u (k) , m; v) = 0 ∀v ∈ V on all loading steps k = 1, 2, . . . , N E . V is the space of admissible virtual displacements.

In vivo extensometer experiments are controlled with prescribed displacement of the pad [START_REF] Jacquet | Ultra-light extensometer for the assessment of the mechanical properties of the human skin in vivo[END_REF]. In this case, the measured force f (k) msr needs to be introduced into the cost functional in a form of a constraint equation. This can be done using the Lagrange multiplier method. By introducing a single multiplier λ the reaction force constraint can be enforced on average over all times. The problem is then to find the stationary points of J(m, λ) = N E k=1 J (k) (m, λ). The optimization problem is defined as :

m = argmin m,λ∈R Nm+1 N E k=1 J (k) (m, λ) (1)
subject to :

F (u (k) , m; v) = 0 ∀v ∈ V (2)
where

J (k) (m, λ) = J(u (k) (m), m, λ) := Γ u msr u (k) (m) -u (k) msr 2 dx + λ Γ f msr f (u (k) (m), m) -f (k) msr dx (3)
with Γ u msr a part of the domain Γ E where the displacement field is measured (i.e. using Digital Image Correlation) and Γ f msr the boundary where the force is measured using deformation gauge. the problem defined by equations ( 1) and ( 2) can be solved using Newton's method.

D mi J (k) = ∂ u J(u (k) , m)[D mi u (k) ] + ∂ mi J(u (k) , m) = 0 ∀i = 1, 2, . . . , N m (4)
Similarly, the Hessian of the model cost can be expressed as :

D 2 mi mj J (k) = ∂ u J[D 2 mi mj u (k) ] + ∂ 2 uu J[D mi u (k) , D mj u (k) ] + ∂ mi u J[D mj u (k) ] + ∂ 2 u mj J[D mi u (k) ] + ∂ 2 mi mj J (5)
Usually, the derivatives D mi u and D 2 mi mj u are not computed explicitly ; instead, the adjoint solution method is involved thereby enabling one to compute D mi J (k) and D 2 mi mj J (k) directly and efficiently.

It is also possible to perform the optimization at each time and determine the locally optimal model parameters. This can be of interest since the local model parameter deviations from the globally optimal values can serve as a good indication of the fitness of a particular hyperelastic model to approximate the biological tissue at hand.

Model parameters sensitivities

Goal of this study is to determine model parameters m minimizing the cost functional J(m) depending on measurements obtained from a DIC analysis. As the precision of the measurements can not be trusted in general, we seek to study the sensitivity of the inferred model parameters to measurement noise. To this end, an artificial noise function can be superposed on the true measurements and the sensitivity of the model parameters to the noise function can be analyzed using automatic differentiation tools within the computational framework. The model parameter sensitivities δm with respect to a noise function v

(k)
msr at pseudo-time k can be computed as :

∂m ∂v (k) msr = - d 2 J dm 2 -1 d 2 J (k) ∂v (k) msr dm (6)
The noise function can be one for the measured displacement field or one for the measured reaction force.

In the case of the displacement field uncertainty, the model parameter sensitivities can be projected on a scalar function space to yield plots indicating where a particular model parameter is most sensitive to the measured displacements. For example, Fig. 4 shows the sensitivity of the generalized Poisson's parameter in the Neo-Hookean hyper-elastic material model. Subsequently, the suitability of a chosen material model can been assessed by quantifying the sensitivity of the model parameters with respect to the noise in the measurement data. Consequently, equation 6 can be used to assess the fitness of a particular material model for the skin-keloid tissue considering the noise in the experimental measurement data.

Application

Figure 5 describes the numerical process from experimental data to outputs like model parameters for skin and surroundig healthy skin and true stress fields. This overall process has been implemented using python packages for data preprocessing (filtering, projection on mesh) and FEniCs tool for the Finite Element simulation. Identification and sensitivity analysis has been implemented in python. ), the second is fixed (u x = 0). External boundaries are such that top and bottom are fixed and right and left are constrained to u y = 0. All the blue domain is meshed but only the rectangular area with apparent elements is used for identification process. After parameters identification from minimization of the cost function, a comparison between raw data reaction force and numerical reaction force calculated has exhibited a satisfactory correlation for the beginning of the extension (deformation of 8% or less) and a lack of correlation for higher deformations.

Discussion and Conclusion

Inverse least square method has been implemented for the uniaxial extension test of a bi-material structure and with simple hyperelastic model. More complex material models able to take into account incompressibility constraint, like Higher order Mooney-Rivlin model should bring more accurate results. The procedure consists in computing sensitivities of model parameters with respect to measurement noise of different models. Furthermore, the sensitivity analysis can be useful in the design of the experiments in order to maximize the information gained from the measurements. Note that there are some problems with high level of noise into experimental data. That's why it seems to be necessary to filter data before inverse analysis. After improvements done to the numerical model, the results can be used for predicting stress fields around the keloid scar and bring the mechanical basis of predicting the preferential directions of its growth. The ultimate goal consists in defining the specifications of a predictive device. Hereafter this concrete example, the focus of the paper has been on the analysis processing able to solve a large variety of problems from raw experimental data to robust mechanical simulations in which measurements that can be obtained are highly variable and several boundary conditions are totally unknown.

  a , A. ELOUNEG a , M. SENSALE a,b , F. CHOULY c , J. CHAMBERT a , A. LEJEUNE a , D. BAROLI b , P. HAUSEUX b , S. BORDAS b , E. JACQUET a, a. Univ. Bourgogne Franche-Comté, FEMTO-ST Institute, Department of Applied Mechanics, Besançon, France. b. University of Luxembourg, Institute of Computational Engineering, Luxembourg. c. Univ. Bourgogne Franche-Comté, Institut de Mathématiques de Bourgogne, Dijon, France.

F 1 -

 1 Keloid scar located on the left-upper arm (X = 15 mm, Y = 47 mm).

F 2 -

 2 Reference displacement fields obtained by digital image correlation and the post-processed strain fields. Fields for healthy skin

F 3 -

 3 Geometry model and boundaries conditions of keloid/skin characterisationFor the sake of brevity suppose the material model parameters and any constraint multipliers can be collectively denoted by mi for i = 1, 2, ..., N m where in the case of a single constraint multiplier N m = N m + 1. In this case the gradient of the cost at load-step k = 1, 2, . . . , N E can be expressed as :

F 4 -

 4 Poisson's parameter sensitivity due to the measured displacement field at the final time-step.

F 5 -

 5 Numerical processFor the estimation of the model parameters, an inverse solution strategy to best fit the model to a set of experimental measurements of the skin-keloid mechanical response due to different loading conditions has been applied. Finite element model of extensometer experiment is illustrated in Figure6. Extension is controlled by the two pads (grey squares). They are modelised like rigid holes. One is uniaxially translated in the longitudinal direction (u x = u prescribed x

F 6 -

 6 Model with external boundaries : fixed top and bottom and u y = 0 Subsequently, the suitability of a chosen material model has been assessed by quantifying the sensitivity of the model parameters with respect to the noise in the measurement data (Fig.7). F 7 -Objective functional derivative values, DJ/Dm, at different observation times, t