D Remache
email: djamel.remache@univ-amu.fr

M Semaan
email: marie.semaan@balamand.edu.lb

J M Rossi

M Pithioux
email: martine.pithioux@univ-amu.fr

J L Milan

Application of the Johnson-Cook plasticity model in the Finite Element simulations of the nanoindentation of the cortical bone

Keywords: Mechanical cortical bone behavior, Nanoindentation test, Inverse optimization approach, Design of experiments (DOE), Johnson-Cook model

, a plastic hardening until failure was observed on the micro-scale of a dry ovine osteonal bone samples subjected to micropillar compression. The purpose of the current study was to apply an isotropic hardening model in the finite element simulations of the nanoindentation of the cortical bone to predict its mechanical behavior. The Johnson-Cook (JC) model was chosen as the constitutive model. The finite element modeling in combination with numerical optimization was used to identify the unknown material constants and then the finite element solutions were compared to the experimental results. A good agreement of the numerical curves with the target loading curves was found and no pile-up was predicted. A Design Of Experiments (DOE) approach was performed to evaluate the linear effects of the material constants on the mechanical response of the material. The strain hardening modulus and the strain hardening exponent were the most influential parameters. While a positive effect was noticed with the Young's modulus, the initial yield stress and the strain hardening modulus, an opposite 1 effect was found with the Poisson's ratio and the strain hardening exponent. Finally, the JC model showed a good capability to describe the elastoplastic behavior of the cortical bone.

Introduction

The instrumented indentation technique is usually performed for local mechanical characterization of the cortical bone [START_REF] Carnelli | A finite element model for direction-dependent mechanical response to nanoindentation of cortical bone allowing for anisotropic post-yield behavior of the tissue[END_REF][START_REF] Carnelli | Nanoindentation testing and finite element simulations of cortical bone allowing for anisotropic elastic and inelastic mechanical response[END_REF][START_REF] Faingold | Nanoindentation of osteonal bone lamellae[END_REF][START_REF] Fan | Effects of viscoelasticity and time-dependent plasticity on nanoindentation measurements of human cortical bone[END_REF][START_REF] Fan | Anisotropic properties of human tibial cortical bone as measured by nanoindentation[END_REF][START_REF] Feng | Mechanical properties of porcine femoral cortical bone measured by nanoindentation[END_REF][START_REF] Franzoso | Elastic anisotropy of human cortical bone secondary osteons measured by nanoindentation[END_REF][START_REF] Isaksson | Precision of nanoindentation protocols for measurement of viscoelasticity in cortical and trabecular bone[END_REF][START_REF] Joo | Mechanical properties on nanoindentation measurements of osteonic lamellae in a human cortical bone[END_REF][START_REF] Rho | Elastic properties of microstructural components of human bone tissue as measured by nanoindentation[END_REF][START_REF] Rho | Elastic properties of human cortical and trabecular lamellar bone measured by nanoindentation[END_REF][START_REF] Spiesz | Elastic anisotropy of uniaxial mineralized collagen fibers measured using two-directional indentation. effects of hydration state and indentation depth[END_REF][START_REF] Vanleene | Influence of nanoindentation test direction on the elastic properties of human cortical bone lamellae[END_REF][START_REF] Wang | Elastic modulus and hardness of cortical and trabecular bovine bone measured by nanoindentation[END_REF]. The well documented Oliver-Pharr (OP) analytical approach [START_REF] Oliver | An improved technique for determining hardness and elastic modulus using load and displacement sensing indentation experiments[END_REF][START_REF] Oliver | Measurement of hardness and elastic modulus by instrumented indentation: Advances in understanding and refinements to methodology[END_REF], available in the software included with most commercial Indentation Testers, has widely been adopted for measuring the elastic modulus and the hardness of materials from indentation load-displacement data obtained during one cycle of loading and unloading. This method has been introduced according to assumptions that the material is homogeneous, linearly elastic, isotropic and incompressible, and during unloading, it is assumed that only the elastic displacements are recovered, i.e., the unloading response is considered as purely elastic [START_REF] Chen | Nanoindentation of hyperelastic polymer layers at finite deformation and parameter re-identification[END_REF][START_REF] Oliver | An improved technique for determining hardness and elastic modulus using load and displacement sensing indentation experiments[END_REF][START_REF] Poon | An analysis of nanoindentation in linearly elastic solids[END_REF]. This method is therefore limited with regard to determining the mechanical properties of the heterogeneous and multi-scale structural materials such as the cortical bone tissue [START_REF] Ghanbari | Nonlinear hierarchical multiscale modeling of cortical bone considering its nanoscale microstructure[END_REF][START_REF] Hamed | Multiscale modeling of elastic properties of cortical bone[END_REF]. Indeed, in addition to the inhomogeneities of the material, such as porosity, defects, grain boundaries, etc..., viscoplastic deformations are added to the elastic deformation [START_REF] Bernard | Accurate measurement of cortical bone elasticity tensor with resonant ultrasound spectroscopy[END_REF][START_REF] Carnelli | A finite element model for direction-dependent mechanical response to nanoindentation of cortical bone allowing for anisotropic post-yield behavior of the tissue[END_REF][START_REF] Carnelli | Nanoindentation testing and finite element simulations of cortical bone allowing for anisotropic elastic and inelastic mechanical response[END_REF][START_REF] Fan | Effects of viscoelasticity and time-dependent plasticity on nanoindentation measurements of human cortical bone[END_REF][START_REF] Fan | Anisotropic properties of human tibial cortical bone as measured by nanoindentation[END_REF][START_REF] Fondrk | Some viscoplastic characteristics of bovine and human cortical bone[END_REF][START_REF] Hoffler | An application of nanoindentation technique to measure bone tissue lamellae properties[END_REF][START_REF] Lang | Elastic coefficients of animal bone[END_REF][START_REF] Peterson | Biomechanics: principles and applications[END_REF][START_REF] Swadener | Effects of anisotropy on elastic moduli measured by nanoindentation in human tibial cortical bone[END_REF][START_REF] Van Buskirk | The elastic moduli of bone[END_REF][START_REF] Zysset | An alternative model for anisotropic elasticity based on fabric tensors[END_REF], making the strain field induced under the indenter very heterogeneous and complex. Furthermore, the mechanical properties of the cortical bone depend also on other factors such as the anatomical zone of the tested sample, sex, age, the hydration state of the material, the probe geometry, etc... [START_REF] Feng | Mechanical properties of porcine femoral cortical bone measured by nanoindentation[END_REF][START_REF] Hoffler | An application of nanoindentation technique to measure bone tissue lamellae properties[END_REF][START_REF] Öchsner | Biomechanics of hard tissues[END_REF][START_REF] Peterson | Biomechanics: principles and applications[END_REF][START_REF] Rodriguez-Florez | Insight into differences in nanoindentation properties of bone[END_REF]. The finite element approach has thus become essential for determining the mechanical properties of materials such as the cortical bone [START_REF] Bolshakov | Influences of stress on the measurement of mechanical properties using nanoindentation: Part ii. finite element simulations[END_REF][START_REF] Carnelli | A finite element model for direction-dependent mechanical response to nanoindentation of cortical bone allowing for anisotropic post-yield behavior of the tissue[END_REF][START_REF] Fan | Three-dimensional finite element analysis of the effects of anisotropy on bone mechanical properties measured by nanoindentation[END_REF], particularly their time-dependent mechanical properties such as creep, viscoelasticity or elastoplasticity [START_REF] Carnelli | A finite element model for direction-dependent mechanical response to nanoindentation of cortical bone allowing for anisotropic post-yield behavior of the tissue[END_REF][START_REF] Chang | Practical assessment of nanoscale indentation techniques for the biomechanical properties of biological materials[END_REF].

One significant problem in the indentation modeling of the elastic-plastic materials is the pile-up of material around the contact impression, due to a plastic flow of material from beneath the indenter toward the sample surface. This phenomenon is frequently observed in the indentation of the elastic-perfectly-plastic materials such as metal materials [START_REF] Fan | Three-dimensional finite element analysis of the effects of anisotropy on bone mechanical properties measured by nanoindentation[END_REF], but it doesn't occur in the indentation of work-hardening materials [START_REF] Bolshakov | Influences of pileup on the measurement of mechanical properties by load and depth sensing indentation techniques[END_REF]. Based on Atomic Force Microscopy [START_REF] Fan | Three-dimensional finite element analysis of the effects of anisotropy on bone mechanical properties measured by nanoindentation[END_REF][START_REF] Mullins | Calibration of a constitutive model for the post-yield behaviour of cortical bone[END_REF][START_REF] Schwiedrzik | Quantitative analysis of imprint shape and its relation to mechanical properties measured by microindentation in bone[END_REF] and Scanning Electron Microscopy [START_REF] Johnson | Characterization of the microindentation of bone via scanning electron microscopy[END_REF] observations and based on a statistical shape modeling of the indentation residual imprint [START_REF] Schwiedrzik | Quantitative analysis of imprint shape and its relation to mechanical properties measured by microindentation in bone[END_REF],

it was observed that pile-up doesn't occur in cortical bone tissue during the nanoindentation testing, which suggests that the material could be considered as a work-hardening material. It would then be appropriate to use a plastic hardening model for the finite element modeling of the nanoindentation of the cortical bone tissue. Indeed, the classical metal plasticity model predicted a significant pile-up of material and an overestimation of the maximum depth by approximately 30% [START_REF] Mullins | Calibration of a constitutive model for the post-yield behaviour of cortical bone[END_REF], however, it has been suggested by [START_REF] Fan | Three-dimensional finite element analysis of the effects of anisotropy on bone mechanical properties measured by nanoindentation[END_REF] that it might be possible to eliminate the pile-up effect by introducing a strain hardening rate to the classical metal plasticity model [START_REF] Bolshakov | Influences of pileup on the measurement of mechanical properties by load and depth sensing indentation techniques[END_REF][START_REF] Fan | Three-dimensional finite element analysis of the effects of anisotropy on bone mechanical properties measured by nanoindentation[END_REF][START_REF] Mullins | Calibration of a constitutive model for the post-yield behaviour of cortical bone[END_REF]. The finite element results of this study have not been validated by comparison with the experimental results. In a recent study [START_REF] Schwiedrzik | In situ micropillar compression reveals superior strength and ductility but an absence of damage in lamellar bone[END_REF], a dry ovine osteonal bone was subjected to micropillar compression, microindentation and macroscopic compression tests to identify the elastic modulus, yield stress, plastic deformation, damage accumulation and failure mechanisms. The results showed that, on the micro-scale a highly ductile behavior with continuous hardening and no damage until failure were exhibited by the samples, whereas, on the macro-scale the samples failed mostly in a quasi-brittle fashion due to the existing of cracks, pores and defects leading to a reduction of strength and ductility. The measured strength on the micro-scale was found to be higher than on the macro-scale. It was assumed that on the micro-scale the indentation results represent the local bulk material properties. In a more recent study [START_REF] Ramezanzadehkoldeh | Nanoindentation response of cortical bone: dependency of subsurface voids[END_REF], nanoindentation was performed on mouse femur and then Focused Ion Beam-Scanning Electron Microscopy (FIB-SEM) was used to map microporosities below the nanoindentation sites. Several cavities such as lacunae were identified, but for the majority of indentation sites, cavity-free material was found, and this corresponded to a reasonable consistency in the maximum indentation depth and Young's modulus. A significantly softer response was measured in the sites with subsurface cavities. Further investigations [START_REF] Ramezanzadehkoldeh | Nanoindentation response of cortical bone: dependency of subsurface voids[END_REF][START_REF] Schwiedrzik | A generalized anisotropic quadric yield criterion and its application to bone tissue at multiple length scales[END_REF][START_REF] Schwiedrzik | Quantitative analysis of imprint shape and its relation to mechanical properties measured by microindentation in bone[END_REF], based on the porous nature of the cortical bone, assumed that a model, which accounts for frictional mechanisms, is required for accurate material behavior prediction [START_REF] Carnelli | A finite element model for direction-dependent mechanical response to nanoindentation of cortical bone allowing for anisotropic post-yield behavior of the tissue[END_REF][START_REF] Carnelli | Nanoindentation testing and finite element simulations of cortical bone allowing for anisotropic elastic and inelastic mechanical response[END_REF][START_REF] Mullins | Calibration of a constitutive model for the post-yield behaviour of cortical bone[END_REF][START_REF] Tai | Nanogranular origins of the strength of bone[END_REF]. The extended Drucker-Prager plasticity model [START_REF] Drucker | Soil mechanics and plastic analysis or limit design[END_REF], intended for geological materials that exhibit pressure-dependent yield, has accurately predicted the post-yield behavior of a cortical bone and no pile-up in the material [START_REF] Mullins | Calibration of a constitutive model for the post-yield behaviour of cortical bone[END_REF].

However, no evidence of fracture events was observed in the experimental testing, as it has been expected.

In the current study, nanoindentation tests were performed on a human femur cortical bone. Thanks to the microscope included with the used commercial Indentation Tester, visible surface defects, holes, and pores were avoided. Hence, it was assumed that the nanoindentation results represent the local bulk material properties. The cortical bone was considered as a work-hardening elastic-plastic material. The Johnson-Cook (JC) isotropic hardening model [START_REF] Johnson | A constitutive model and data for metals subjected to large strains, high strain rates and high temperatures[END_REF] was chosen as the constitutive model. Unlike the most previous studies, we were interested in the elastic-plastic response during the loading phase of a cortical bone sample subjected to nanoindentation. To minimize errors that might result from underestimating of the contact area, the elastic constants, corresponding to the initial elastic contact, and the JC material constants were identified from the non-linear curve-fit of the force-displacement raw data using an inverse analysis. The elastic constants were predicted independently and simultaneously unlike in the OP method, which requires a priori knowledge of one constant to compute the second one. The finite element method in combination with numerical optimization were used to identify the unknown material parameters. This indentation data analysis method has already been used to investigate some classical materials such as polymer layers [START_REF] Chen | Nanoindentation of hyperelastic polymer layers at finite deformation and parameter re-identification[END_REF][START_REF] Chen | Parameter re-identification in nanoindentation problems of viscoelastic polymer layers: small deformation[END_REF] and steel [START_REF] Kang | Determination of elastic and viscoplastic material properties obtained from indentation tests using a combined finite element analysis and optimization approach[END_REF].

Materials and Methods

Experimental tests

A cortical bone sample was obtained from the femur of a human female cadaver (82 years old) without bone disease. The sample was cut into a parallelepiped shape of a few millimeters by a water-cooled low-speed diamond saw (Buehler Isomet 4000, Buehler, Lake Bluff, IL, USA) and was stored at -20°C. No chemicals were used for cleaning and dehydration in order to preserve the mechanical properties of the material microstructural features [START_REF] Reisinger | Principal stiffness orientation and degree of anisotropy of human osteons based on nanoindentation in three distinct planes[END_REF]. Also, resin was avoided since its infiltration will alter the measured elastic properties [START_REF] Bushby | Nanoindentation of bone: Comparison of specimens tested in liquid and embedded in polymethylmethacrylate[END_REF]. The sample sides were polished using silicon carbide abrasive papers of decreasing grit size (400, 800 and 1200 grit) under deionized water [START_REF] Rho | Elastic properties of human cortical and trabecular lamellar bone measured by nanoindentation[END_REF][START_REF] Rho | Elastic properties of microstructural components of human bone tissue as measured by nanoindentation[END_REF] then, with a diamond suspension (1 µm). The sample was cleaned of polishing dust with water.

Nanoindentation tests were conducted at room temperature (∼23°C) using an indentation Tester of Anton Paar (Tester NHT 2 , Switzerland and Austria) with the modified Berkovich diamond tip. The tip area function was calibrated by a fused silica reference sample. The tests were performed in the load control mode in the two principal bone directions, i.e., the axial direction (corresponding to the long bone axis) and the transverse direction (perpendicular to the long bone axis). In the axial direction, the tests were performed at two locations; osteons and interstitial lamellae observed by optical microscopy (Fig. 1). In the transverse direction, it was difficult to distinguish between osteons and interstitial lamellae. The nanoindentation tests were conducted up to a maximum load of 20 mN, at a loading/unloading rate of 40 mN/min. A total of 26 indents were made: 8 in osteons locations, 8 in interstitial lamellae locations, and 10 in the transversal side.

The elastic modulus and the hardness constant of the cortical bone sample were measured by the OP technique [START_REF] Oliver | An improved technique for determining hardness and elastic modulus using load and displacement sensing indentation experiments[END_REF][START_REF] Oliver | Measurement of hardness and elastic modulus by instrumented indentation: Advances in understanding and refinements to methodology[END_REF], available in the software included with the commercial Indentation Tester. The elastic modulus is given by:

E = 1 -ν 2 1 Er - 1-ν 2 i Ei (1)
where E i is Young's modulus of the indenter (E i =1141 GPa), ν i (ν i =0.07) and ν are Poisson's ratio of the indenter and the material, respectively, and E r is the reduced elastic modulus given by:

E r = S. √ π 2.β. A (h c ) (2)
where S = dF/dh is the contact stiffness (defined as the slope of the initial part of the unloading curve during the first stage of unloading), A is the projected area of the contact impression, β is a constant (a correction factor) which depends on the indenter geometry [START_REF] Rho | Elastic properties of human cortical and trabecular lamellar bone measured by nanoindentation[END_REF] (for the Berkovich indenter β = 1.034), h c is the contact depth of the indenter with the sample at the maximum indentation force. A value of 0.3 was chosen for ν.

The hardness, H, is given by the following relationship:

H = P max A (3)
where P max is the maximum indentation force.

Finite Element analysis 2.2.1. Finite element model

The numerical simulation of the Berkovich indentation experiments was modeled in two-dimensional axisymmetric finite element model using the commercial finite element code Abaqus ® 6.14-2. A rectangular region of dimensions 40×40 µm was assumed as the problem domain. The Berkovich indenter was defined as an analytical rigid surface part with an angle of 65.27°between the vertical axis and the indenter face, and the cortical bone was assumed to be a deformable body. Asymmetric linear triangular elements of type CAX3 were used for meshing the sample geometry (total number of elements: 5202 and total number of nodes: 2704). After mesh convergence study, finer elements were used in the vicinity of the indenter surface and coarse elements farther away. A reference point was defined on the indenter geometry, on which a total indentation load of 20 mN was applied under load-control conditions. The geometry and boundary conditions of the finite element model are illustrated in Fig. 2. The indenter was defined as master surface, while the nodes at the top of the sample part were defined as slave surface. A surface-to-surface (standard) contact using the finite-sliding formulation was defined as a type of interaction, and to avoid the penetration of the slave surface into the master surface at the constraint locations, the "hard" contact was chosen as a contact property.

A rigid body constraint was applied to the indenter and a frictionless hard contact was assumed. The elastic-plastic behavior of the material was governed by the Johnson-Cook behavior law presented in the following section.

Johnson-Cook constitutive model

The nanoindentation tests were simulated using the isotropic hardening model proposed by Johnson-Cook (JC) [START_REF] Johnson | A constitutive model and data for metals subjected to large strains, high strain rates and high temperatures[END_REF]. The JC model is a phenomenological model for materials subjected to large strains, high strain rates and high temperatures. The flow stress of this particular type of Mises plasticity model is expressed as:

σ = A + Bε n p 1 + C ln ε ε0 [1 -T * m ] (4)
T * = T -T 0 T m -T 0
where A , B, n, C, and m are the material constants; ε p is the equivalent plastic strain, ε and ε0 are the plastic strain rate and the reference plastic strain rate (i.e., which corresponds to the quasi-static test), respectively. T * is the normalized temperature; T 0 is a reference (or transition) temperature (i.e., the room temperature in a quasi-static test), and T m is the melting temperature. A is defined as the initial static yield stress (i.e., σ = A when ε p = 0 , ε= ε0 and T = T 0 ); B is the strain hardening modulus and n is the strain hardening exponent.

The expression in the first set of brackets gives the stress as a function of strain for ε = ε0 and T * = 0. The expressions in the second and third sets of brackets represent the effects of elevated strain rates and temperature, respectively. Since the experimental indentation was a quasi-static test without significant thermal changes, the second and third factor in the constitutive relationship were neglected. Therefore, we had to determine the constants in the first set of brackets. Finally, the model was used in conjunction with a linear elastic material model. The choice of this model was based on the fact that it is able to describe the strain hardening effect in materials and because of its good prediction and simplicity [START_REF] Dey | On the influence of constitutive relation in projectile impact of steel plates[END_REF]. This model was used in bone cutting simulations [START_REF] Alam | Finite element analysis of forces of plane cutting of cortical bone[END_REF][START_REF] Alam | Thermal analysis of orthogonal cutting of cortical bone using finite element simulations[END_REF][START_REF] Alam | 3d finite-element modelling of drilling cortical bone: Temperature analysis[END_REF][START_REF] Zhang | On the selection of johnson-cook constitutive model parameters for ti-6al-4v using three types of numerical models of orthogonal cutting[END_REF] and in simulations of impact and penetration related problems [START_REF] Banerjee | Determination of johnson cook material and failure model constants and numerical modelling of charpy impact test of armour steel[END_REF].

Calibration procedure

Numerical optimization is a mathematical formulation that allows to minimize or maximize a particular objective function subjected to constraints on its variables [44].

The objective function depends on certain variables called unknowns or parameters (frequency, stress, displacement, material parameter...). So, in the optimization process, the objective, the variables, and the constraint functions are to be defined. In this work, our goal was to determine the cortical bone constants so that the numerical forcedisplacement curve matches the experimental one. The force-displacement raw data were used in the inverse analysis as input data in order to minimize errors that might result from underestimating of the contact area. It was assumed that the measurement indentation errors accumulated during the loading phase such as, e.g. errors derived from friction, surface roughness or the material inhomogeneities are involved in the experimental raw data [START_REF] Chen | Nanoindentation of hyperelastic polymer layers at finite deformation and parameter re-identification[END_REF][START_REF] Chen | Parameter re-identification in nanoindentation problems of viscoelastic polymer layers: small deformation[END_REF].

The parameters identification procedure was developed combining the FEM code Abaqus ® with HyperStudy ® software, which is a multi-disciplinary design exploration software in Altair ® HyperWorks ® . The software uses mathematical approaches such as design of experiments (DOE), optimization, and stochastic studies. Each approach serves a different purpose in the design study. In this study, the optimization approach was used for calibration and the DOE was used for sensitivity analysis.

The calibration procedure is illustrated in Fig. 3. The procedure can be described in the following steps: (1) the finite element model was developed in Abaqus ® (See section 2.2.1) and, then the Abaqus ® input file was imported in HyperStudy ® as well as the experimental force-displacement data. The numerical force-displacement history is extracted from the Abaqus ® output database file. Since the indenter was load-controlled, the displacement is predicted. Abaqus ® was specified as a solver in HyperStudy ® . [START_REF] Alam | Finite element analysis of forces of plane cutting of cortical bone[END_REF] the material constants were defined as design variables and the parameterization of the FE model was performed in HyperStudy ® . The initial guess values assigned to the unknown design variables and their lower and upper boundaries are shown in Table 1.

The upper bounds of the elastic modulus, the Poison's ratio and the initial yield stress were set according to their maximum values reported in the literature [START_REF] Bayraktar | Comparison of the elastic and yield properties of human femoral trabecular and cortical bone tissue[END_REF]11,[START_REF] Cowin | Bone mechanics handbook[END_REF][START_REF] Fan | Effects of viscoelasticity and time-dependent plasticity on nanoindentation measurements of human cortical bone[END_REF][START_REF] Rho | Elastic properties of microstructural components of human bone tissue as measured by nanoindentation[END_REF].

The upper bound of B was arbitrary chosen and the one of n is its maximum theoretical value. The initial values were taken from [START_REF] Alam | Finite element analysis of forces of plane cutting of cortical bone[END_REF] (see Table 3), which is the only study that used the Johnson Cook model for modeling, except for the initial value of n, an average value with respect to its theoretical range of values [0;1] has been chosen because we found that the reported value (n=0.08) is very close to that of the perfectly plastic solid (i.e., n=0) such as metal. The lower bounds were arbitrary chosen. (3) for the calibration process, a response has to be defined as an objective function required to be minimized in the optimization process. The curve difference function was defined as the objective function. The curve difference function is an HyperMath function that calculates the integral of the absolute value of the difference between two curves over the common domain of the supplied functions. The curve difference function is expressed as:

ˆb a |f (x) -g(x)| dx (5)
(4) in the optimization approach, the system identification procedure was used as the optimization problem formulation to minimize the objective function. This procedure consists to fit a set of responses to a set of target values by modifying a set of variables.

In the case of the current study, the numerical curve was fitted to the experimental curve by modifying the material constants using a least squares formula described as:

min f i -fi fi 2 (6)
where fi response is the target value of the i th response, and where the target values are the experimental data.

In For further details about this approach, see [44].

Design Of Experiments (DOE)

DOE is a mathematical approach developed to investigate the relationships between the different input variables of a design problem which influence the output responses [START_REF] Giunta | Overview of modern design of experiments methods for computational simulations[END_REF][START_REF] Tanco | Implementation of design of experiments projects in industry[END_REF]. The DOE process available in Altair HyperStudy ® was used to determine which material constants are most influential on the material responses to nanoindentation stress. The elastic (E and ν) and JC material constants were selected as design variables in the DOE process and the maximum equivalent stress as the material response. A DEO of type full factorial was used in this study. The full factorial DOE evaluates all 1 for the current study) and then plotted by drawing a straight line between the two average values corresponding to each design variable. The linear effect value of each design variable corresponds to the difference between its average responses. In the plot format, the effect is larger as the line becomes more vertical. A positive slope indicates a positive effect, ie., if a variable increases the response will also increase. In the full factorial DOE, the responses are evaluated for all combinations of design variable levels. The total number of runs of the full factorial DOE is given by this relationship: L v , where L is the number of level studies, which corresponds to the number of bounds, and v the number of variables. For the current study, we have five variables and two-level studies (ie., for the lower bound and for the upper bound). The total number of runs of the full factorial DOE is then 32 runs, 16 runs for each bound value.

Results and Discussion

Experimental results

The average values of the experimental indentation moduli and hardnesses are summarized in Table 2. These values, obtained with the OP approach available in the software included with the commercial Indentation Tester, are within the range of values reported in previous studies for the same approach [START_REF] Hoffler | An application of nanoindentation technique to measure bone tissue lamellae properties[END_REF][START_REF] Rho | Elastic properties of microstructural components of human bone tissue as measured by nanoindentation[END_REF][START_REF] Spiesz | Computational and experimental methodology for site-matched investigations of the influence of mineral mass fraction and collagen orientation on the axial indentation modulus of lamellar bone[END_REF][START_REF] Szabó | Anisotropy of bovine cortical bone tissue damage properties[END_REF]. A minimum difference of 4.6 GPa was found between the average values of the elastic modulus in the axial and transverse directions. Fig. 4 shows an example of the force-displacement curves corresponding to the axial and transverse directions. The curves show that the mechanical response of the cortical bone to nanoindentation depends on the loading direction, as found in the literature [START_REF] Carnelli | A finite element model for direction-dependent mechanical response to nanoindentation of cortical bone allowing for anisotropic post-yield behavior of the tissue[END_REF][START_REF] Ghanbari | Nonlinear hierarchical multiscale modeling of cortical bone considering its nanoscale microstructure[END_REF][START_REF] Vanleene | Influence of nanoindentation test direction on the elastic properties of human cortical bone lamellae[END_REF][START_REF] Szabó | Anisotropy of bovine cortical bone tissue damage properties[END_REF]. The axial direction was found to be the stiffest bone direction while the interstitial lamellae were found to be stiffer than the osteons lamellae, as reported by [START_REF] Fan | Anisotropic properties of human tibial cortical bone as measured by nanoindentation[END_REF][START_REF] Hoffler | An application of nanoindentation technique to measure bone tissue lamellae properties[END_REF][START_REF] Rho | Elastic properties of microstructural components of human bone tissue as measured by nanoindentation[END_REF][START_REF] Rho | Elastic properties of human cortical and trabecular lamellar bone measured by nanoindentation[END_REF][START_REF] Rho | Variations in the individual thick lamellar properties within osteons by nanoindentation[END_REF][START_REF] Vanleene | Influence of nanoindentation test direction on the elastic properties of human cortical bone lamellae[END_REF][START_REF] Zysset | Elastic modulus and hardness of cortical and trabecular bone lamellae measured by nanoindentation in the human femur[END_REF]. For the same site, a scatter in nanoindentation results was observed [START_REF] Faingold | Nanoindentation of osteonal bone lamellae[END_REF][START_REF] Rho | Variations in the individual thick lamellar properties within osteons by nanoindentation[END_REF]. Fig. 5 provides an illustration of scatter in the nanoindentation results of the interstitial bone site. A deviation of 2.36 GPa in elastic modulus was found between its minimum value and maximum value. This deviation corresponds to a deviation of 122 nm in depth. This difference was 1.48 GPa in osteons site.

Numerical results

As mentioned in the introduction, one significant problem in nanoindentation modeling of the elastic-plastic materials is pile-up deformations. This phenomenon doesn't occur in work-hardening materials [START_REF] Bolshakov | Influences of pileup on the measurement of mechanical properties by load and depth sensing indentation techniques[END_REF] and is not produced in cortical bone, as an elasticplastic material [START_REF] Fan | Three-dimensional finite element analysis of the effects of anisotropy on bone mechanical properties measured by nanoindentation[END_REF][START_REF] Johnson | Characterization of the microindentation of bone via scanning electron microscopy[END_REF][START_REF] Mullins | Calibration of a constitutive model for the post-yield behaviour of cortical bone[END_REF]. The cortical bone was then considered as a work-hardening elastic-plastic material. The purpose of this study was to apply a plastic hardening model to predict the mechanical behavior of the cortical bone in nanoindentation. The Johnson-Cook isotropic hardening model was chosen as the constitutive model. The finite element analysis was performed to solve the nanoindentation mechanical problem.

With this approach, it was expected to predict simultaneously: no pile-up of material, the force-displacement curve of the nanoindentation test and realistic yield stress values, ie, which should be in agreement with those reported in the literature. To the best of our knowledge, this study is the first one that used an isotropic hardening model in simulations of the nanoindentation testing of bone.

The mean values of the elastic (E and ν) and the JC cortical bone constants determined by the inverse analysis are summarized in Table 3. To illustrate the results, the finite element modeling of the nanoindentation test was performed using the identified constants. Fig. 6 shows a representative example of the fit between the experimental and the predicted curves. As we can see from this figure, the numerical curve matches well with the experimental one. Fig. 7 shows an example of the finite element modeling of the nanoindentation testing. The equivalent plastic strain (PEEQ) contours within the cortical bone sample indicates important plastic strains beneath the indenter and shows that the stress region is located far away from the boundaries, which suggests that the considered problem domain was enough for the numerical simulation of the nanoindentation testing. The figure also shows that there is no obvious pile-up of material around the contact impression, which is characterized by an upward vertical displacement of the outer boundary of the contact area for Berkovich indenters [START_REF] Bolshakov | Influences of pileup on the measurement of mechanical properties by load and depth sensing indentation techniques[END_REF][START_REF] Carnelli | A finite element model for direction-dependent mechanical response to nanoindentation of cortical bone allowing for anisotropic post-yield behavior of the tissue[END_REF][START_REF] Oliver | Measurement of hardness and elastic modulus by instrumented indentation: Advances in understanding and refinements to methodology[END_REF]. This agrees with the bone experimental results [START_REF] Fan | Three-dimensional finite element analysis of the effects of anisotropy on bone mechanical properties measured by nanoindentation[END_REF][START_REF] Johnson | Characterization of the microindentation of bone via scanning electron microscopy[END_REF][START_REF] Mullins | Calibration of a constitutive model for the post-yield behaviour of cortical bone[END_REF]. Another result, the predicted mean values of the initial yield stress, A, lie within the range of values [0.04-0.12] GPa, reported in the literature [START_REF] Alam | Finite element analysis of forces of plane cutting of cortical bone[END_REF][START_REF] Bayraktar | Comparison of the elastic and yield properties of human femoral trabecular and cortical bone tissue[END_REF][START_REF] Cowin | Bone mechanics handbook[END_REF]. The results presented thus far provide evidence that the JC model can predict accurately the elastoplastic behavior of the cortical bone. However, it would be noted that Table 3 shows a difference of about 0.5 GPa between the predicted and the corresponding measured mean values of the elastic modulus (Table 2). This difference can be due to the applied data analysis method [START_REF] Rodriguez-Florez | Insight into differences in nanoindentation properties of bone[END_REF], i.e., the proposed analysis method and the OP method. Indeed, the predicted elastic modulus was determined from the initial elastic contact and the experimental one was determined from the elastic response of the maximum contact area, beneath which a plastic straining was produced. On the other hand, unlike the OP method, which requires a priori knowledge of Poisson's ratio to compute the elastic modulus (equation 1), the proposed method allows to predict both the elastic constants independently and simultaneously. Previous investigation has proposed to decouple the elastic constants (E and ν) by using two different indenters to achieve accurate measurements [START_REF] Poon | An analysis of nanoindentation in linearly elastic solids[END_REF]. The reported difference in elastic modulus can also be related to errors in Poisson's ratio fixed in the experimental testing [START_REF] Poon | An analysis of nanoindentation in linearly elastic solids[END_REF][START_REF] Hay | A critical examination of the fundamental relations used in the analysis of nanoindentation data[END_REF], Poisson's ratio was fixed at 0.3 in the experimental testing while its predicted mean values varied from 0.28 to 0.48, which is in good agreement with other studies [START_REF] Cowin | Bone mechanics handbook[END_REF]. Fig. 8 provides an illustration of the variation in the force as a function of the displacement for two different values of Poisson's ratio. This variation (of about of 2.5 mN) corresponded to an overestimation of the elastic modulus of about 1 GPa. Finally, the predicted JC constants were found to be elevated compared to those obtained from the tensile testing of the cortical bone by [START_REF] Alam | Finite element analysis of forces of plane cutting of cortical bone[END_REF] (Table 3).

Linear effects of the elastic-plastic constants on the material response were evaluated using the full factorial DOE. The maximum equivalent stress was considered as the material response. Linear effects are shown in plot format in Fig. 9. For each material constant, the average response value is calculated for 16 runs of the full factorial DOE when the constant is at its lower bound value and 16 runs when the constant is at its upper bound value (see Table 1). The average values are indicated by the ends of the lines. The difference between the two average values, that is the linear effect value, is indicated in the legend. For example, for the lower bound of Young's modulus (15 GPa), an average maximum equivalent stress of 1.835 GPa was calculated for 16 results, whose minimal value was 0.015 GPa and the maximum value was 5.116 GPa. For its upper bound (30 GPa), an average maximum equivalent stress of 1.905 GPa was calculated for 16 results, whose minimal value was 0.0163 GPa and maximum value was 5.122 GPa.

The deducted linear effect value was then 0.069 GPa. As the figure shows, the Young's modulus, E, the strain hardening modulus, B, and the initial yield stress, A, have a positive effect and the Poisson's ratio, ν, and the strain hardening exponent, n, have a negative effect, which means an inverse effect on results, i.e., the maximum equivalent stress decreases when n increases [START_REF] Milani | An improved multi-objective identification of johnson-cook material parameters[END_REF]. Comparing the linear effect values displayed in the legends, and going from the least influence to the most influence on the material response within the current lower and upper boundaries, the parameters can be sorted as: E, A, ν, n, B. It is important to note that the order of the influence of the constants depends on the lower and upper boundaries.

As the full factorial DOE analysis has shown, the strain hardening modulus and the strain hardening exponent were the most influential parameters on the material response to nanoindentation stress. Thus, the stiffest bone region corresponded to the lowest magnitude of B and n. The strain hardening exponent measures the strain hardening ability of materials [11]. A perfectly elastic-plastic solid has a value of n=0. A mean values of n lying between 0.3 to 0.6 were predicted for the cortical bone. This ability of the cortical bone to strain-harden is probably due to its complex hierarchical structure, in particular its inorganic phase (33-43% bone volume) [START_REF] Hamed | Multiscale modeling of elastic properties of cortical bone[END_REF], whose the volume fraction varies with the indentation site and even with the single lamella. Indeed, the cortical bone is a dense connective tissue consisting of dense and hard concentric layers called lamellae [START_REF] Hamed | Experimentally-based multiscale model of the elastic moduli of bovine trabecular bone and its constituents[END_REF][START_REF] Hambli | Multiscale Computer Modeling in Biomechanics and Biomedical Engineering, 1st Edition[END_REF]. Each single lamella consists of mineralized collagen fibrils, a soft and ductile collagen fibrils reinforced with stiff hydroxyapatite crystals [START_REF] Jasiuk | Modeling of bone at a single lamella level[END_REF][START_REF] Lai | Mechanical behaviour of staggered array of mineralised collagen fibrils in protein matrix: Effects of fibril dimensions and failure energy in protein matrix[END_REF][START_REF] Vercher-Martinez | Influence of the mineral staggering on the elastic properties of the mineralized collagen fibril in lamellar bone[END_REF]. In each single lamella, the mineralized collagen fibrils, considered as the primary building block of bone [START_REF] Hamed | Multiscale modeling of elastic properties of cortical bone[END_REF], are oriented in one direction and rotate relative to the adjacent lamella [START_REF] Yoon | The estimated elastic constants for a single bone osteonal lamella[END_REF]. So, it has been suggested that the local material response depends mainly on the mineral volume fraction and the fibrils arrangement [START_REF] Evans | Microhardness and young's modulus in cortical bone exhibiting a wide range of mineral volume fractions, and in a bone analogue[END_REF][START_REF] Ghanbari | Nonlinear hierarchical multiscale modeling of cortical bone considering its nanoscale microstructure[END_REF][START_REF] Hamed | Multiscale modeling of elastic properties of cortical bone[END_REF][START_REF] Öchsner | Biomechanics of hard tissues[END_REF][START_REF] Porter | Pragmatic multiscale modelling of bone as a natural hybrid nanocomposite[END_REF]. A strengthening of the material could occur at the mineralized collagen fibrils level as a function of the tip penetration. Indeed, nanoindentation tests performed on hydroxyapatite single crystal revealed a nanoscale plasticity in the form of pile-up along the edges of the indenter [START_REF] Viswanath | Mechanical properties and anisotropy in hydroxyapatite single crystals[END_REF],

which may induce a uniform increase of the yield stress as the indenter penetrates the material. In a recent study, it was shown that a highly ductile behavior with continuous hardening until failure was observed on the micro-scale of a dry ovine osteonal bone samples subjected to micropillar compression [START_REF] Schwiedrzik | In situ micropillar compression reveals superior strength and ductility but an absence of damage in lamellar bone[END_REF]. The material strengthening could also occur because of lamellae dislocation through slip [START_REF] Schwiedrzik | In situ micropillar compression reveals superior strength and ductility but an absence of damage in lamellar bone[END_REF] as it has been found in metals (e.g., dislocations pile-up at grain boundaries [11]). In contrast to these findings, a further study suggested that other parameters could affect the indentation measurements such as nano-porosity, micro-damage, collagencross-links and non-collagenous proteins [START_REF] Spiesz | Computational and experimental methodology for site-matched investigations of the influence of mineral mass fraction and collagen orientation on the axial indentation modulus of lamellar bone[END_REF]. In a recent study [START_REF] Ramezanzadehkoldeh | Nanoindentation response of cortical bone: dependency of subsurface voids[END_REF], several cavities such as lacunae were identified below the nanoindentation sites of a mouse femur using FIB-SEM. However, a cavity-free material was found in the majority of the indentation sites.

Finally, the limitation of this work is that the bone tissue was assumed as an homogeneous material, i.e., its microstructure was not accounted for and the viscoelasticity and the anisotropy of the material [START_REF] Bekker | Introducing objective power law rate dependence into a visco-elastic material model of bovine cortical bone[END_REF][START_REF] Feng | Mechanical properties of porcine femoral cortical bone measured by nanoindentation[END_REF][START_REF] Olesiak | Viscous-elastic-plastic behavior of bone using berkovich nanoindentation[END_REF][START_REF] Zhang | Mechanical property determination of bone through nano-and micro-indentation testing and finite element simulation[END_REF] were not taken into the count in this study. Indeed, the strain rate effect on the mechanical behavior of the cortical bone has not been investigated, since the use of the JC strain rate dependence is required for high-strain-rate deformation of materials. For quasi-static tests, as is the case with the current study, the second factor in the constitutive relationship is neglected. To determine its anisotropic elastic properties, further tests need to be performed along multiple orientations since the cortical bone has a hierarchical and heterogeneous structure. Another limitation of this work is that the sensitivity of the JC material constants to the penetration depth and the indenter geometry was not investigated.

Summary and conclusions

Nanoindentation tests were performed on a cortical bone of a cadaver femur ;

The experimental results showed that the mechanical response of the material depends on the loading direction and for the same indentation site, a scatter in nanoindentation results was observed ;

The finite element analysis was used for modeling the nanoindentation testing ;

The material was considered as work-hardening elastic-plastic material ;

The JC isotropic hardening model was used to predict the mechanical behavior of the material ;

The inverse optimization approach was used for determining the unknown material constants from the force-displacement raw data in loading phase, as input data in the optimization process ;

The predicted elastic constants were determined from the initial elastic contact ;

The unknown material constants were predicted independently and simultaneously ;

The predicted values of the yield stress and the elastic constants were in agreement with the existing results in the literature ;

The numerical curves were well fitted to the target curves ;

No pile-up of material was predicted around the contact impression by the JC model ;

The DOE approach was used to investigate the linear effects of the material elasticplastic constants on the mechanical response of the material ;

It has been shown that the strain hardening modulus and the strain hardening exponent have the most significant effect on the material response ;

Finally, the results demonstrated that the JC model can predict accurately the elastoplastic behavior of the cortical bone in nanoindentation.

  the optimization step, several iterations are performed. During each iteration, the optimization algorithm updates the parameter values of the FE model until the optimal values are obtained. Several optimization methods are used by HyperStudy ® such as Adaptive Response Surface Method (ARSM), Global Response Surface Method (GRSM), Sequential Quadratic Programming (SQP), etc....In this study, SQP was used for the objective optimization. SQP is a gradient-based iterative optimization method.

Fig. 1 :Fig. 2 :

 12 Fig. 1: Optical microscopy image of human cortical bone obtained from a cadaver femur.

Fig. 3 :

 3 Fig. 3: Optimization process diagram.

Fig. 4 :

 4 Fig. 4: The force-displacement curves corresponding to the maximum depth (a) and the minimum depth (b) of the indenter tip for the axial and transverse directions.

Fig. 5 :

 5 Fig. 5: Loading curves obtained from nanoindentation testing performed on the interstitial bone. Eight indents were performed.

Fig. 6 :

 6 Fig. 6: Representative fit between experiment and JC model of a nanoindentation test performed on the interstitial lamellae of a cortical bone sample.

Fig. 7 :Fig. 8 :

 78 Fig.7: Finite element modeling of the nanoindentation test performed on a cortical bone. The contour map of the equivalent plastic strain (PEEQ) corresponds to the maximum indenter depth.

Fig. 9 :

 9 Fig.9: A Design Of Experiments (DOE) was performed to evaluate the linear effects of the elasticplastic Johnson-Cook constants on the mechanical response (here, the maximum equivalent stress) of a cortical bone to nanoindentation. The linear effect is evaluated by calculating the averaged responses at the lower and upper bounds of each constant fixed in the optimization process, and then plotted by drawing a straight line between the two average values. The difference between the two average values, that is the linear effect value, is indicated in the legend. For each material constant, the linear effect is calculated for 32 runs of the full factorial DOE. A positive slope indicates a positive effect.

Table 1 :

 1 Initial guess values and the lower and upper boundaries of the design variables in the optimization process. E 0 is the measured elastic modulus corresponding to each experimental test.

	variable name Lower Bound Initial Upper Bound
	E (GPa)	15	20	30
	ν	0.01	0.3	0.5
	A(GPa)	0.01	0.05	0.12
	B (GPa)	0.01	0.1	5
	n	0.01	0.5	1
	linear effects and interactions between the different design variables. The linear effect is
	evaluated by calculating the averaged responses at the lower and upper bounds of the
	design variables (given by Table

Table 2 :

 2 The average values of the experimental nanoindentation moduli and hardnesses and their standard mean deviations of human cortical bone compared to some values found in the literature. In parentheses, is shown the number of indentations.

		Direction	Axial (in O.R)	Axial (in I.R)	Transverse
	Current Study E (GPa)	22.13±0.74 (8)	23.9±1.18 (8)	17.45± 0.47 (10)
		H(GPa)	0.72±0.05 (8)	0.83±0.10 (8)	0.63±0.05 (10)
	[55] (1)	E (GPa)	22.4±1.20 (72)	25.7±1.00 (58)	16.6±0.74 (58)
		H(GPa)	0.617±0.039 (72) 0.736±0.044 (58) 0.564±0.034 (58)
	[33] (2)	E (GPa)	18.6±4.2	20.3±5.1	NR
		H(GPa)	0.52±0.15	0.59±0.2	NR
	[77] (3)	E (GPa)	19.1 ±5.4	21.2 ±5.3

NR

O.R: Osteonal Region. I.R: Interstitial Region. NR = not reported. (1) tibiae cortical bone samples of two male humans (ages of 60 and 66). (

2

) femoral diaphysis cortical bone samples of four human subjects (males: age 53 and 77; females: age 77 and 93). (3) femoral diaphysis cortical bone samples of four females and four males (age: from 53 to 93 yr old).

Table 3 :

 3 Mean values and standard mean deviations of the elastic and the JC cortical bone constants identified by the inverse analysis.

	Direction	E (GPa)	ν	A (GPa)	B (GPa)	n
	Axial (in O.R) 21.90±0.45 0.29±0,10 0.091±0,025 1.63±0.48 0.43±0.04
	Axial (in I.R) 23.07±1.98 0.38±0,06 0.094±0.007 1.29±0.19 0.34±0.13
	Transverse	18.27±1.28 0.44±0,04 0.043±0.033 1.77±0.30 0.52±0.07
	[2]	20	-	0.05	0.1	0.08
	O.R: Osteonal Region. I.R: Interstitial Region.

Acknowledgments

This work was supported by QUICKMOLD project.

Conflict of interest statement

None of the authors have any conflicts of interest.