Simon Le Floc'h

Guy Cloutier

Yoshifumi Saijo

Gérard Finet

K Saami

Flavien Yazdani

Gilles Deleaval

Roderic I Rioufol

Jacques Pettigrew

Ohayon

Saami K Yazdani

Flavien Deleaval

Gilles Rioufol

Roderic I Pettigrew

Professor Jacques Ohayon
email: jacques.ohayon@imag.fr.

A four-criterion selection procedure for atherosclerotic plaque elasticity reconstruction based on in vivo coronary intravascular

Keywords: Atherosclerosis, vulnerable plaques, coronary arteries, radial strain elastography, linear elasticity, inverse problem

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

INTRODUCTION

Atherosclerotic plaque rupture is the major cause of acute coronary syndrome, myocardial infarction and stroke in the western world [START_REF] Virmani | Lessons from sudden coronary death: a comprehensive morphological classification scheme for atherosclerotic lesions[END_REF]. The thin-cap fibroatheroma (TCFA) is the precursor lesion that once ruptured, may lead to the formation of a thrombus causing an acute syndrome and possibly death [START_REF] Virmani | Pathology of the vulnerable plaque[END_REF]. Typically, an unstable vulnerable coronary plaque (VP) has a large extracellular necrotic core thickness, a fibrous cap thickness less than 65 µm, and activated macrophages [START_REF] Virmani | Lessons from sudden coronary death: a comprehensive morphological classification scheme for atherosclerotic lesions[END_REF][START_REF] Ohayon | Necrotic core thickness and positive arterial remodeling index: emergent biomechanical factors for evaluating the risk of plaque rupture[END_REF]. Cellular mechanisms associated with plaque instability include inflammation, reduced collagen synthesis, local expression of collagenase, and smooth muscle cell apoptosis [START_REF] Falk | Coronary plaque disruption[END_REF]. Such biological processes alter mechanical properties of the plaque surface [START_REF] Arroyo | Mechanisms of plaque rupture: mechanical and biologic interactions[END_REF].

Coronary VPs can be detected clinically by various intravascular techniques, including ultrasound (IVUS) [START_REF] Rioufol | Multiple atherosclerotic plaque rupture in acute coronary syndrome: a three-vessel intravascular ultrasound study[END_REF][START_REF] Carlier | Studying coronary plaque regression with IVUS: a critical review of recent studies[END_REF], optical coherence tomography (OCT) [START_REF] Jang | Visualization of coronary atherosclerotic plaques in patients using optical coherence tomography: comparison with intravascular ultrasound[END_REF][START_REF] Tearney | Three-dimensional coronary artery microscopy by intracoronary optical frequency domain imaging[END_REF]) and magnetic resonance imaging (IV-MRI) [START_REF] Larose | Characterization of human atherosclerotic plaques by intravascular magnetic resonance imaging[END_REF][START_REF] Briley-Saebo | Magnetic resonance imaging of vulnerable atherosclerotic plaques: current imaging strategies and molecular imaging probes[END_REF]. Predicting coronary plaque rupture, on the other hand, is still imprecise, as imaging features of plaque morphology and composition are insufficient predictors of risk [START_REF] Loree | Effects of fibrous cap thickness on peak circumferential stress in model atherosclerotic vessels[END_REF][START_REF] Ohayon | Necrotic core thickness and positive arterial remodeling index: emergent biomechanical factors for evaluating the risk of plaque rupture[END_REF]. Prediction of rupture also requires knowledge of the plaque elasticity [START_REF] Ohayon | In-vivo prediction of human coronary plaque rupture location using intravascular ultrasound and the finite element method[END_REF][START_REF] Finet | Biomechanical interaction between cap thickness, lipid core composition and blood pressure in vulnerable coronary plaque: impact on stability or instability[END_REF], thus numerically estimating the peak stress amplitude of the fibrous cap, which appears to be a good biomechanical predictor of plaque rupture [START_REF] Loree | Effects of fibrous cap thickness on peak circumferential stress in model atherosclerotic vessels[END_REF][START_REF] Ohayon | Necrotic core thickness and positive arterial remodeling index: emergent biomechanical factors for evaluating the risk of plaque rupture[END_REF]. Therefore, the plaque elasticity reconstruction (i.e. modulogram) is an essential step for a reliable computation of intraplaque stresses.

Computation of such elasticity maps, based on the estimation of the intraplaque strain field (i.e. elastogram) obtained from various ultrasound-based techniques [START_REF] Wan | Strain imaging and elasticity reconstruction of arteries based on intravascular ultrasound video images[END_REF][START_REF] De Korte | Morphological and mechanical information of coronary arteries obtained with intravascular elastography; feasibility study in vivo[END_REF][START_REF] Kim | Vascular intramural strain imaging using arterial pressure equalization[END_REF][START_REF] Maurice | On the potential of the lagrangian estimator for endovascular ultrasound elastography: in vivo human coronary artery study[END_REF][START_REF] Zhang | Automatic segmentation of calcifications in intravascular ultrasound images using snakes and the contourlet transform[END_REF][START_REF] Hu | Intravascular ultrasound area strain imaging used to characterize tissue components and assess vulnerability of atherosclerotic plaques in a rabbit model[END_REF] and OCT [START_REF] Chan | OCT-based arterial elastography: robust estimation exploiting tissue biomechanics[END_REF][START_REF] Rogowska | Optical coherence tomographic elastography technique for measuring deformation and strain of atherosclerotic tissues[END_REF][START_REF] Van Soest | Robust intravascular optical coherence elastography by line correlations[END_REF]) remains a difficult task that has been tackled by a rather large diversity of approaches [START_REF] Vorp | Identification of elastic properties of homogeneous, orthotropic vascular segments in distension[END_REF]Soualmi et al. 1997;[START_REF] Vonesh | Regional vascular mechanical properties by 3-D intravascular ultrasound with finite-element analysis[END_REF][START_REF] Beattie | Mechanical analysis of heterogeneous, atherosclerotic human aorta[END_REF][START_REF] Kim | Vascular intramural strain imaging using arterial pressure equalization[END_REF][START_REF] Fehrenbach | Detection of small inclusions by elastography[END_REF][START_REF] Baldewsing | An inverse method for imaging the local elasticity of atherosclerotic coronary plaques[END_REF][START_REF] Cimrman | On the identification of the arterial model parameters from experiments applicable •in vivo[END_REF][START_REF] Pazos | Mechanical Characterization of Atherosclerotic Arteries using Finite Element Modeling: Feasibility Study on Mock Arteries[END_REF][START_REF] Franquet | Identification of heterogeneous elastic properties in stenosed arteries: a numerical plane strain study[END_REF][START_REF] Doyley | Model-based elastography: a survey of approaches to the inverse elasticity problem[END_REF]). Our group recently proposed a robust IVUS elasticity modulus approach for plaque imaging that combined a dynamic watershed segmentation (DWS) method with an optimization procedure (Le [START_REF] Floc | Vulnerable atherosclerotic plaque elasticity reconstruction based on a segmentation-driven optimization procedure using strain measurements: theoretical framework[END_REF]. This Young's modulus reconstruction approach 'iMOD' was successfully tested in vitro on polyvinyl alcohol cryogel vessel phantoms with controlled intraparietal stiffness inclusions (Le [START_REF] Floc | On the potential of a new IVUS elasticity modulus imaging approach for detecting vulnerable atherosclerotic coronary plaques: in vitro vessel phantom study[END_REF].

Plaque elasticity reconstruction models based on IVUS strain sequences in patients, however, are difficult due to motions created by the contraction and relaxation of the myocardium and the pulsatile nature of the coronary flow, preventing optimum strain field measurements and precise tracking of a same vessel cross-section during the whole cardiac cycle. Although original local sine-wave approach [START_REF] Arts | Mapping Displacement and Deformation of the Heart With Local Sine-Wave Modeling[END_REF], spatial priors techniques [START_REF] Richards | Investigating the impact of spatial priors on the performance of modelbased IVUS elastography[END_REF]Barbone et al. 2002) or more appropriate regularization technique [START_REF] Doyley | Model-based elastography: a survey of approaches to the inverse elasticity problem[END_REF]) could potentially overcome these hurdles, an alternative strategy would be to "precondition" the algorithm by selecting reliable input plaque strain and displacement fields, which is the main goal of the current study. We therefore designed a technique to extract accurate strain fields and modulograms from recorded IVUS sequences. We identified a set of four criteria based on (1) tissue overlapping, (2) RF-correlation coefficient between two successive frames, (3) performance of the elasticity reconstruction method to recover the measured radial strain, and (4) reproducibility of the computed set of modulograms over the cardiac cycle. This four-criterion selection procedure was successfully tested on IVUS sequences obtained in twelve patients referred for a directional coronary atherectomy (DCA) intervention.

METHODS

IVUS image analysis and plaque strain reconstruction

The first step of the IVUS data pre-processing consisted of segmenting IVUS images to detect vessel boundaries. The segmentation algorithm to detect the lumen and adventitia boundaries was based on a fast-marching model combining region and contour information (Roy [START_REF] Cardinal | Intravascular ultrasound image segmentation: a three-dimensional fast-marching method based on gray level distributions[END_REF]. Contours were validated by a cardiologist before the next processing step. The artery wall (area between the detected lumen and adventitia boundary) defined the region of interest (ROI) and was used for image registration before the radial strain estimation. The second step consisted of a rigid motion compensation to remove the rotation observed between consecutive images. Each IVUS frame made of 256 radial lines in polar coordinates was unwrapped to obtain a matrix where rows correspond to propagation depths and columns to angles. A lateral translation in this representation thus corresponded to a rotation centered on the middle of the catheter in the Cartesian system. The first 256 columns were selected as the initial image of the sequence. The second image was searched by lateral 2D correlations of regions of interest (ROIs) and the rotation artifact was determined as the lateral shift with maximum correlation. Next, the new image compensated for rotation served as the initial image to align the third one and so on. A complete sequence was reconstructed with all images compensated for rotation. Finally, the radial strain was computed using the Lagrangian speckle model estimator (LSME) [START_REF] Maurice | Adapting the Lagrangian speckle model estimator for endovascular elastography: theory and validation with simulated radio-frequency data[END_REF]). The first step of this algorithm consisted of a local rigid registration on overlapping sub-windows (measurement windows -MWs) within the ROI that allowed compensating for potential residual translation movement using a 2D cross-correlation analysis. Then, for each MW, displacement field and RF-correlation coefficient were computed. Notice that the LSME algorithm was formulated as a nonlinear minimization problem based on the optical flow equations solved for each MW. This technique allows ) [START_REF] Maurice | Adapting the Lagrangian speckle model estimator for endovascular elastography: theory and validation with simulated radio-frequency data[END_REF].

Plaque elasticity reconstruction method

The theoretical framework of this method has been described earlier [START_REF] Floc | Vulnerable atherosclerotic plaque elasticity reconstruction based on a segmentation-driven optimization procedure using strain measurements: theoretical framework[END_REF]. Briefly, considering two successive frames (at time t j-1 and t j ) of the IVUS sequence, we applied the imaging modulography technique (iMOD), which involves three successive steps: i) the computation of a pseudo-gradient elasticity map, ii) the DWS segmentation procedure that makes use of the previous step results to extract the inclusions' contours, and finally iii) the mathematical optimization procedure that provides the estimated Young's moduli of detected inclusions and surrounding tissue. The last step iii relies on using the adjoint method to estimate the gradient of the root mean squared error (RMSe rror ) between radial strains:

(1)

where n is the total number of nodes in the plaque mesh and x i the position of the plaque node i. In this process, a Young's modulus set solution was found acceptable when the gradient-based optimization procedure reached either a tolerance termination value lower than 10 -6 or a maximum number of iterations equal to 200. During the optimization procedure, local Young's moduli were constrained to remain positive and between 1 kPa and 10 4 kPa.

Four-criterion selection procedure used to extract plaque modulograms

For each patient's atherosclerotic lesion, approximately 30 displacement maps were reconstructed based on the IVUS sequence recorded during one second. We used four criteria to extract accurate displacement, strain and elasticity maps. The first criterion estimates the quality of the measured displacement field (displacement quality index Q displa ); the second criterion estimates the accuracy of the spatial distribution of the RFcorrelation coefficient between two successive frames (tracking quality index Q track ); the third criterion estimates the propensity of the Young's modulus reconstruction method to recover the measured radial strain (strain quality index Q strain ) and finally, the last criterion estimates the reproducibility of resulting modulograms over the cardiac cycle (reproducibility quality index Q reprod ).

Displacement quality index (Q displa )-The displacement field resulting from the LSME process was used to detect sites with inaccurate measurements corresponding to regions with tissue overlap, which induces strain artifacts. An overlapping region was defined as a site where the regular polar mesh is strongly distorted so that an element overlap its neighbouring elements. Such condition are not physically admissible. The positions of the nodes of the elements in the deformed configuration were obtained by using the displacement vector estimated at each node of the regular polar mesh thanks to the LSME algorithm. Therefore, to select only frames with few tissue overlap sites, a global displacement quality index (0% ≤ Q displa ≤ 100%) was defined at each time-frame t j as:

(2) where q displa is a local value associated to the position x i and equal either to 0 or 1 whether an overlapping is observed or not, respectively. The parameter n is the total number of pixels representing the cross-section image of the pathological artery. To minimize the effects of such inaccurate displacement sites on plaque elasticity reconstruction, only frames for which more than 60% of the plaque area was found without tissue overlaps (i.e. with local q displa value equal to 1) were kept for our analysis (i.e. with Q displa ≥60%).

Tracking quality index (Q track )-The tracking quality index results, obtained from the LSME process, highlight the robustness of the MW tracking method. This criterion Q track was defined for each timeframe t j as:

(3) where q correl is the local RF-correlation coefficient measured at position x i (0 ≤ q correl ≤ 1) [START_REF] Maurice | Adapting the Lagrangian speckle model estimator for endovascular elastography: theory and validation with simulated radio-frequency data[END_REF]). The amplitude of the criterion Q track varies from 0% ≤ Q track ≤ 100% ; low values indicate that almost all MWs were poorly tracked due to the out-of-plane motion of either the catheter or the arterial wall. At the opposite, high values of Q track point out accurate tracking of all MWs. Only frames for which Q track exceeded 60% were kept for our analysis.

Strain quality index (Q strain )-To estimate the reliability of computed strain fields obtained with the iMOD algorithm, we defined for each time-frame t j the following strain quality index:

(4) where 0% ≤ Q strain ≤ 100%. A high value of Q strain indicates a good match between the LSME-measured and iMOD-estimated strain fields. Only frames for which Q strain exceeded 60% were kept.

Reproducibility quality index of modulograms (Q reprod )-For a given IVUS sequence, only N* frames -for which all three quality indexes Q displa , Q track and Q strain exceeded 60% -were kept. The reproducible quality index determine the similarity between modulograms. Such criterion rejects modulograms that deviate from a specific mean modulus elastogram. This specific mean elasticity map was computed by averaging the best combination of the extracted K reproducible modulograms (2≤ K ≤ N*). Since we are looking simultaneously for the best reproducible combination of K modulograms and the resulting specific mean modulogram, extensive simulations were performed. Thus, for each combinations of K modulograms were selected and reproducibility quality indexes Q reprod were computed. This index -for a given j th modulogram (1≤ j ≤ K) -was estimated by computing the relative change between its value at time t j (i.e. the spatial Young's modulus map E(t j ,x i )) and the mean reference Young's modulus map Ẽ K (x i ) averaged over all K modulograms:

(5) with ( 6)

The best set of K reproducible modulograms was defined as the combination with the higher Q reprod mean value. From this best set, the final modulogram was chosen as the one with the higher Q strain value. In situations where Q reprod mean value remained less than 60% (i.e. when the modulograms are not reproducible), the final modulogram was chosen to be the one with the higher Q strain amplitude.

Validation of the four-criterion selection procedure-The elasticity reconstruction method was conducted by using a four-criterion selection procedure (4-CSP) in which the four criteria Q displa , Q track , Q strain and Q reprod were applied successively. However, to investigate the impact of each criterion, we selected modulograms by applying them individually on the whole database of 12 IVUS sequences of patients referred for DCA intervention.

Clinical data and histological study on excised DCA lesions

Selected patients (including two no re-flow cases and one perforation case) gave their written informed consent before being recruited. The research protocol was approved by the Review Ethical Committee of Sendai University hospital. All 12 patients were scheduled to undergo a DCA intervention. Routine IVUS observations were performed and RF signals were acquired at the same moment. IVUS scans were conducted with a Galaxy II® echograph (Boston Scientific, Natick, MA, USA) equipped with 40 MHz mechanically rotating probes. One sequence of approximately 30 images was acquired for each patient with a fixed longitudinal position of the catheter within the maximum stenosis, at a frame rate of 30 images/sec. Blood pressure increment between acquired RF images was assessed by assuming a linear variation between normal systemic diastolic and systolic phases.

The DCA procedure consisted in removing part of the atherosclerotic plaque with a catheter designed for radial lesion cutting and excision (Flexicut, Guidant Corporation, Santa Clara, California). The excised specimens were fixed in 10% formalin. The samples were then embedded in paraffin using standard protocols, and 4 µm transversal slices were obtained using a microtome for histological staining with Elastica-Masson's thrichrome. Such staining gave information on plaque constituents only. Indeed, due to the DCA technique, only the excised "deteriorated" portions of the atherosclerotic plaques and not the intact excised lesions were collected making it impossible to correlate between the modulagram image and the plaque morphology. We therefore investigated only the atherosclerotic plaque composition by performing the histological study.

The DCA procedure does not remove the whole lesion. A post DCA IVUS scan was performed for all patients. But unfortunately we were unable to identify the extracted part of the lesion by using this post DCA IVUS scan. After DCA procedure, the cross-section of the coronary (at the level of the treated lesion) becomes smaller (due to a global contraction process) making the superposition between the pre and post DCA IVUS images impossible. Therefore, we asked the cardiologist to roughly define on the pre DCA IVUS image the boundary of the lesion that was excised during the DCA procedure.

Image analysis of the histological samples-Histological stained samples were observed using a microscope, with images acquired concomitantly. Each image was then digitized and saved to a computer for postprocessing.

Atherosclerotic plaque composition and characterization-The elastica-Masson's trichrome (EMT) staining was used to subdivide atherosclerotic excised lesions into 7 distinct constituents (Fig. 1). The thrombotic (Thb) and the fibrin (Fi) regions consisted of high density of red blood cells and fibrin, respectively; the cellular fibrous (CeFb) region consisted of smooth muscle cells or other cells admixed with a low collagen content or elastic fiber; the hypocellular fibrous (HyFb) region contained extracellular connective tissue matrix with collagen and few cells; the collagen (Co) region was defined as site with high density of collagen fibers.

The intensity of the blue staining color reveals the amount of collagen content (strong blue indicates sites with high collagen content while light blue or no blue indicates sites with low collagen content). Constituent proportions were determined using the ImageJ software (ImageJ, NIH, Bethesda, MD, USA) and were expressed as a percentage of the total excised lesion area. We were unable to quantify the lipid-rich (LpRi) region since such constituent did not appear clearly in lesion samples. However, the presence of macrophage-derived foam-cells and lipid crystals was analyzed. From this histological analysis performed by biologists, plaque composition was subdivided in two groups. The first group named "low collagen content" included constituents without or with few collagen (i.e. Thrombosis + Fibrin + Cellular fibrosis). The second group included constituents with high collagen content (i.e. Hypocellular fibrosis + Collagen). Then both low and high collagen areas were expressed as a percentage of the total excised lesion area.

Correlation between lesion contents and plaque stiffness

Based on the study of [START_REF] Lee | Structure-dependent dynamic mechanical behavior of fibrous caps from human atherosclerotic plaques[END_REF] conducted on mechanical properties of human atherosclerotic plaques, two groups of stiffness were defined. The soft medium consisted of plaque regions with Young's moduli ranging between 5 kPa and 700 kPa, and the stiff medium was defined as lesion sites with Young's moduli > 700 kPa. These areas were expressed as a percentage of the total excised lesion area. Correlations between these two families with distinct rigidities (determined from computed modulograms) and plaque composition (obtained by histology) were analyzed using a commercially available software package (SigmaStat 3.5, Systat Software, Point Richmond, CA). Regressions with P values <0.05 were considered statistically significant. Since the plaque was partially removed during DCA intervention, the comparison study was conducted on a specific region of each modulogram by assuming that at post DCA, the coronary artery became concentric with a uniform wall thickness close to 1 mm. This assumption allowed us to define the specific modulogram's region that was correlated with the histological measurements of the excised lesion.

RESULTS

Performance of the four-criterion selection procedure (4-CSP)

A total of 337 measured displacement, RF-correlation coefficient and measured strain maps, and 254 computed strain and modulograms were computed. Table 1 provides details for each analyzed plaque, regarding the number of computed displacement fields, RFcorrelation coefficient maps, strain fields and modulograms at each step of the iterative 4-CSP. Using this procedure, we could select at least one modulogram over all available frames for 7 of the 12 patients (i.e. for 58% of the population). For one of these 7 patients (patient # 8), only one modulogram fulfilling the 4-CSP was obtained while more than 2 reproducible modulograms were computed for the 6 other patients.

Performance of the plaque elasticity reconstruction algorithm

Reproducible modulograms within a cardiac cycle

In 5 patients (patients # 2, 5, 9, 10 and 11, Table 1), we successfully extracted 3 reproducible modulograms from the IVUS scan of the one cardiac cycle available. Illustrations of such sets of reproducible modulograms for patients # 2 and 10 are presented in figures 4 and 5, respectively. The data revealed that approximately 2/3 of selected modulograms were computed during the decrease of the coronary pressure.

Severity of each individual criterion

Plaque elasticity reconstructions were also conducted by considering one-criterion (1-CSP) at a time (i.e. either Q displa , Q track , Q strain or Q reprod ). A quality index threshold of 60% was also selected for these computations. The most selective criterion was the strain quality index Q strain since the plaque elasticity reconstruction was possible for only 20% of the total number of frames. The other indexes Q displa , Q track and Q reprod were less restrictive since modulograms could be reconstructed for 79%, 77% and 70% of total frames, respectively (Fig. 6).

Importance of the threshold value for the 4-CSP

Interestingly, by reducing the threshold values of the four quality indexes Q displa , Q track , Q strain and Q reprod to 50%, we found more than 2 reproducible modulograms for 11 of the 12 patients (i.e. for 92% of the population, see Table 1).

Histological study of the excised plaque samples

Histological analysis was performed on the stained samples of the seven excised lesions for which modulograms could be obtained (Fig. 7). All lesions were found with significant densities of cellular fibrosis (CeFb), collagen (Co) and fibrin (Fi) (mean values: 37.5%± 21.7%, 25.3%± 17.5% and 22.2%± 16.2% of the total excised lesion area, respectively). Few hypocellular fibrosis (HyFb) and thrombosis (Thb) were found (mean values: 10.4%± 10.7% and 4.6%± 5.9%, respectively). Foam-cells and lipid crystals were present in 7/7 and 6/7 samples, respectively. Table 2 summarizes our findings.

Correlation between excised lesion content and tissue stiffness

A strong correlation (p < 0.035 and R=0.79) was found between the percentage of plaque area of distinct stiffness (namely [5 kPa-700 kPa] and > 700 kPa, derived from computed modulograms) and the percentage of plaque area with low collagen content (i.e. Cellular fibrosis + Thrombosis + Fibrin) and high collagen content (i.e. Hypocellular fibrosis + Collagen) obtained from the histological analysis, respectively (Fig. 8). Moreover, since stiffness range and plaque composition were divided into two groups, a similar correlation was found between the percentage of stiff plaque area (stiffness > 700 kPa, derived from the computed modulograms) and the percentage of plaque area with high collagen content (i.e. Hypocellular fibrosis + Collagen).

DISCUSSION

A critical key in prediction of plaque rupture is the accurate quantification of both the morphology and mechanical properties of the diseased artery [START_REF] Cheng | Distribution of circumferential stress in ruptured and stable atherosclerotic lesions. A structural analysis with histopathological correlation[END_REF][START_REF] Finet | Biomechanical interaction between cap thickness, lipid core composition and blood pressure in vulnerable coronary plaque: impact on stability or instability[END_REF]. Whereas several theories were developed to reconstruct the plaque elasticity map based on the strain field measurement [START_REF] Doyley | Model-based elastography: a survey of approaches to the inverse elasticity problem[END_REF], none of them preconditioned the algorithm by selecting reliable input strain and displacement fields to prevent several sources of artifacts that can be encountered in vivo. Therefore, in this work, we designed a four-criterion selection procedure for this purpose.

Our data demonstrate that performing the 4-CSP procedure is an essential step to obtain in vivo coronary IVUS strain sequences using reconstructed elasticity methods. The results showed strong correlations between the computed and the measured radial strain fields and between plaque constituents and Young's modulus amplitudes. Moreover, the reproducibility of the computed modulograms highlighted the robustness of the 4-CSP to extract accurate elasticity maps.

Is the 4-CSP sufficient for an accurate plaque elasticity reconstruction?

Due to the motion of the heart, frames acquired at the different blood pressures may be misaligned preventing optimum strain field measurements and precise tracking of a same vessel cross-section during the whole cardiac cycle. Therefore, in their IVUS strain elastography study, de [START_REF] De Korte | Characterization of plaque components and vulerability with intravascular ultrasound elastography[END_REF] used a similarity index between sequential echo frames which excludes such misaligned frames. We did not used the similarity index proposed by de Korte since our 4-CSP was already able to detect such misaligned frames. Indeed, our tracking quality index Q track allowed to exclude frames that were poorly tracked due to the out of plane motion of the catheter.

As suggested by a few authors [START_REF] Li | Elasticity reconstruction from displacement and confidence measures of a multi-compressed ultrasound RF sequence[END_REF][START_REF] Sumi | Regularization of tissue shear modulus reconstruction using strain variance[END_REF], the local RF-correlation coefficient q correl and the displacement index q displa can be used to increase the performance of the optimization scheme. Indeed, since our methodology allowed to clearly identify plaque sites with tissue overlaps and poor correlation coefficients, it is legitimate to disregard such locations from the optimization process for an optimum elasticity map reconstruction. To investigate this point, we modified the initial RMSe rr or cost function given by Eq. ( 1), by including the local RF-correlation coefficient q correl as a weighting term: (7) where n* is the total number of nodes in the plaque mesh excluding overlap sites (i.e. sites where q displa = 0). For each patient, we investigated the performance of this alternative cost by considering only frames for which the quality indexes Q displa , Q track and Q strain were above 60%. Interestingly, the modulograms computed with the modified cost function (Eq. 7) were found similar to those computed with the initial cost function (Eq. 1) since the mean reproducibility index was found larger than 90 % (Fig. 9). This result confirms the stability of the proposed method. Moreover, [START_REF] Doyley | Model-based elastography: a survey of approaches to the inverse elasticity problem[END_REF] in his elegant review article mentioned the uniqueness of the mathematical solution for such hard-prior elasticity reconstruction method which does not require any type of regularization since the reconstruction problem is well conditioned. However, such technique is prone to errors because the DWS preconditioning step may generate some errors during the segmentation procedure [START_REF] Doyley | Model-based elastography: a survey of approaches to the inverse elasticity problem[END_REF].

What is the most appropriate cardiac phase to estimate modulograms? Baldewsing and co-workers (2008) observed that the decrease of the arterial pressure was a more appropriate cardiac phase to extract reliable strain fields during in vivo IVUS sequence of a vulnerable plaque. Our results are in agreement, as 66% of our extracted reliable modulograms were obtained during the decrease of the coronary blood pressure.

Is the elastogram sufficient to detect vulnerable plaques?

Knowing that the rate of deformation (strain) of a tissue is directly related to the stress applied on it and to its mechanical properties, Ophir and colleagues [START_REF] Ophir | Elastography: a quantitative method for imaging the elasticity of biological tissues[END_REF][START_REF] Cespedes | Elastography: elasticity imaging using ultrasound with application to muscle and breast in vivo[END_REF] were the pioneers on developing imaging techniques based on the strain field. In this context, several elegant IVUS strain reconstructions (de [START_REF] De Korte | Morphological and mechanical information of coronary arteries obtained with intravascular elastography; feasibility study in vivo[END_REF][START_REF] Maurice | On the potential of the lagrangian estimator for endovascular ultrasound elastography: in vivo human coronary artery study[END_REF]) and palpography methods [START_REF] Schaar | Intravascular palpography for high-risk vulnerable plaque assessment[END_REF] were developed to investigate atherosclerotic coronary artery plaques and to predict their vulnerability to rupture. Techniques based on the optical flow [START_REF] Maurice | Adapting the Lagrangian speckle model estimator for endovascular elastography: theory and validation with simulated radio-frequency data[END_REF] or time-delay correlation estimation (de [START_REF] De Korte | Morphological and mechanical information of coronary arteries obtained with intravascular elastography; feasibility study in vivo[END_REF] allowed the calculations of intraplaque strain images during the cardiac cycle. Their main approach to detect VP was based on the idea that the deformation of soft plaque heterogeneities (i.e. the lipid core) should be much larger than those of stiff inclusions (i.e. calcified and fibrosis areas) in response to blood pressure changes (de [START_REF] De Korte | Morphological and mechanical information of coronary arteries obtained with intravascular elastography; feasibility study in vivo[END_REF]. However, these methods did not overcome a main limitation (de [START_REF] De Korte | The vulnerable atherosclerotic plaque: strategies for diagnosis and management[END_REF]) related to the complex geometries of atherosclerotic plaques, which alter the intraplaque strain fields and inhibit direct translation into Young's modulus maps. Figure 3 indeed demonstrates clearly such limitation. The strain map of patient # 8 (Fig. 3) highlights a high strain site located at the endoluminal layer between 6 and 8 o'clock. The strain field alone seems to indicate that the plaque is vulnerable and unstable with the possibility of a large necrotic core close to the lumen. However, the associated modulogram reveals a small soft region in the deeper plaque layer with a thick and stiff fibrous cap (cap thickness > 800 µm and mean Young modulus close to 335 kPa).

Study limitations

Several limitations deserve to be pointed out, even if the results include findings which may help to increase the performance of the modulography technique: First, despite our elasticity reconstruction maps being based on real in vivo IVUS data, correlations between patient's modulogram and the histological analysis of the associated excised lesion presented in this study were relevant, but not sufficient to broadly validate our computational findings. However, our findings are supported by a pre-clinical study showing that the plaque morphology and mechanical properties of the atherosclerotic lesion can be accurately derived by using the proposed modulography approach (Le [START_REF] Floc | On the potential of a new IVUS elasticity modulus imaging approach for detecting vulnerable atherosclerotic coronary plaques: in vitro vessel phantom study[END_REF].

Second, an important study limitation remains the assumption made on blood pressure. Indeed, since blood pressure was not measured during the IVUS acquisition sequence, the modulograms were obtained by assuming a bilinear pressure-time evolution during cardiac cycle. Such assumption may bias the resulting modulogram. However, by using a linear elasticity theory to solve the inverse problem, the values of the Young's modulus solutions appear to be proportional to the pressure gradient amplitude. Therefore one option to overcome such limitation, would have been to normalize the Young's modulus by the pressure gradient amplitude to obtain a normalized modulogram. Such option was not chosen since the mean relative error made by assuming such simplified linear temporal pressure distribution was found lower than 26% (results not shown). We quantified such error by comparing the pressure gradient amplitude obtained with the bilinear curve with the one obtained by considering a more realistic pressure-time curve.

Third, the proposed technique for rigid motion compensation between consecutive frames does not correct for the non uniform rotation of the single ultrasound probe element that may affect the resulting strain elastogram.

Finally, the influence of residual stresses [START_REF] Matsumoto | Residual stress and strain in the lamellar unit of the porcine aorta: experiment and analysis[END_REF]), generated by remodeling and plaque growth processes, has been ignored in this study. Residual stress patterns have been studied ex vivo in several human vulnerable coronary plaque samples [START_REF] Ohayon | Influence of residual stress/strain on the biomechanical stability of vulnerable coronary plaques: potential impact for evaluating the risk of plaque rupture[END_REF]. Neglecting residual stresses in structural analysis could bias the mechanical stress quantification but not the characterization of the mechanical properties, which was the main goal of this study.

Clinical implications

Increased arterial stiffness is associated with markers of cardiovascular risk [START_REF] Glasser | Vascular compliance and cardiovascular disease: a risk factor or a marker?[END_REF][START_REF] Van Bortel | Pulse pressure, arterial stiffness, and drug treatment of hypertension[END_REF]. Our group recently showed that local increase of vessel stiffness resulting from wall-strain stiffening phenomenon is critical in the evolution of atherosclerotic lesions [START_REF] Ohayon | Is arterial wall-strain stiffening an additional process responsible for atherosclerosis in coronary bifurcations? an in vivo study based on dynamic CT and MRI[END_REF]. The proposed LSME-iMOD imaging method may help in localizing such sensitive zones.

Several animal and clinical studies conducted to analyze the structural variation in the fibrous cap and necrotic core demonstrated that statin treatment enhances plaque stability [START_REF] Libby | Inflammation and atherosclerosis[END_REF]. Additionally, as suggested in a previous study [START_REF] Finet | Biomechanical interaction between cap thickness, lipid core composition and blood pressure in vulnerable coronary plaque: impact on stability or instability[END_REF]), a very slight increase in the mechanical properties of plaque constituents, namely the hardening of the necrotic core, can tilt a VP from instability to stability. The proposed LSME-iMOD imaging method may be used to explain such findings [START_REF] Nozue | Statin treatment for coronary artery plaque composition based on intravascular ultrasound radiofrequency data analysis[END_REF]) by investigating the evolution of mechanical properties of atherosclerotic plaques during the statin therapy. Performance of the four-criterion selection procedure to extract morphology and elasticity map from the IVUS sequence of patient # 9. Notice that regions with tissue overlap, observed in the displacement map, correlate well with low RF-correlation coefficient sites. The high value of the strain quality index found (Q strain =67%) indicates that the method reproduced reasonably the measured radial strain field. The green contour in the modulogram corresponds to the boundary of the lesion that was excised during the DCA procedure. Performance of the four-criterion selection procedure to extract morphologies and elasticity maps from IVUS sequences of patients # 8, 11 and 12. Plaques of patients # 11 and 12 have large necrotic cores (dark regions in the modulograms). The vulnerable plaque of patient # 8 appears to be likely stable since the small necrotic core at 7 o'clock is far from the lumen. The green contours in the modulograms correspond to the boundaries of the lesions that were excised during the DCA procedure. Performance of the four-criterion selection procedure to extract reproducible modulograms from the IVUS sequence of patient # 10. Approximately 2/3 of the computed modulograms were obtained during the decrease of the coronary blood pressure (corresponding to an increase of the cumulated radial strain). The green contours in the modulograms correspond to the boundaries of the lesion that was excised during the DCA procedure. Performance of the four-criterion selection procedure to extract reproducible modulograms from the IVUS sequence of patient # 2. The strong correlations found between the computed and the measured radial strain fields (Q strain > 70%) and the significant reproducibility of the three resulting modulograms highlight the robustness of the method. The green contours in the modulograms correspond to the boundaries of the lesion that was excised during the DCA procedure. Results regarding the severity of each criterion. The strain quality index Q strain was found to be the most selective since modulograms were reconstructed for 20% of the total recorded frames only. Microscopic images of the histological stained samples of the seven excised lesions for which modulograms were obtained. The elastica-Masson's trichrome (EMT) staining was used to characterize the composition of the atherosclerotic excised lesions. These stained samples were used for the histological study. A good agreement was found when comparing the percentage of soft plaque area (i.e. 5 kPa < stiffness < 700 kPa) values derived from the computed modulograms with the percentage of plaque area with low collagen content (i.e. Cellular fibrosis + Thrombosis + Fibrin) obtained from the histological analysis. The linear regression line (dark line) and the confidence of the prediction interval domains, which are the domain limited by the two internal and external dashed lines, respectively, are also given. Performance of the considered cost function (see Eq. 1) in the optimization algorithm. Improving the cost function (see Eq. 7) by restricting the optimization process to only regions with no tissue overlap and with high RF-correlation coefficients did not modify the resulting modulograms. These simulations were performed on patient # 5. This result demonstrates the stability and robustness of the proposed four-criterion selection procedure. The green contours in the modulograms correspond to the boundaries of the lesion that was excised during the DCA procedure. Number of computed displacement fields, RF-correlation coefficient maps, strain fields and modulograms at each step of the iterative four-criterion selection procedure. Simulations were performed with two criteria thresholds (namely 50% and 60%). The four criteria Q displa , Q track , Q strain and Q reprod were used to extract accurate modulograms with regard to the quality of the displacement, RF-correlation coefficient, strain and reproducible elasticity maps, respectively.

Figure 2

 2 Figure 2 displays the best computed modulogram and the associated displacement, RF-) maps obtained for patient # 9. The modulogram reveals an extensive endoluminal soft region between 8 and 12 o'clock (mean Young's modulus < 25 kPa). This soft region correlates well to the high strain site and corresponds to a region with almost no tissue overlap and with high values of the local RF-correlation coefficients. Figure 3 presents modulograms computed for plaques of patients # 8, 11 and 12. The strong similarities found between the computed and measured strain fields (in both Figs. 2 and 3) point out the performance of the Young's modulus reconstruction technique combined with the 4-CSP.

Figure 1 .

 1 Figure 1. Images showing morphologically distinguishable subtypes of atherosclerotic plaque constituents. A) Thrombotic (Thb) region. B) Fibrin (Fi) region. C) Cellular fibrous (CeFb) region. D) Hypocellular fibrous (HyFb) region. E) Collagen (Co) region. F) Lipid crystals (arrows). G) Macrophage-derived foam-cells (arrow). See text for complete description.

  Figure 2.

  Figure 3.

  Figure 4.

  Figure 5.

  Figure 6.

  Figure 7.

  Figure 8.

  Figure 9.

  Histological measurements. Excised lesions composition, expressed as % total plaque area. Thb, Thrombosis; Fi, Fibrin. CeFb, Cellular fibrosis; HyFb, Hypocellular fibrosis; Co, Collagen.

Table 1

 1

Acknowledgments

The authors thank Arnold Fertin, PhD. (TIMC Laboratory, Grenoble, France) for his help on image analysis processing. Acknowledgements are also addressed to Dr Younes Majdouline and Dr Marie-Hélène Roy Cardinal of the Laboratory of Biorheology and Medical Ultrasonics, University of Montreal Hospital Research Center, for their computational supports regarding segmentation and radial strain elastograms.

Grants

Grant supports were provided by the Agence Nationale de la Recherche (ANR), France (ATHEBIOMECH project), and by the collaborative health research joint program of the Natural Sciences and Engineering Research Council of Canada (NSERC #323405-06) and Canadian Institutes of Health Research (CIHR #CPG-80085). This research is now supported by a joint international program of the ANR (MELANII project # 09-BLANC-0423) and NSERC strategic grant #STPGP-381136-09. Financial support was also provided by a post-doctoral appointment (S. Le Floc'h) from a valorization research program of the University Joseph Fourier (UJF), France in partnership with Floralis-France.


C A Plastiras

E Berlioz

G Merceron

F Guy

D S Kostopoulos

Exploring the feeding ecology and habitat of the primate from Dafnero-3

Introduction

The locality of Dafnero is an Early Pleistocene locality situated in Kozani, Northwestern Greece, that was discovered in the early 90's [START_REF] Koufos | A new mammalian locality in the Villafranchian of western Macedonia (Greece)[END_REF]. The aim of our study is to characterize the environmental context and feeding habits of the fossil genus Procynocephalus/ Paradolichopithecus (P/P) from the Early Pleistocene site of Dafnero-3 (MNQ-17) [START_REF] Kostopoulos | A 2Ma-old baboon-like monkey from Northern Greece: Procynocephalus aff. arvernensis (Primates:Cercopithecidae)[END_REF]. To achieve that, we explore the dental microwear texture of Eucladoceros ctenoides (n=3) and Metacervoceros rhenanus (n=6) from Dafnero, because the microwear texture of herbivores, as direct plant consumers, is a good proxy to characterize paleo-environments. To do so, our results are compared with those of two contrasted reference populations of extant Cervus elaphus (n=116), Lugar Nuevo, Southern Spain and Białowieża (n=23), Eastern Poland. Also the dental microwear texture data and 3D enamel thickness (micro CT scan imaging) of the Procynocephalus/Paradolichopithecus (DFN3-150) are explored and two indices that describe the molar flare are measured: molar flare index (MF) and dentine flare index (DFI) (F.Guy). The results are compared with a set of fossil Papio specimens (n=75) and 6 primate species with different dietary habits.

Materials and Methods

Results

Fig. 2 3D elevation model of DFN3-150 dental facet 9 microwear texture.

Dental Microwear Textural Analysis

(DMTA) i. Eucladoceros ctenoides (n=3), Metacervoceros rhenanus (n=6) from Dafnero. ii. Cervus elaphus, Lugar Nuevo (n=116), Białowieża (n=23). iii.P/P (DFN3-150). iv.Fossil Papio (n=75). v. 6 extant forms [START_REF] Scott | Dental Microwear Texture and Anthropoid Diets[END_REF].

3D Enamel Thickness, Molar Flare

and Dentine Flare. i. Papio anubis (n=3), Papio hamadryas (n=1) and P/P from Dafnero-3. ii. Micro CT, CT-scans using EasyScan (Rx solutions) (Centre de Microtomographie at the University of Poitiers). iii.© Avizo v7 ,© Geomagic studio 2013.

Dental capabilities Dietary habits

Environmental context of Dafnero-3 P/P Fig. 3 Schematic representation of MF index (left) , which is the intercuspal distance divided by the distance at the base of the molar (red lines).

The DFI index (right) is the area of the dentine basin divided by the area of the whole dentine, both projected in two dimensions.

Fig. 4 The 3D average enamel thickness (left) shows that DFN3-150 possessed thicker average enamel than all the other taxa compared. However this measurement is biased by the size effect. If we remove the size effect (right), still DFN3-150 show the higher median value of 3D relative scale free enamel thickness. Fig. 7 The position of the taxa on the plot derives from the percentage of the specimens that surpassed the cutpoint values for Asfc and epLsar that we have set (1.5 and 0.003 respectively), the bubble width is proportional of the sample size. Even though P/P is represented only by two specimens, based on its position on the plot we can presume that hard food objects was a part of its diet and the dietary habits are similar to fossil and modern day baboons.

Concluding remarks

• In our study, the fossil deer from Dafnero are characterized by a diet plasticity that is not surprising if compared with the plasticity in dietary habits of modern Cervus elaphus but is at odds with traditional paleoenvironmental interpretations that consider the presence of fossil deer in a locality as a sufficient paleoenvironmental proxy to draw conclusions in terms of the existence of a tree cover in the habitat. The dental microwear textures of deer, rather than their occurrence, constitutes an adequate paleo-habitat proxy.

• The dental microwear textures of Dafnero deer are in favor of the availability of both browse and abrasive herbaceous monocots in the habitat, a result which is coherent with previous interpretations drawn on the basis of the faunal assemblage [START_REF] Kostopoulos | Palaeoecological remarks on Plio-Pleistocene mammalian faunas, in: Interim Colloquium RCMNS[END_REF].The lack of frank signal is indubitably influenced by the really small sample size for Eucladoceros ctenoides and Metacervoceros rhenanus, that has to be taken into consideration in order not to overvalue the interpretations.

• Comparisons of enamel thickness indicate that DFN3-150 from Dafnero-3 possessed thick enamel.

• The comparisons of the MF and DFI indices show that P/P from Dafnero-3 possessed pronounced flare.

• The dietary strategy and habits of P/P from Dafnero-3 resemble the ones of modern and fossil baboons and possibly was dependent on seasonal and spatial availability of food resources.

Fig. 1

 1 Fig. 1 The locality of Dafnero, Northwestern Greece (left) and the fossiliferous site DFN3 (right) (photos by D.S Kostopoulos).

Fig. 5

 5 Fig. 5 Comparisons of the molar flare index (MF) (left) and dentine flare index (DFI) (right). The molar flare index value places P/ P from Dafnero-3 with the hominoids with diet consisted of hard food objects, but this comparison is not in the appropriate phylogenetic and temporal context (fide Singleton, 2003). The dentine flare index places P/P from Dafnero-3 higher than the Papio sample indicating more pronounced flare, and this comparison seems more appropriate because it is in the same phylogenetic and relatively same temporal context.

Fig. 6

 6 Fig.6The wide range of dental microwear textures for both Eucladoceros ctenoides and Metacervoceros rhenanus from Dafnero, occupying the entire ecospace between the reference localities, reflects the high plasticity of their diet. While there is no noticeable tendency in the textural signal of Eucladoceros ctenoides (due to a small sample size; n=3), most Metacervoceros rhenanus (n=6) occupy a grazing ecospace similar to the extant red deer from Lugar Nuevo.


Djafar Iabbaden
email: djafar.iabbaden@univ-st-etienne.fr

Jonathan Amodeo

Claudio Fusco

Florence Garrelie

Jean-Philippe Colombier
email: jean.philippe.colombier@univ-st-etienne.fr

Molecular dynamics simulation of structural evolution in crystalline and amorphous CuZr alloys upon ultrafast laser irradiation

A promising route for manufacturing emergent metamaterials is the use of ultrafast laser pulses in complex alloy processing. Laser-induced structural transformations in Zr based alloys for crystalline and glassy states are investigated here. The ultrafast thermomechanical response is compared between relevant stable crystalline structures B2-Cu50Zr50, C11 b -Cu33.3Zr66.7 and the amorphous structures α-Cu50Zr50, α-Cu33.3Zr66.7. The subsurface modification resulting from ultrafast laser irradiation is investigated by a hybrid simulation to capture the phenomenon occurring at a picosecond time scale. This combines a two temperature model and molecular dynamics approaches to simulate laser matter interaction at the mesoscale. Our results indicate that the involved structural dynamics strongly depend on the initial atomic composition and phase structure. In particular, a martensite phase transition is unveiled for B2 crystalline alloy, and defects are induced in the irradiated C11 b phase, whereas the amorphous state of photoexcited metallic glasses remains remarkably preserved.

I.

INTRODUCTION

Upon ultrafast laser structuring process, light is temporally confined allowing photoexcitation of band gap materials and a quasi-isochoric heating of metallic or semi-metallic materials [START_REF] Gamaly | The physics of ultra-short laser interaction with solids at non-relativistic intensities[END_REF]. At the same time, beam focusing and polarization manipulation allow the energy concentration at the nano-scale [START_REF] Stoian | Advances in ultrafast laser structuring of materials at the nanoscale[END_REF] and significant results evidenced that phase transformation can be triggered in material alloys, implying both glasses and crystalline materials [START_REF] Danz | Ultrafast nanoimaging of the order parameter in a structural phase transition[END_REF][START_REF] Takeda | Ultrafast crystalline-to-amorphous phase transition in ge2sb2te5 chalcogenide alloy thin film using single-shot imaging spectroscopy[END_REF][START_REF] Buividas | Phase transformation in laser-induced micro-explosion in olivine (fe, mg) 2sio4[END_REF][START_REF] Bressel | Femtosecond laser induced density changes in geo2 and sio2 glasses: fictive temperature effect[END_REF]. Subpicosecond laser irradiation shorter than electron-phonon relaxation time is able to cause structural distortion in solid lattices, inducing orderdisorder changes in covalently-bonded materials [START_REF] Sokolowski-Tinten | Thermal and nonthermal melting of gallium arsenide after femtosecond laser excitation[END_REF][START_REF] Sundaram | Inducing and probing nonthermal transitions in semiconductors using femtosecond laser pulses[END_REF], metal-insulator transitions [START_REF] Cavalleri | Femtosecond structural dynamics in vo 2 during an ultrafast solid-solid phase transition[END_REF][START_REF] Wall | Ultrafast changes in lattice symmetry probed by coherent phonons[END_REF]. The typical martensitic transformation is the face-centered-cubic to bodycentred-cubic transformation in steels [START_REF] Urbassek | Molecular dynamics modeling of martensitic transformations in steels[END_REF]. Martensite structural changes occur in shape memory alloys, and also in pure elements such as tantalum, and titanium [START_REF] Lu | Phase transformation in tantalum under extreme laser deformation[END_REF][START_REF] Petry | Phonon dispersion of the bcc phase of group-iv metals. i. bcc titanium[END_REF], where a body-centred-cubic to hexagonal-closepacked transformation is observed [START_REF] Li | Structural phase transformation in ni-hf and ni-ti systems studied by molecular dynamics simulation[END_REF]. Besides, femtosecond (fs) irradiation of titanium alloys revealed the opportunity to tune structural states of their surface, through a laser-driven α to β phase transformation related to the combined effect of ultrafast heating and cooling or strong shock waves [START_REF] Kolobov | Structural transformation and residual stresses in surface layers of α+ β titanium alloys nanotextured by femtosecond laser pulses[END_REF]. Recent theoretical predictions highlighted that direct solid-solid transition can be transiently triggered in pure metals by ultrafast laser due to the strong coupling between electronic and ionic structures but remains to be experimentally confirmed for permanent crystalline states [START_REF] Murphy | Dynamical simulations of an electronically induced solid-solid phase transformation in tungsten[END_REF][START_REF] Ben-Mahfoud | Structural instability of transition metals upon ultrafast laser irradiation[END_REF]. Metallic glasses, usually formed by liquid quenching, have also shown capability for undergoing crystallization during molding on the micrometer (µm) to nanometer (nm) scale owing to the low viscosity in the supercooled liquid state [START_REF] Kumar | Nanomoulding with amorphous metals[END_REF]. This behaviour remains unexplored in response to ultrashort laser irradiation.

In this work, we employ atomistic modeling to investigate laser-induced architecturation of CuZr alloys. This reveals how mechanical stress induced by shock wave affects the structural features, and the stability of certain compositions. The energy dynamics of the photoexcited solid is investigated by employing the Two Temperature Model (TTM) [START_REF] Anisimov | Electron emission from metal surfaces exposed to ultrashort laser pulses[END_REF]. This approach is used to describe the thermal diffusion and the energy transfer from the excited electrons to the colder lattice under an electron-phonon non-equilibrium state. In order to understand the response to this swift heating, it is necessary to couple the TTM to classic hydrodynamics or Molecular Dynamics (MD). TTM-MD has successfully been implemented to provide insights on phase transition [START_REF] Duffy | Including the effects of electronic stopping and electron-ion interactions in radiation damage simulations[END_REF][START_REF] Rutherford | The effect of electron-ion interactions on radiation damage simulations[END_REF], defect formation [START_REF] Lin | Computational study of the generation of crystal defects in a bcc metal target irradiated by short laser pulses[END_REF][START_REF] Sedao | Growth twinning and generation of high-frequency surface nanostructures in ultrafast laserinduced transient melting and resolidification[END_REF] or laser-ablation mechanisms [START_REF] Ivanov | Combined atomisticcontinuum modeling of short-pulse laser melting and disintegration of metal films[END_REF]. Most of the reported works using MD approaches explored the interaction of an ultrafast laser with pure metals as Al [START_REF] Wu | Microscopic mechanisms of laser spallation and ablation of metal targets from largescale molecular dynamics simulations[END_REF], Cr [START_REF] Karim | Molecular dynamics study of femtosecond laser interactions with cr targets[END_REF], Cu [START_REF] Karim | A computational study on the microstructural evolution in nearsurface copper grain boundary structures due to femtosecond laser processing[END_REF] and Au [START_REF] Stegailov | Atomistic simulation of laser ablation of gold: The effect of electronic pressure[END_REF], semiconductors as Si [START_REF] Lee | Applying tersoff-potential and bond-softening models in a molecular dynamics study of femtosecond laser processing[END_REF] or even graphite [START_REF] Jeschke | Theory for the ultrafast ablation of graphite films[END_REF]. A first attempt to understand the behavior of a CuZr under a ultrafast laser shot was done in [START_REF] Marinier | Femtosecond laser ablation of cux zr1-x bulk metallic glasses: A molecular dynamics study[END_REF] but no phase transitions were detected especially in the crystalline states. The choice of a CuZr system is not fortuitous, since in its amorphous state this specific Metallic Glass (MG) plays a fundamental role in advanced material engineering, from biomedical applications for antimicrobial properties [START_REF] Chen | Antimicrobial properties of zr-cu-al-ag thin film metallic glass[END_REF] to metallurgy and anticorrosion properties [START_REF] Yu | Corrosion behavior of bulk (zr58nb3cu16ni13al10) 100-xyx (x= 0, 0.5, 2.5 at.%) metallic glasses in sulfuric acid[END_REF]. Recently, this kind of amorphous system was irradiated by ultrafast laser for surface nanostructuring applications [START_REF] Li | Initial cumulative effects in femtosecond pulsed laser-induced periodic surface structures on bulk metallic glasses[END_REF][START_REF] Prudent | Initial morphology and feedback effects on laser-induced periodic nanostructuring of thin-film metallic glasses[END_REF]. Single shot ultrafast laser interaction with metallic alloys has been investigated in nanostructuring conditions far from equilibrium conditions. Unusual phase transformation emerges at the surface or in the subsurface region of the ultrafast-irradiated target. Phase transitions are expected to originate from atomic rearrangements either due to thermoelastic stresses [START_REF] Titus | Dislocation injection in strontium titanate by femtosecond laser pulses[END_REF][START_REF] Echlin | Materials response to glancing incidence femtosecond laser ablation[END_REF][START_REF] He | Computational study of short-pulse laserinduced generation of crystal defects in ni-based singlephase binary solid-solution alloys[END_REF] in response to the propagation of the pressure wave throughout the system relaxation or during the resolidification process that proceeds under high undercooling conditions [START_REF] Ivanov | Nanocrystalline structure of nanobump generated by localized photoexcitation of metal film[END_REF][START_REF] Wu | Generation of subsurface voids and a nanocrystalline surface layer in femtosecond laser irradiation of a single-crystal ag target[END_REF][START_REF] Sedao | Growth twinning and generation of high-frequency surface nanostructures in ultrafast laser-induced transient melting and resolidification[END_REF]. This paper aims to settle the relative role of both processes by focusing on the CuZr binary alloys with different compositions in both crystal and amorphous structures, in particular for α-Cu 50 Zr 50 and α-Cu 33.3 Zr 66.7 . Our results demonstrate that phase stability upon laser pulse and microstructural transformations depend on the initial crystallographic structure. The semiempirical potential used to model the vitrification of the Cu-Zr alloys confirms that the initial concentrations (50 % and 66.7 % of Zr) are stable and leads to a reasonable agreement with the experimental and phase diagram observations [START_REF] Zhou | Phase stability for the cu-zr system: first-principles, experiments and solution-based modeling[END_REF]. In contrast, thermodynamical response is substantially different for amorphous systems irradiated by similar laser conditions. No structural transformation was observed for the evaluated stoichiometry as MG accommodates the distortion even in the lack of long-range translational order maintaining mechanical stability.

II. COMPUTATIONAL MODELLING

The numerical simulations are carried out using the LAMMPS software [START_REF] Plimpton | Fast parallel algorithms for short-range molecular dynamics[END_REF] in which the CuZr systems are computed using the Embedded Atom Model (EAM) of Mendelev et al. [START_REF] Mendelev | Development of a semi-empirical potential suitable for molecular dynamics simulation of vitrification in cuzr alloys[END_REF]. A timestep of ∆t = 1 fs is used throughout the simulations in which the local atomic structures are investigated using the Polyhedral Template Matching Algorithm (PTM) [START_REF] Larsen | Robust structural identification via polyhedral template matching[END_REF] as implemented in the OVITO software [START_REF] Stukowski | Visualization and analysis of atomistic simulation data with ovitothe open visualization tool[END_REF].

A. TTM-MD Theoretical Framework

The TTM model is executed using the USER-MISC package [47] derived from cascade algorithm simulations [START_REF] Jay | Simulation of single particle displacement damage in siliconpart ii: Generation and long-time relaxation of damage structure[END_REF]. This LAMMPS package accounts for electronic effects and provides the implementation of the TTM-MD hybrid method algorithm as described in [START_REF] Duffy | Including the effects of electronic stopping and electron-ion interactions in radiation damage simulations[END_REF][START_REF] Rutherford | The effect of electron-ion interactions on radiation damage simulations[END_REF][START_REF] Chen | A semiclassical twotemperature model for ultrafast laser heating[END_REF][START_REF] Norman | Atomistic modeling of warm dense matter in the two-temperature state[END_REF][START_REF] Pisarev | Atomistic simulation of ion track formation in uo2[END_REF][START_REF] Seaton | Domain decomposition of the two-temperature model in dl poly 4[END_REF][START_REF] Daraszewicz | Determination of the electron-phonon coupling constant in tungsten[END_REF]. This combines two schemes where electrons are modeled as a continuum on a regular grid, and ions are tracked using a regular classical MD as:

m j ∂v j ∂t = -∇ j U (r 1 , . . . , r n )+F lang j (T e -T a )- ∇P e n i (1)
where m j and v j are respectively the mass and velocity of the j th atom, ∇ j U (r 1 , . . . , r n ) is the force produced by the empirical potential. P e is the electronic pressure, and n i is the ionic density. T a , T e represent the ionic and electronic temperature respectively. The t, r represent the considered time and position, where the F lang j is the random force due to the electronphonon coupling defined as:

F lang j (t) = -γ p m i v i + F (t) (2)
γ p represents the electron-phonon contribution coefficient. The F (t) is a random force giving the effect of noise satisfying the fluctuation dissipation theorem [START_REF] Phillips | A twotemperature model of radiation damage in -quartz[END_REF][START_REF] Mcclure | Role of electronic thermal transport in amorphous metal recrystallization: A molecular dynamics study[END_REF].

F (t) = 0 F (t ). F (t) = 2k B T e m i γ p δ(t -t) (3)
where k B is the Boltzmann constant, m i is the atomic mass of the used species, and δ(t -t) is the Dirac function.

The electronic subsystem is modeled by a one dimensional heat diffusion equation in the x direction with an external heat source:

C e (T e )
∂T e ∂t = ∇(κ e (T e )∇T e ) -g p (T e -T a ) + I(t)e -x/lp l p (4) where κ e is the electronic thermal conductivity, C e (T e ) the electronic specific heat, l p the penetration depth, g p electron-phonon coupling strength coefficient, and I(t) = I 0 Π(t). Where Π(t) is a rectangular function with Π(t) = 1 if 0 < t < 100 fs and 0 otherwise. The I 0 is the initial laser intensity. The equivalent absorbed laser fluence is defined as F abs = I.τ with τ the pulse duration. Furthermore, the friction due to the electron-ion coefficient is given by:

g p = 3N k B γ p ∆V m i (5)
∆V and N represent the volume of a unitary electronic grid and the total number of atoms inside each grid element respectively.

The C e (T e ) is considered electronic temperature dependent, this quantity is either obtained by using available experimental data (cf. Tab I) or by performing first principle calculation (cf. Tab S1) and fit the result according to Eq. [START_REF] Bressel | Femtosecond laser induced density changes in geo2 and sio2 glasses: fictive temperature effect[END_REF].

C e (X) = C 0 + (a 0 + a 1 X + a 2 X 2 + a 3 X 3 + a 4 X 4 )e -(AX) 2
(6) where C 0 , a 0 , a 1 , a 2 , a 3 , a 4 , A are the corresponding fitting parameters and X = T e /1000.

An electron pressure contribution P e (T e ) = BC e (T e )T e is expected to induce a blast force in the computation according to Eq. ( 1), with B a material-dependent constant. In the present study, we have neglected this surface effect and fix B = 0. For the crystalline structures the C e (T e ) are computed using density functional theory. However, for the amorphous structures, a free electron model is used as C e (T e ) = γT e [START_REF] Kittel | Introduction to Solid State Physics[END_REF]. The electronic specific heat coefficient γ of each composition, reliable experimental data mentioned in [START_REF] Choy | Thermal conductivity of amorphous alloys above room temperature[END_REF][START_REF] Garoche | Comparison between amorphous and crystalline phases of copper-zirconium alloys by specific-heat measurements[END_REF] are available. The values are summarized in table I.

The electron-phonon coupling factor g p is considered as constant in this model. However, a linear response approximation [START_REF] Marinier | Femtosecond laser ablation of cux zr1-x bulk metallic glasses: A molecular dynamics study[END_REF] relation is used to compute the electron-phonon coupling factor in both crystalline and amorphous states by the relation:

g p = xg Cu p +(1-x)g Zr p .
Where x is the percentage of Cu atoms present in a given Cu x Zr 1-x structure.

A proper descriptions of the samples thermodynamics are ruled by κ e (T e ). In fact, to feed the simulation model, experimental data are taken at 300 K. The thermal diffusion is related to the previous quantities by:

D e = κ e (T e ) C e (T e ) (7)

B. Samples elaboration

Two crystalline CuZr supercells are prepared using ATOMSK [START_REF] Hirel | Atomsk: A tool for manipulating and converting atomic data files[END_REF], the B 2 -Cu 50 Zr 50 , C11 b -Cu 33.3 Zr 66.7 of dimension 9.69×3.55×3.55 nm 3 and 3.56×3.56×12.35 nm 3 respectively. The energy of the B2 and C11 b are minimized and relaxed using the conjugate gradient algorithm [START_REF] Yuan | A modified polakribirepolyak conjugate gradient algorithm for large-scale optimization problems[END_REF]. The system is then replicated along the three directions of space up to a size of about 1092 × 15 × 15 nm 3 . Periodic Boundary Conditions (PBCs) are used during the whole preparation process. To cancel residual pressure in the samples during the relaxation, a force norm criterion of 10 -8 eV/ Å is used. Then, the system is thermalized up to 300 K using the NPT ensemble (5 ps) and further equilibrated at constant PT conditions for another 5 ps. To prepare both MGs samples namely α-Cu 50 Zr 50 and α-Cu 33.3 Zr 66.7 , the B 2 and C11 b structures are first equilibrated at 300 K in the NVE ensemble for 100 ps and then heated up to 2000 K above the theoretical melting point T m = 1207 K, 1401 K (±50 K) for Cu 50 Zr 50 , and Cu 33.3 Zr 66.7 respectively using NPT during 170 ps as described in the supplemental material [63] (see, also, references [START_REF] Kokalj | Xcrysdena new program for displaying crystalline structures and electron densities[END_REF][START_REF] Porter | Phase transformations in metals and alloys[END_REF][START_REF] Zhang | Molecular dynamics simulations of the melting curve of nial alloy under pressure[END_REF][START_REF] Gonze | The abinit project: Impact, environment and recent developments[END_REF][START_REF] Torrent | Implementation of the projector augmentedwave method in the abinit code: Application to the study of iron under pressure[END_REF] therein). Annealing is performed at 2000 K for 100 ps to erase the crystalline memory of the samples. This process is followed by a fast quench run from 2000 K to 300 K at a cooling rate of 10 12 K.s -1 . The obtained glassy states are finally equilibrated at ambient temperature of 300 K for extra 100 ps to remove internal stress within the final structure. The MG are replicated along the 3 directions of space before a last 300 K equilibration run of 10 ps. The damping coefficients of 0.1 ps and 1.0 ps in Nose-Hoover (thermostatbarostat) style are used respectively for temperature and pressure [START_REF] Evans | The nose-hoover thermostat[END_REF]. The MG samples final size is about 874 × 7 × 13 nm 3 . Finally, for both MGs and crystalline structures the half-part along x of the samples is suppressed (i.e., the matter stands in the [L x /2,L x ] with L x the box length along x) to introduce the empty space necessary to accommodate the laser-matter interaction and possible matter ablation.

C. Laser-matter interaction setup

The heat diffusion Eq. ( 4) is solved numerically on a 702×1×1 regular electronic cubic mesh as shown in figure 1. The laser shot attenuation length are set according to the value shown in table I. They are computed using the formula described in [START_REF] Marinier | Femtosecond laser ablation of cux zr1-x bulk metallic glasses: A molecular dynamics study[END_REF] as l p = xl Cu p + (1 -x)l Zr p , where l Cu p = 12.2 nm, l Zr p = 16.4 nm, are penetration depth of Cu, Zr, and x is the percentage of Cu atoms present in each composition. A laser wavelength of 800 nm and a step function dx ≈ 3 nm for the T e variation is used. Here, PBCs are set in y and z lateral directions and free boundaries are imposed following the laser shot x = [100]. Ultrafast laser irradiation induces a strong temperature gradient normal to the surface which leads to the generation of a pressure wave (P-wave) propagating to the bulk of the target [START_REF] Shugaev | Thermoelastic modeling of laser-induced generation of strong surface acoustic waves[END_REF][START_REF] Shugaev | Laser-induced thermal processes: Heat transfer, generation of stresses, melting and solidification, vaporization, and phase explosion[END_REF]. Hence, a trial-and-error approach focused on the gradual extinction of the reflected P-wave is used. For this purpose, we define non reflecting boundary conditions (NRBC), made of a viscous damping region [74-76] with dimensions of L NRBC = 100 nm. A damping factor γ damp = 2.5 10 -3 eV ps Å-2 is used in this work in all simulations. The ultrafast process is modeled as a homogeneous dose energy deposition on the surface with a pulse duration of a 100 fs. The TTM-MD stage pushes the system to an extreme non equilibrium state which is held for 30 ps. The second stage of the simulation is a regular MD run in the NVE ensemble for (i) 170 ps for crystals and (ii) 970 ps for MGs.

III. RESULTS AND DISCUSSION

Here we investigate the influence of a single ultrafast laser pulse on B 2 and C11 b bulk crystalline CuZr as well as on two amorphous phases using atomistic simulations. The pulse duration is 100 fs long and the absorbed fluence is 160 mJ/cm 2 . After the energy deposition on the surface, the samples thermalization is achieved by a strong electron-phonon coupling where both electrons and the ions reach the same equilibrium temperature. The final temperature gradient defines the position of the melting front below the surface.

A. Irradiation of Cu50Zr50: B2 crystal and amorphous

Figure 2(a) shows the evolution of the B 2 microstructure during the simulation. The local crystalline structure is investigated by applying the PTM method to the alloyed system without atomic species differentiation which leads to characterize the B 2 perfect crystal 4). There is no energy transfer with the vacant space beyond the top surface. Additionally, the rear region of the sample represent the NRBC region.

Table I. Physical properties of crystalline and amorphous CuZr binary systems, namely specific heat constant γ and electron-phonon factor gp, electronic thermal conductivity κe, and electronic thermal diffusion De. The transformed domain takes the shape of two nanoblades of ∼ 8.5 nm thickness each that propagate across the sample due to PBCs. The deformation gradient in the x direction displayed in figure 2(b) shows that the nano-blades follow a distribution with a spatial pseudoperiodicity separated by a thin B 2 inter-layer. A closer look to the transformed structure allows the identification of two simple cubic sub-lattices attributed to Cu and Zr as in the original B 2 phase. However, the atomic bonding changes from aligned to wavy or zigzag configurations along certain directions (see figure 2(d)) typical of the B 19 crystalline structure (e.g., AuCd). This phase transition is similar to the martensitic transforma-tion observed in B 2 NiTi alloys [START_REF] Zhang | Crystallography of the b2 r b19 phase transformations in niti[END_REF][START_REF] Ko | Development and application of a ni-ti interatomic potential with high predictive accuracy of the martensitic phase transition[END_REF][START_REF] Bowles | The crystallography of martensite transformationsi[END_REF]. The relative orientation between the two nano-blades allows to confirm the nucleation and propagation of two {111}-oriented B 19 nanotwins.

Property ρ0 gp γ κe De lp Unit [g cm -3 ] [10 17 Wm -3 K -1 ] [J m -3 K -2 ] [W m -1 K -1 ] [10 -5 m 2 s -1 ] [nm] B2-
The cooperative displacement of atoms requires a significant amount of energy that here relies on the laser pulse. Indeed, the pulse emission is responsible of local temperature gradients which induce the shock wave (P-wave) and then the aforementioned phase transformation. The pressure and temperature profiles timeevolution in the transformed zone are shown in figure 3. The pressure profile is characterized by a loading (pressure wave) and unloading (rarefaction wave) P-waves. At short times (t = 0-60 ps), the front of the loading P-wave reaches the emphasized region with a maximal value of ∼ 7 GPa but no transformation is noticed. The B 2 to B 19 phase transformation that starts at t = 80 ps exactly coincides with the appearance of the rarefaction wave (∼ -5 GPa extremum) as the loading wave progressively reflects to the colder material. Hence, the phase transition is not due to the propagation of the pressure wave but more precisely to the unloading one that has exceeded a critical level to activate the transformation [START_REF] Sutrakar | Single and multistep phase transformation in cuzr nanowire under compressive/tensile loading[END_REF]. Meanwhile, the temperature profile shows an increase up to ∼ 1000 K below the melting temperature T m = 1207 K. A mean T of ∼ 600 K is noticed on the whole transformed zone, confirming that the solid-solid phase transformation observed is driven by mechanical effects as observed in figure 2(c).

An equivalent protocol was applied to the Cu-Zr amorphous (α-Cu 50 Zr 50 ) sample having the same stoichiometry. For comparison purposes with the B 2 crystal, identical laser conditions are used while the total simulation time is extended up to 1 ns. Figure 4(a) shows the structural evolution of the α-Cu 50 Zr 50 sample after the pulse. As in the B 2 case, results show surface liquefaction concomitant to the nucleation of voids preceding an ablation event. However, no phase transformation is observed. Regarding the CuZr phase diagram [START_REF] Arias | Cu-zr (copper-zirconium)[END_REF], local transformation into B 2 is possible and segregation could drive the formation of single crystalline inclusions (with FCC or HCP crystal structure). However, the nucleation of more complex phases is completely excluded. We conclude here that α-Cu 50 Zr 50 laser-induced local phase transition can not be obtained using the same laser pulse conditions (duration and fluence) than the one used to phase-transform B 2 CuZr crystals. Figure 5 shows the structural evolution of the C11 b -Cu 33.3 Zr 66.7 sample. In particular, it provides the atomistic picture of the point defects generated in this case by the ultrafast pulse. First in the t = 0 -100 ps range, a liquid front is formed characterized by the coalescence of small bubbles formation process that generates voids below the surface in the t = 100 -160 ps range. However, the ablation regime only starts at t ∼ 160 ps causing a surface expansion of about 60 nm thickness which expands freely in the vacuum region. Figure 5(b) emphasizes a small area of subsurface material oriented along x =[100]. Using this projection, the classical C11 b stacking along z =[001] shows some disturbances and the original [100]-oriented Cu and Zr columns of atoms observed at 0 K can now be schematized by pseudo-columnar domains separated by a thinner layer of atoms attributed to projected elasto-plastic deformation (that has not to be confused with sub-lattice degradation or melting). In addition, point defects are identified up to x ∼ 740 nm.

Various kinds of point defects including interstitial as well as substitutional Cu and Zr atoms are noticed. Atoms in defective environment (point defects and highly-sheared regions) are identified using a criterion on the atomic displacement along y and z, δy and δz respectively (taking t = 0 as a reference). The displacement criterion is chosen as about half the Cu-Zr or Zr-Zr blond length along y, ensuring to avoid counting atoms fluctuating around their equilibrium position. We choose the y and z directions since they are less impacted by the P-wave induced homogeneous matter displacement merely observed along the x direction. While the amount of atoms with significant δz are only few (0.03 % of the atoms in the emphasized zone at t = 200 ps), the number of atoms with δy > 1 Å continuously increases during the simulation (Figure 5(c)). The anisotropy observed between δy and δz corroborates the observed shear zones in the z = [001] plane. Indeed, no defects are detected in the emphasized domain of the material from t = 0 -25 ps, the fraction of atoms in defective configuration grows up to a maximum 1.12 % at t = 175 ps before decreasing to 0.93 % at the end of the simulation due to relaxation. The defects concentration were computed by evaluating N def ects /N total in the total system volume region (241984 nm 3 ), where N def ects (69985 atoms) and N total (6248652 atoms) are the number of defects and the number of atoms, respectively.

To better understand this phenomenon, we again look at the PT profiles as shown in Figure 6. Figure 6(a) shows that a strong P-wave is generated and reaches a maximum value of ∼ 8 GPa at t = 20 ps. At t = 60 ps the rarefaction P-wave propagates with a negative extremum front of ∼ -6 GPa appears. The rarefaction wave and the defects population are highly-correlated in terms of space and time. The lattice temperature shown in figure 6(b) confirms that at the liquid-solid interface the temperature is ∼ 2000 K, higher than the melting point T m = 1401 K. We have to mention that below the solid-liquid interface the system is highly perturbed due to heat diffusion phenomenon but not melted. Here we show that the electron-phonon coupling plays an important role in the laser-matter interaction, and since the coupling factor is stronger for C11 b -Cu 33.3 Zr 66.7 than for I) the temperature and P-wave maximum are higher in the case of C11 b -Cu 33.3 Zr 66.7 .

Finally, the amorphous α-Cu 33.3 Zr 66.7 sample is subjected to a 1 ns simulation with same fluence and pulse duration. Results are shown in figure 4(b). A spallation process is noticed at t = 100 ps through the growth of subsurface voids inducing the ejection of a top liquid layer but again without any sign of crystal nucleation. One can hypothesize that the absence of crystallization is due (i) either to the tensile rarefaction wave magnitude sufficient enough to induce spallation in the molten part but not to reorganize atoms in an ordered manner or the (ii) a non-favorable stoichiometry.

IV. CONCLUSION

The local subsurface structural phase transition and defects formation process attributed to ultrafast laser shot in crystalline CuZr systems are predicted using TTM-MD atomistic calculations. In summary, using the same laser conditions, two remarkable outcomes are reported here depending on the initial crystallographic structures. The original B 2 crystalline sample phase-transform into twinned B 19 while C11 b shows a huge deformation pattern including points defects and shear zones. Both transformations are driven by the rarefaction pressure wave. They are permanent as the materials reach metastable states in our calculations and localized within the subsurface. Therefore, they can be captured experimentally. Concerning the B 2 structure, the simulation corroborates that the phase transformation and the twinning process emerge in a quasi-ordered manner, exhibiting no spatial variability across the crystalline structure. In contrast, the response of the C11 b structure is different as the defects generation is noticed instead of a martensite phase transition. On the other hand, the response of the irradiated MGs counterparts show a completely different structural response. Indeed, no crystallization is noticed for the same laser conditions. Beside the surface ablation in the MGs, we might have expected the laser shock to induce recrystallization by playing a similar role as a nanoindenter [START_REF] Sort | Out-of-plane magnetic patterning based on indentation-induced nanocrystallization of a metallic glass[END_REF][START_REF] Fornell | Enhanced mechanical properties in a zr-based metallic glass caused by deformation-induced nanocrystallization[END_REF][START_REF] Ma | Cyclic stress induced surface nanocrystallization adjacent to indentation edge of zr-based bulk metallic glass at room temperature[END_REF].

Targeting innovative technological applications requires the control of local crystalline structures properties at the nanoscale that can be promoted by lasermatter interactions. Dislocations and twins that appear at the nanoscale as a strain relaxation mechanism can be frozen within a solid alloy. As predicted in the B 2 structure, this leads to the emergence of advantageous mechanical properties [START_REF] Zhang | Dislocation mechanisms and 3d twin architectures generate exceptional strength-ductility-toughness combination in crconi medium-entropy alloy[END_REF]. Nano-sized twins are known to enhance hardness [START_REF] Tian | Ultrahard nanotwinned cubic boron nitride[END_REF] and nano-twinned materials become more resilient to fracture and fatigue crack initiation during plastic deformation [START_REF] Uttam | Nan-otwinning: Generation, properties, and application[END_REF]. Nanoarchitecturing of the C11 b structure with defects (points defects, dislocations, stacking faults) could be could be exploited for battery applications [START_REF] Pan | Point defect of titanium sesquioxide ti2o3 as the application of next generation li-ion batteries[END_REF], elaboration of high-performance thermoelectric materials [START_REF] Hu | Point defect engineering of high-performance bismuth-telluride-based thermoelectric materials[END_REF] or improving the fracture toughness [START_REF] Liu | Physical metallurgy and mechanical properties of transition-metal laves phase alloys[END_REF]. The originality of this work is that ultrafast laser irradiation could provide an effective method to tune a new generation of nanoarchitectured CuZr alloys which might open new routes future industrial applications.

Figure 1 .

 1 Figure 1. Simulation setup geometry used to model ultrafast laser pulse interaction in CuZr alloys using TTM-MD simulations. Atoms are located inside the rectangular box and the overlying red grid represents the electronic cell domain for the finite difference solution of Eq. (4). There is no energy transfer with the vacant space beyond the top surface. Additionally, the rear region of the sample represent the NRBC region.

  structure and transformed domains environments (liquid, amorphous or new phase). Just after the laser pulse energy deposition, the upper part of the samples instantaneously melts as shown by the large amount of atoms in defective configurations. At t = 80 ps, several voids and bubbles form below the CuZr surface at position x = 572 nm and expand rapidly reaching their critical size at t = 200 ps. In the mean time, a new crystalline phase nucleates and grows at position x ∼ 700 nm far deep in the sample starting from t = 60 ps as shown in figures 2(c) and (d).

Figure 2 .

 2 Figure 2. Snapshots of the atomic configurations evolution in figure a) at several times : 0 ps, 20 ps, 60 ps, 80 ps, 100 ps, 120 ps, 160 ps, and 200 ps in B2-Cu50Zr50 target during ultrafast laser irradiation of τ = 100 fs pulse duration and F abs = 160 mJ/cm 2 absorbed fluence. The atoms are colored according to their local atomic structure of the binary alloy by using the PTM, the double dashed arrow in the upper part of target shows the liquefied surface. b) Deformation gradient in the transformed zone following the x direction at different times computed by the method described in [77]. c) Zoom on the dashed region where the PTM algorithm is applied on each sub-lattice monoatomic lattice showing twinned martensite transformation. d) Enlarged view of the chemical species of the new twinned-B19 phase emergence after the irradiation.

Figure 3 .

 3 Figure 3. Spatiotemporal evolution of thermodynamic quantities in the B2-Cu50Zr50 transformed region located in subsurface (zoom region) area target: a) pressure, and b) lattice temperature. The initial surface position is given by x surface = 556.68 nm following the x axis direction.

Figure 4 .

 4 Figure 4. Snapshots of the structural evolution in the amorphous structures: a) α-Cu50Zr50, and b) α-Cu33.3Zr66.7 using TTM-MD method at various time steps: 0 ps, 100 ps, 300 ps, 400 ps, 700 ps, and 1 ns. The atoms are colored according to their local atomic structure using PTM algorithm, two atomic contrasts are shown, the zone delimited by a double dashed arrow in the upper part of targets correspond to liquid phase and the inner part is the amorphous phase. The irradiation conditions are τ = 100 fs pulse duration and F abs = 160 mJ/cm 2 of absorbed fluence.

Figure 5 .B 2 -

 52 Figure 5. Snapshots of atomic configurations evolution in figure a) at several times: 0 ps, 20 ps, 60 ps, 80 ps, 100 ps, 120 ps, 160 ps, and 200 ps in the C11 b -Cu33.3Zr66.7 sample target during ultrafast laser irradiation of τ = 100 fs pulse duration and F abs = 160 mJ/cm 2 of absorbed fluence, the double dashed arrow in the upper part of target shows the liquefied surface. The atoms are colored according to their local atomic structure using PTM algorithm. b) Enlarged view of interstitial and substitutional defects formation in subsurface within C11 b -Cu33.3Zr66.7 target after the irradiation at different times. c) Temporal defects ratio evolution in Cu33.3Zr66.7 (C11 b ) target sorted according to the displacement vector δy > 1 Å. Atoms satisfying this condition are expected to be in defective positions, this is achieved using displacement vectors analysis [83].

Figure 6 .

 6 Figure 6. Spatiotemporal evolution of thermodynamic quantities in the subsurface transformed region (zoom region) : a) pressure, and b) lattice temperature evolution in (C11 b -Cu33.3Zr66.7 sample target. The initial surface position is given by x surface = 547.30 nm following the x axis direction.

ACKNOWLEDGEMENTS

This work is supported by the ANR project MEGALIT (ANR-18-CE08) and the project FORMEL of the Pack Ambition Research program of the Auvergne Rhne-Alpes Region. Numerical calculations have been performed using HPC resources from GENCI-TGCC, CINES (project gen7041) and P2HPD computing center facility of Lyon university.


E Berlioz
email: emilie.berlioz@univ-poi�ers.fr

D S Kostopoulos

C Blondel

G Merceron

St-Vallier Senèze

Krimni Gerakarou

C Ramosus

A Gallicus

E Ctenoides

M Rhenanus

Exploring the dietary plas�city of Eucladoceros ctenoides: what does its occurrence in the fossil record mean?

F. U .

St-Vallier

Gerakarou

Results

The dental microwear texture of E. ctenoides show an appreciable varia�on similar to C. elaphus [4], reflec�ng the diversity of the food categories this taxon is able to consume (fig. 3).

Varia�ons of the dental microwear textures of E. ctenoides through �me reflects structural changes in ecosystems that miror early Pleistocene clima�c oscilla�ons (fig. 4).

Introduc�on

Cervidae are o�en used as paleoenvironmental proxies. In par�cular they are used for the early Pleistocene paleoenvironmental reconstruc-�ons as they were par�cularly abundant and diversified during this period of clima�c oscilla�ons and Homo dispersal in Europe [1]. This family has been used for long in the framework of environmental reconstruc�ons but the presence of deer in a fossil locality is tradi�onally interpreted as witnessing a significant tree cover in the habitat. Among cervids, the large deer Eucladoceros ctenoides is present throughout the whole Europe since ca. 2.5 Ma (middle Villafranchian) and presents a body mass similar to that of extant European Cervus elaphus. Considering the large body mass of E. ctenoides we expected this taxon to present a significant dietary plas�city and therefore to cons�tute an efficient paleoenvironmental indicator.

Conclusions

It is the first comprehensive explora�on through �me and space of the dietary plas�city of Eucladoceros ctenoides, an emblema�c deer species of the early Pleistocene. E. ctenoides presents a significa�ve trophic plas�city. Such a result is not surprising when compared with the plas�city in dietary habits of modern C. elaphus. However, there are at odds with tradi�onal paleoenvironmental interpreta�ons. In our study, we show that the feeding ecology of E. ctenoides, rather than its occurence in a fossil assemblage, is an efficient proxy for paleoenvironmental reconstruc�ons.

Method

We performed a polyvinylsiloxane mold of the dental facet (fig. 2-a). The topology of the surface is obtained by scanning the dental mold with a DCM8 white light confocal surface profilometer (obj. x100; NA=0.90; working distance = 0.9 mm)(fig. 2-b). Dirty surfaces are discarded. (fig. 2-c) We pre-processed the textural surfaces (fig. 2-d) before the acquisi�on of the textural parameters (SSFA).

  Material fossils : we analysed the teeth of 145 E. ctenoides specimens from 8 early Pleistocene European locali�es (fig. 1-a) (sample sizes: St-Vallier=22, Chilhac=78, Varshets=21, Dafnero=3, Senèze=11, Slivnitsa=1, Gerakarou=9, Krimni=1). reference : we compared our results with 2 modern Cervus elaphus popula�ons : -grazing C. elaphus from Lugar-Nuevo, south of Spain [2] (fig. 1-b), -browsing C. elaphus from Białowieża, polish primeval forest [3] (fig. 1 c).

Fig. 1 :Fig. 2 :

 12 Fig. 1: fossil locali�es with Eucladoceros ctenoides. (a): geographic localiza�on of locali�es. Lugar-Nuevo (b) and Białowieża (c) are reference locali�es with extant Cervus elaphus of known diet and environment (illustrated at the bo�om of the map). All locali�es belong to the early Pleistocene (d) and to middle or late Villafranchian. Fig. 2: Acquisi�on of textural parameters.

Fig. 3 : 4 :

 34 Fig. 3: Distribu�on (mean and standard error of mean) of E. ctenoides depending on the complexity and anisotropy of dental facets. The distribu�on of all reference specimens is illustrated by transparent green crosses (Białowieża) and red circles (Lugar-Nuevo). Fig. 4: Evolu�on through �me of the complexity, heterogeneity (9 cells) of the complexity and textural fill volume of E. ctenoides. Grey lines are for western european locali�es and black lines for eastern ones. Significa�ve differences are illustrated by red stars.

  Fig. 3