Wafa Khrouf
email: wafa.khrouf@supcom.tn

Zeineb Hraiech
email: zeineb.hraiech@supcom.tn

Fatma Abdelkefi
email: fatma.abdelkefi@supcom.tn

Mohamed Siala
email: mohamed.siala@supcom.tn

Matthieu Crussière
email: matthieu.crussiere@insa-rennes.fr

On the Joint Use of Time Reversal and POPS-OFDM for 5G Systems

Keywords: Waveform Optimization, Waveform Design, Pulse Shaping, POPS-OFDM, Signal to Interference plus Noise Ratio (SINR), Time Reversal

This paper investigates the efficiency of the combination of the Ping-pong Optimized Pulse Shaping-Orthogonal Frequency Division Multiplexing (POPS-OFDM) algorithm with the Time Reversal (TR) technique. This algorithm optimizes the transmit and receive OFDM waveforms with a significant reduction in the system Inter-Carrier Interference (ICI)/Inter-Symbol Interference (ISI) and guarantees maximal Signal to Interference plus Noise Ratio (SINR) for realistic mobile radio channels in 5G Systems. To this end, we characterize the scattering function of the TR channel and we derive the closedform expression of the SINR as a Generalized Rayleigh Quotient. Numerical analysis reveals a significant gain in SINR and Out-Of-Band (OOB) emissions, brought by the proposed TR-POPS-OFDM approach.

I. INTRODUCTION

OFDM modulation has witnessed a considerable interest from both academia and industry. This is due to many advantages such as low complexity receiver, simplicity and efficient equalization structure. Indeed, it has been adopted in various wired and wireless standards, such as ADSL, DVB-T, Wimax, WiFi, and LTE. Nevertheless, in its present form, OFDM presents several shortcomings and as such it is not capable of guaranteeing the quality of service of new and innovative applications and services that will be brought by 5G systems. In fact, it has a high spectral leakage and it requires strict frequency synchronization because it uses a rectangular waveform in time, leading to significant sidelobes in frequency [START_REF] Yunzheng | A Survey: Several Technologies of Non-Orthogonal Transmission for 5G[END_REF]. As a consequence, any lack of perfect frequency synchronization, to be expected from most of the innovative Machine Type Communications (MTC) in 5G, causes important Inter-Carrier Interferences (ICI). In addition to that, a variety of services and new applications will be provided by 5G systems, such as high data-rate wireless connectivity, which requires large spectral and energy efficiency, and Internet of Things (IoT), requiring robustness to time synchronization errors [START_REF] Luo | Signal Processing for 5G: Algorithms and Implementations[END_REF].

In order to overcome OFDM limitations and meet 5G requirements, various modulations have been suggested in the literature such as Generalized Frequency Division Multiplexing (GFDM), Universal Filtered Multi-Carrier (UFMC) and Filter Bank Multi-Carrier (FBMC), which are proposed in 5GNOW project [START_REF] Wunder | 5GNOW: Non-Orthogonal, Asynchronous Waveforms for Future Mobile Applications[END_REF]. It is shown at [START_REF] Wunder | 5GNOW: Non-Orthogonal, Asynchronous Waveforms for Future Mobile Applications[END_REF] that GFDM offers high flexibility for access to fragmented spectrum and low Out-Of-Band (OOB) emissions. However, in contrast to UFMC, GFDM has a low robustness to frequency synchronization errors in the presence of Doppler spread. Moreover, like UFMC [START_REF] Schaich | Waveform Contenders for 5G -OFDM vs. FBMC vs. UFMC[END_REF], FBMC has a high spectral efficiency and a good robustness to ICI. Nevertheless, FBMC, because of its long shaping filters, cannot be used in the case of low latency, sporadic traffic and small data packets transmission. Furthermore, authors in [START_REF] Siala | Novel Algorithms for Optimal Waveforms Design in Multicarrier Systems[END_REF] and [START_REF] Hraiech | POPS-OFDM: Ping-pong Optimized Pulse Shaping-OFDM for 5G Systems[END_REF] propose new class of waveforms, namely POPS-OFDM, which iteratively maximize the SINR in order to create optimal waveforms at the transmitter (Tx) and the receiver (Rx) sides. The obtained waveforms are well localized in time and frequency domains and they are able to reduce the ISI and the ICI as they are not sensitive to time and frequency synchronization errors. Another alternative which aims to reduce interference especially in time, caused by highly dispersive channels, is the time reversal (TR) technique which has been recently proposed for wireless communications systems. Its time and space focusing properties make it an attractive candidate for green and multiuser communications [START_REF] Dubois | Performance of Time Reversal Precoding Technique for MISO-OFDM Systems[END_REF]. In fact, it reduces ISI at the Rx side and it mitigates the channel delay spread.

This paper aims to design new waveforms for 5G systems by combining the benefits of POPS and TR techniques in terms of interference resilience. For this end, we analyze the corresponding system and we derive the SINR expression. We also evaluate the performances of the proposed approach in terms of SINR and OOB emissions.

The remaining of this paper is organized as follows. In Section III, we present the system model. In Section IV, we focus on the derivation of the SINR expression for TR systems and describe the TR-POPS-OFDM algorithm for waveforms design. Section V is dedicated to the illustration of the obtained optimization results and shed light on the efficiency of the proposed TR-POPS-OFDM approach. Finally, Section VI presents conclusion and perspectives to our work.

II. NOTATIONS

Boldface lower and upper case letters refer to vectors and matrices, respectively. The superscripts . * and . T denote the element-wise conjugation and the transpose of a vector or matrix, respectively. We denote by v = (. . . , v -2 , v -1 , v 0 , v 1 , v 2 , . . .) T = (v q ) q∈Z = (v q ) q the infinite vector v. The last notation, (v q ) q , with the set of values taken by q is not explicitly specified, that means that q spans Z.

Let M = (M pq ) p∈Z,q∈Z = (M pq ) pq refers to the infinite matrix M. The matrix shift operator Σ k (•) shifts all matrix entries by k parallel to the main diagonal of the matrix, i.e. if M = (M pq ) pq is a matrix with (p, q)-th entry M pq , then Σ k (M) = (M p-k,q-k ) pq . The symbol ⊗ is the convolution operator of two vectors and the symbol is the componentwise product of two vectors or matrices. We denote by E the expectation operator and by |.| the absolute value.

III. SYSTEM MODEL

In this Section, we first present the TR principal. Then, we describe the channel and system models in which we will apply our approach.

A. Time Reversal Principle

The Time Reversal (TR) principle [START_REF] Fink | Time Reversal of Ultrasonic Fields -Part I: Basic Principles[END_REF], [START_REF] Lerosey | Time Reversal of Electromagnetic Waves[END_REF], comes from the acoustic research field and allows a wave to be localized in time and space. Such a technique can be exploited to separate users, addressed simultaneously on the same frequency band, by their different positions in space.

The use of TR in transmission systems has generated a particular excitement as it allows, from a very high temporal dispersion channel, to obtain an ideal pulse in time and in space. This property has several useful advantages in wireless communications, among which we cite the followings:

• Negligible or null ISI brought a nearly "no memory" equivalent channel. • Minimum inter-user interference thanks to space power localization, with received negligible power outside a focal spot targeted to a given Rx. • Physical-layer-secured data transmission towards a desired user, as other users located outside of the focal spot of the targeted user will receive only a few power. TR integration into a telecommunication system is very simple. It consists in applying a filter on the transmitted signal. We suppose that we have a perfect knowledge of the transmission channel and that it is invariant between the instants of its measurements and the application of the TR at the Tx side. This filter is made up of the Channel Impulse Response (CIR) reversed in time and conjugated. It has the form of a matched filter to the propagation channel which guarantees optimal reception in terms of Signal to Noise Ratio (SNR). Then, the transmitted signal will cross an equivalent filter equal to the convolution between the channel and its time reversed version.

B. Channel Model

We consider a Wide Sense Stationary Uncorrelated Scattering (WSSUS) channel in order to have more insights on the TR-POPS-OFDM performances in the general case. To simplify the derivations, we consider a discrete time system. We denote by T s the sampling period and by R s = 1 Ts the sampling rate. We suppose that the channel is composed of K paths and that the Tx has a perfect knowledge of the channel state at any time. Note that this hypothesis is realistic in the case of low Doppler spread. Let

h (p) = (h (p) 0 , h (p) 1 , . . . , h (p) K-1 ) T be the channel discrete version at instant p such as h (p) l = M -1 m=0
h lm e j2πν lm pTs is the path corresponding to a delay lT s , where M is the number of Doppler rays, h lm and ν lm denote respectively the amplitude and the Doppler frequency of the l th path and the m th Doppler ray. The ray amplitudes, h lm , are supposed to be centered, independent and identically distributed complex Gaussian variables with average powers

π lm = E[|h lm | 2 ]. We denote by π l = M -1 m=0 π lm , where K-1 l=0 π l = 1.
The channel time reversed version at the instant p can be written as:

g (p) = h (p) * K-1 , . . . , h (p) * 1 , h (p) * 0 T . (1)
When we apply the TR technique at the Tx in a Single Input Single Output (SISO) system, the equivalent channel, experienced by the transmission at time instant pT s at the Rx, could be seen as the convolution between the channel and its time reversed version, as follows:

H (p) = h (p) ⊗ g (p) (2) = H (p) -(K-1) , . . . , H (p) -1 , H (p) 0 , H (p) 1 , . . . , H (p) K-1 T ,
where

H (p) k = K-1-|k| l=0 M -1 m,m =0 f (k, l, m, m ), such as: f (k, l, m, m ) = h * lm h l+k,m e j2π(ν l+k,m -ν lm )pTs , if k ≥ 0 h lm h * l-k,m e -j2π(ν l-k,m -ν lm )pTs , else.
It should be noted that

H (p) -k = H (p) * k
, which means that the channel is Hermitian symmetric, and that the equivalent aggregate channel coefficients, H (p) k , are still decorrelated, as in the actual channel.

C. OFDM System

In this paper, we consider a discrete time version of the waveforms to simplify the theoretical derivations that will be investigated.

Let T and F refer to the OFDM symbol duration and the frequency separation between two adjacent subcarriers respectively. The sampling period is equal to T s = T N where N ∈ N. We denote by δ = 1 F T = Q N the time-frequency lattice density, where Q = 1 TsF ≤ N is the number of subcarriers. We denote by e = (e q ) q the sampled version of the transmitted signal at time qT s , with a sampling rate R s = 1 Ts , expressed as: e = m,n a mn ϕ ϕ ϕ mn ,

where ϕ ϕ ϕ mn = (ϕ q-nN ) q (e j2πmq/Q ) q is the time and frequency shifted version of the OFDM transmit prototype waveform, ϕ ϕ ϕ = (ϕ q ) q , used to transmit the symbol a mn . We suppose that the transmitted symbols are decorrelated, with zero mean and energy equal to

E = E[|a mn | 2 ] ϕ ϕ ϕ 2 .
The received signal is expressed as:

r = mn a mn φ φ φmn + n, (4)
where

[ φ φ φmn ] q = K-1 k=-(K-1) H (q)
k [ϕ ϕ ϕ mn ] q-k is the channel distorted version of ϕ ϕ ϕ mn and n = (n q ) q is a discrete complex Additive White Gaussian Noise (AWGN), with zero mean and variance N 0 .

The decision variable, denoted Λ kl , on the transmitted symbol a kl is obtained by projecting r on the receive pulse ψ ψ ψ kl , such as:

Λ kl = ψ ψ ψ kl , r = ψ ψ ψ H kl r, (5)
where ψ ψ ψ kl = (ψ q-lN ) q (e j2πkq/Q ) q is the time and frequency shifted version of the OFDM receive prototype waveform ψ ψ ψ = (ψ q ) q and •, • is the Hermitian scalar product over the space of square-summable vectors.

IV. TR-POPS ALGORITHM

The main objective of this part is to optimize the waveforms at the Tx/Rx sides in our system based on TR technique. To this end, we adopt the POPS-OFDM principal [START_REF] Siala | Novel Algorithms for Optimal Waveforms Design in Multicarrier Systems[END_REF]. This algorithm consists in maximizing the SINR for fixed synchronization imperfections and propagation channel.

Without loss of generality, we will focus on the SINR evaluation for the symbol a 00 . Referring to (5), the decision variable on a 00 can be written as:

Λ 00 = a 00 ψ ψ ψ 00 , φ φ φ00 + (m,n) =(0,0)
a mn ψ ψ ψ 00 , φ φ φmn + ψ ψ ψ 00 , n and it is composed of three terms. The first term is the useful part, the second term is the ISI and the last term presents the noise term. Their respective powers represent useful signal, interference and noise powers in the SINR and which we derive their closed form expressions in the sequel. This SINR will be the same for all other transmitted symbols.

A. Average Useful, Interference and Noise Powers

The useful term is denoted U 00 = a 00 ψ ψ ψ 00 , φ φ φ00 . For a given realization of the channel, the average power of the useful term can be written as:

P h S = E ϕ ϕ ϕ 2 | ψ ψ ψ 00 , φ φ φ00 | 2 .
Thus, the useful power average over channel realizations is given by:

P S = E P h S = E ϕ ϕ ϕ 2 E | ψ ψ ψ 00 , φ φ φ00 | 2 . ( 6
)
The interference term, I 00 = (m,n) =(0,0) a mn ψ ψ ψ 00 , φ φ φmn , results from the contribution of all other transmitted symbols a mn , such as (m, n) = (0, 0).

For a given realization of the channel, the average power of the interference term can be written as:

P h I = E ϕ ϕ ϕ 2 (m,n) =(0,0) | ψ ψ ψ 00 , φ φ φmn | 2 .
Therefore, the interference power average over channel realizations has the following expression:

P I = E P h I = E ϕ ϕ ϕ 2 (m,n) =(0,0) E | ψ ψ ψ 00 , φ φ φmn | 2 , ( 7
)
where

E | ψ ψ ψ 00 , φ φ φmn | 2 = ψ ψ ψ H E φ φ φmn φ φ φH mn ψ ψ ψ. (8)
In the sequel, we consider a diffuse scattering function in the frequency domain, with a classical Doppler spectral density, decoupled from the dispersion in the time domain. So,

E [ φ φ φmn ] p φ φ φH mn q = K-1 k=-(K-1) Π k J 2 0 (πB d T s (p -q)) ×[ϕ ϕ ϕ mn ] p-k [ϕ ϕ ϕ mn ] * q-k
, where B d is the Doppler Spread, J 0 (•) is the Bessel function of the first kind of order zero and

Π k =          K-1 l=0 π 2 l + K-1 l=0 π l 2 , if k = 0 K-1-|k| l=0 π l π l+|k| , else
is the average power of the global channel. Then the average useful and interference powers have the following expressions:

P S = E ϕ ϕ ϕ 2 ψ ψ ψ H KS ϕ ϕ ϕ ψ ψ ψ and P I = E ϕ ϕ ϕ 2 ψ ψ ψ H KI ϕ ϕ ϕ ψ ψ ψ, (9)
where KS ϕ ϕ ϕ and KI ϕ ϕ ϕ are Hermitian, symmetric, positive and semidefinite matrices:

KS ϕ ϕ ϕ =   K-1 k=-(K-1) Π k Σ k ϕ ϕ ϕϕ ϕ ϕ H   Λ (10)
and

KI ϕ ϕ ϕ =     n ΣnN   K-1 k=-(K-1) Π k Σ k ϕ ϕ ϕϕ ϕ ϕ H     Ω   -KS ϕ ϕ ϕ . (11)
The entries of matrices Λ and Ω are defined as:

Λ pq = J 2 0 (πB d T s (p -q))
and

Ω pq = QJ 2 0 (πB d T s (p -q)) , if (p -q) mod Q = 0, 0, else, with p, q ∈ Z.
The noise term is given by N 00 = ψ ψ ψ 00 , n n n . Thus, the noise power average is the following:

P N = E | ψ ψ ψ 00 , n | 2 = ψ ψ ψ H E nn H ψ ψ ψ.
As the noise is supposed to be white, its covariance matrix is equal to R nn = E nn H = N 0 I, where I is the identity matrix. Consequently,

P N = N 0 ψ ψ ψ 2 . (12)

B. Optimization Technique

The SINR expression is the following:

SIN R = P S P I + P N = ψ ψ ψ H KS ϕ ϕ ϕ ψ ψ ψ ψ ψ ψ H KIN ϕ ϕ ϕ ψ ψ ψ , (13)
where

KIN ϕ ϕ ϕ = KI ϕ ϕ ϕ + N0 E ϕ ϕ ϕ 2 I.
Our optimization technique is an iterative algorithm where we maximize alternately the Rx waveform ψ ψ ψ, for a given Tx waveform ϕ ϕ ϕ, and the Tx waveform ϕ ϕ ϕ, for a given Rx waveform ψ ψ ψ.

Note that (13) can also be written as:

SIN R = ϕ ϕ ϕ H KS ψ ψ ψ ϕ ϕ ϕ ϕ ϕ ϕ H KIN ψ ψ ψ ϕ ϕ ϕ , (14)
where KS ψ ψ ψ and KIN ψ ψ ψ are expressed as:

KS ψ ψ ψ =   K-1 k=-(K-1) Π k Σ k ψ ψ ψψ ψ ψ H   Λ (15)
and

KIN ψ ψ ψ = KI ψ ψ ψ + N 0 E ψ ψ ψ 2 I (16)
with

KI ψ ψ ψ =   n Σ nN   K-1 k=-(K-1) Π k Σ k ψ ψ ψψ ψ ψ H     Ω-KS ψ ψ ψ .
(17) Thus, the optimization problem is equivalent to maximizing a generalized Rayleigh quotient. As N0 E ϕ ϕ ϕ 2 > 0, we can affirm that KIN ϕ ϕ ϕ is always invertible and relatively wellconditioned.

The main steps of the proposed algorithm, presented by Figure 1, are the following:

• Step 1: We initialize the algorithm with ϕ ϕ ϕ (0) , • Step 2: For the iteration (i), we compute ψ ψ ψ (i) as the eigenvector of (KIN ϕ ϕ ϕ (i) ) -1 KS ϕ ϕ ϕ (i) with maximum eigenvalue,

• Step 3: For the obtained ψ ψ ψ (i) , we determine ϕ ϕ ϕ (i+1) as the eigenvector of (KIN ψ ψ ψ (i) ) -1 KS ψ ψ ψ (i) with maximum eigenvalue,

• Step 4: We proceed to the next iteration, (i + 1),

• Step 5: We stop the iterations when we obtain a negligible variation of SINR. We note that eig, used in Figure 1, is a function that returns the eigenvector of a square matrix with the largest eigenvalue.

V. SIMULATION RESULTS

In this section, the performances of the proposed TR-POPS technique are evaluated. To show the gain in terms of SINR and Power Spectral Density (PSD), a comparison with POPS-OFDM and conventional OFDM with TR is also realized.

The results of POPS-OFDM algorithm applied to our system based on TR technique are carried out for a discrete timefrequency lattice. The optimal Tx/Rx waveform couple maximizing the SINR, ϕ ϕ ϕ opt , ψ ψ ψ opt , is evaluated for a Gaussian initialization waveform ϕ ϕ ϕ (0) .We presume having an exponential truncated decaying model. Figure 2 presents the evolution of the SINR versus the normalized Doppler spread B d /F for a normalized channel delay spread T m /T where Q = 128, N = 144, a lattice density is equal to 8/9 and for a waveform support duration D = 3T . The obtained results demonstrate that TR-POPS-OFDM approach improves the SINR with a gain of 2.3 dB for B d /F = 0.1 compared with POPS-OFDM and a gain that can reach 5.2 dB for B d /F = 0.02 compared with conventional OFDM with TR. Moreover, this figure is a mean to find the adequate couple (T, F ) of an envisaged application to insure the desired transmission quality. Figure 4 shows that, thanks to the TR technique, the obtained optimal transmit waveform, ϕ ϕ ϕ opt , reduces the OOB emissions by about 40 dB compared to the POPS-OFDM system without TR. We present in Figure 5 the Tx/Rx waveforms, ϕ ϕ ϕ opt and ψ ψ ψ opt , corresponding to the optimal SINR for Q = 128, N = 144, F T = 1 + 16 128 , B d T m = 0.001 and D = 3T . Since the channel is characterized by an Hermitian and symmetric response thanks to the TR effect, we obtain identical Tx/Rx waveforms as it is illustrated in this figure.

  0 Initialize φ         0 0 1 0 eig         φ φ ψ KIN KS         0 0 1 1 eig         ψ ψ K N KS φ I         1 1 1 1 eig         φ φ ψ KIN KS        

VI. CONCLUSION

In this paper, we studied the association of POPS-OFDM algorithm with TR precoding technique to design novel waveforms for 5G systems. To this end, we presented the corresponding system model and we derived the analytical SINR expression. Despite the additional complexity of applying the combination process, simulation results showed that the proposed approach offers a highly flexible behavior and better performances in terms of maximization of the SINR and reduction of the ISI/ICI. Another possible challenging research axis consists in applying this combination in MIMO-OFDM and FBMC/OQAM systems.

Figure 1 :

 1 Figure 1: Optimization philosophy.

Figure 2 :

 2 Figure 2: Optimized SINR as a function of B d F for Q = 128, SN R = 30 dB, B d T m = 10 -3 and D = 3T .

Figure 3

 3 Figure3illustrates the effect of TR by showing the evolution of the SINR with respect to the time-frequency parameter F T . As in Figure2, our proposed system outperforms the POPS-OFDM system and conventional OFDM with TR. The presented results reveal an increase in the obtained SINR that can reach 1.45 dB for F T = 1 + 8 128 compared with POPS-OFDM and an increase of 4.5 dB for F T = 1+48 128 compared to conventional OFDM with TR.Figure4shows that, thanks to the TR technique, the obtained optimal transmit waveform, ϕ ϕ ϕ opt , reduces the OOB

Figure 3 :

 3 Figure 3: SINR versus F T for Q = 128, SN R = 30 dB, B d T m = 10 -2 and D = 3T .

Figure 4 :

 4 Figure 4: PSD of the optimized transmit waveform for Q = 128, SN R = 30 dB, B d T m = 10 -3 , F T = 1.25 and D = 3T .

Figure 5 :

 5 Figure 5: Tx/Rx optimized waveforms for D = 3T .


M Verezhak

E F Rauch

M Véron

C Lancelon-Pin

J.-L Putaux

M Plazanet

A Gourrier

Ultrafine heat-induced structural perturbations of bone mineral at the individual nanocrystal level

Keywords: Bone, mineral nanocrystals, hydroxyapatite, TEM, electron diffraction, heating effects

The nanoscale characteristics of the mineral phase in bone tissue such as nanocrystal size, organization, structure and composition have been identified as potential markers of bone quality. However, such characterization remains challenging since it requires combining structural analysis and imaging modalities with nanoscale precision. In this paper, we report the first application of automated crystal orientation mapping using transmission electron microscopy (ACOM-TEM) to the structural analysis of bone mineral at the individual nanocrystal level. By controlling the nanocrystal growth of a cortical bovine bone model artificially heated up to 1000 ºC, we highlight the potential of this technique. We thus show that the combination of sample mapping by scanning and the crystallographic information derived from the collected electron diffraction patterns provides a more rigorous analysis of the mineral nanostructure than standard TEM. In particular, we demonstrate that nanocrystal orientation maps provide valuable information for dimensional analysis. Furthermore, we show that ACOM-TEM has sufficient sensitivity to distinguish between phases with close crystal structures and we address unresolved questions regarding the existence of a hexagonal to monoclinic phase transition induced by heating. This first study therefore opens new perspectives in bone characterization at the nanoscale, a daunting challenge in the biomedical and archaeological fields, which could also prove particularly useful to study the mineral characteristics of tissue grown at the interface with biomaterials implants.

Introduction.

Bone tissue is a biological nanocomposite material essentially composed of hydrated collagen fibrils of ~100 nm in diameter and up to several microns in length, reinforced by platelet-shaped nanocrystals of calcium phosphate apatite of ~ 4×25×50 nm 3 in size [START_REF] Weiner | The Material Bone: Structure-Mechanical Function Relations[END_REF]. These mineralized fibrils constitute the building blocks of bone tissue, and their specific arrangement is known to depend primarily on the dynamics of the formation and repair processes. Since these cellular processes can occur asynchronously in space and time, the mineralized fibrils adopt a complex hierarchical organization [START_REF] Weiner | The Material Bone: Structure-Mechanical Function Relations[END_REF], which was shown to be a major determinant of the macroscopic biomechanical properties [START_REF] Zimmermann | Intrinsic mechanical behavior of femoral cortical bone in young, osteoporotic and bisphosphonate-treated individuals in low-and high energy fracture conditions[END_REF]. Extensive research programs are therefore currently focused on bone ultrastructure for biomedical diagnoses or tissue engineering applications.

However, structural studies at the most fundamental scales remain challenging due to the technical difficulties imposed by nanoscale measurements and by the tissue heterogeneity. Nevertheless, as a natural extension of bone mineral density (BMD) analysis, an important marker in current clinical studies, the following key characteristics of the mineral nanocrystals have been identified as potential markers of age and diseases: chemical composition, crystallinity vs disorder, crystal structure, size, shape and orientation [START_REF] Matsushima | Age changes in the crystallinity of bone mineral and in the disorder of its crystal[END_REF][START_REF] Boskey | Variations in bone mineral properties with age and disease[END_REF]. Recent progress in the field showed that in order to obtain a deeper medical insight into the mechanisms of bone function, several such parameters need to be combined and correlated to properties at larger length scales [START_REF] Granke | Microfibril Orientation Dominates the Microelastic Properties of Human Bone Tissue at the Lamellar Length Scale[END_REF]. Interestingly, from a totally different point of view, the archaeological community has drawn very similar conclusions concerning nanoscale studies for the identification, conservation and restoration of bone remains and artifacts [START_REF] Chadefaux | Archaeological Bone from Macro-to Nanoscale: Heat-Induced Modifications at Low Temperatures[END_REF].

From a materials science perspective, this is a well identified challenge in the analysis of heterogeneous nanostructured materials. Yet, technically, a major difficulty stems from the fact that most of the identified nanostructural bone markers require individual measurements on dedicated instruments which are generally difficult to combine in an integrative approach.

One 'gold standard' in nanoscale bone characterization is X-ray diffraction (XRD), which allows determining atomic-scale parameters averaged over the total volume illuminated by the X-ray beam. An important result from XRD studies conducted with laboratory instruments is that the bone mineral phase has, on average, a poorly crystalline apatite structure which, to a certain extent, is induced by a high fraction of carbonate substitutions [START_REF] De | Le substance minerale dans le os[END_REF]. Such studies enabled localization of a substantial number of elements other than calcium and phosphorus present in bone via ionic substitutions [START_REF] Wopenka | A mineralogical perspective on the apatite in bone[END_REF], which can lead to serious pathological conditions, e.g. skeletal fluorosis [START_REF] Boivin | Fluoride content in human iliac bone: Results in controls, patients with fluorosis, and osteoporotic patients treated with fluoride[END_REF]. When an average description of bone properties is insufficient, synchrotron X-ray beams focused to a typical diameter of 0.1 -10 µm [START_REF] Schroer | Hard x-ray nanoprobe based on refractive x-ray lenses[END_REF][START_REF] Fratzl | Position-Resolved Small-Angle X-ray Scattering of Complex Biological Materials[END_REF][START_REF] Paris | From diffraction to imaging: New avenues in studying hierarchical biological tissues with x-ray microbeams (Review)[END_REF] operated in scanning mode allow mapping the microstructural heterogeneities. However, this remains intrinsically an average measurement and the current instrumentation limits prevent any analysis at the single mineral nanocrystal level.

Transmission electron microscopy (TEM) is a second 'gold standard' in bone characterization at the nanoscale. In high resolution mode, it allows reaching sub-angstrom resolution [START_REF] Xin | HRTEM Study of the Mineral Phases in Human Cortical Bone[END_REF] and therefore provides atomic details of the crystals. This increased resolution comes at the cost of the image field of view, which may not provide representative results due to the tissue heterogeneity. This limitation can partly be alleviated in scanning mode, which is more adapted to the collection of a large amount of data for statistical usage. In particular, for the process known as Automated Crystal Orientation Mapping (ACOM-TEM) [START_REF] Rauch | Automated crystal orientation and phase mapping in TEM[END_REF], diffraction patterns are systematically acquired while the electron beam is scanning micron-sized areas, such that the structural parameters of hundreds of individual nanocrystals may be characterized and used to reconstruct orientation maps with nanometer spatial resolution.

To our best knowledge, the present study is the first reported use of the ACOM-TEM method to analyze mineral nanocrystals in bone tissue. To demonstrate the potential of this technique for bone studies, a test object is required which structure should be as close as possible to native bone while offering a wide range of nanocrystal dimensions. Heated bone provides an ideal model for such purposes, ensuring a tight control over the nanocrystal size by adjusting the temperature.

This system was extensively studied in archeological and forensic contexts. Upon heating to 100-150 °C, bone is progressively dehydrated [START_REF] Legeros | Types of 'H2O' in human enamel and in precipitated apatites[END_REF] and collagen is considered to be fully degraded at ~ 400 °C [START_REF] Kubisz | Differential scanning calorimetry and temperature dependence of electric conductivity in studies on denaturation process of bone collagen[END_REF][START_REF] Etok | Structural and chemical changes of thermally treated bone apatite[END_REF]. Most X-ray studies concluded an absence of mineral crystal structure modifications before 400 °C, while a rapid crystal growth has been reported at ~ 750 °C [START_REF] Rogers | An X-ray diffraction study of the effects of heat treatment on bone mineral microstructure[END_REF][START_REF] Hiller | Bone mineral change during experimental heating: An X-ray scattering investigation[END_REF][START_REF] Piga | A new calibration of the XRD technique for the study of archaeological burned human remains[END_REF]. In a recent study we provided evidence that the mineral nanocrystals increase in size and become more disorganized at temperatures as low as 100 °C [START_REF] Gourrier | Nanoscale modifications in the early heating stages of bone are heterogeneous at the microstructural scale[END_REF]. In addition, many debates remain open concerning the nature of a postulated high temperature phase transition, the coexistence of different crystallographic phases, as well as the presence of ionic defects above and below the critical temperature of Tcr = 750 °C [START_REF] Greenwood | Initial observations of dynamically heated bone[END_REF]. The heated bovine cortical bone model therefore presents two main advantages to assess the potential of ACOM-TEM: 1) the possibility to fine-tune the mineral nanocrystal size upon heating and 2) the existence of a phase transition at high temperatures.

Using a set of bovine cortical bone samples in a control state and heated at eight temperatures ranging from 100 to 1000 °C, we show that ACOM-TEM provides enough sensitivity to probe fine crystalline modifications induced by heating; in particular, nanocrystal growth, subtle changes in stoichiometry and space group. Those results provide new insight into the detailed effects of heating on bone and validate the use of ACOM-TEM for fundamental studies of the nanoscale organization of bone tissue in different contexts.

Materials and methods.

Sample preparation: A bovine femur was obtained from the local slaughterhouse (ABAG, Fontanil-Cornillon, France). The medial cortical quadrant of a femoral section from the mid-diaphysis was extracted with a high precision circular diamond saw (Mecatome T210, PRESI) and fixed in ethanol 70 % for 10 days (supplementary information, Fig. S1). Nine 2×2×10 mm 3 blocks were cut in the longitudinal direction and subsequently dehydrated (48 hours in ethanol 70 % and 100 %) and slowly dried in a desiccator. One block was used as a control, while the others were heated to eight temperatures: 100, 200, 300, 400, 600, 700, 800 and 1000 °C for 10 min in vacuum (10 -2 mbar) inside quartz tube and cooled in air. The temperature precision of the thermocouple was ~ 2-3 °C and the heating rate was ~ 30-40 °C/min. The heating process resulted in color change, as shown in Fig. S2 of supplementary information. The samples were then embedded in poly-methyl methacrylate (PMMA) resin following the subsequent steps: impregnation, inclusion and solidification. For impregnation, a solution of methyl methacrylate (MMA) was purified by aluminum oxide and a solution of MMA was prepared with dibutyl phthalate in a 4:1 proportion (MMA1). The samples were kept at 4 °C in MMA1 for 5 days. For inclusion, the samples were stored in MMA1 solution with 1 w% of benzoyl peroxide for 3 days and in MMA1 solution with 2 w% of benzoyl peroxide for 3 days. The solidification took place in PTFE flat embedding molds covered by ACLAR film at 32 °C for 48 h. The resin-embedded blocks were then trimmed and cut with a diamond knife in a Leica UC6 ultramicrotome. The 50-nm-thick transverse sections (i.e., normal to the long axis of the femur) were deposited on 200 mesh Cu TEM grids coated with lacey carbon.

TEM data acquisition:

The measurements were performed using a JEOL 2100F FEG-TEM (Schottky ZrO/W field emission gun) operating at an accelerating voltage of 200 kV and providing an electron beam focused to 2 nm in diameter at sample position. A camera was positioned in front of the TEM front window to collect diffraction patterns as a function of scanning position with a frame rate of 100 Hz. The regions of interest were first selected in standard bright-field illumination (supplementary Fig. S3). The field of view for ACOM acquisition was chosen to be 400×400 nm 2 with a 10 ms acquisition time and a 2 nm step size. A sample-tocamera distance of 30 cm was chosen for all samples, except for the larger crystals treated at 800 and 1000 °C, for which a camera length of 40 cm was used, applying a precession angle of 1.2 ° at a frequency of 100 Hz in order to minimize dynamical effects [START_REF] Rauch | Automated crystal orientation and phase mapping in TEM[END_REF]. Following distortion and camera length corrections, a virtual brightfield image was reconstructed numerically by selecting only the transmitted beam intensities.

Radiation damage assessment:

No severe radiation damage was observed during ACOM-TEM data acquisition. This was assessed by independent bright-field acquisitions in the region close to the one scanned by ACOM-TEM for each heat-treatment temperature. These measurements were performed under the same conditions but with a smaller spot size of 0.7 nm (i.e. with a higher radiation dose) and were repeated 25 times to emphasize potential damage. Examples of bright-field images before the first and last (25th) frames are shown in Fig. S5 (supplementary information), showing very limited radiation damage.

ACOM-TEM analysis:

The data analysis relies on the comparison between the electron diffraction patterns collected at every scan position and simulated patterns (templates) calculated for a given crystal structure in all possible orientations [START_REF] Rauch | Automated crystal orientation and phase mapping in TEM[END_REF][START_REF] Rauch | Rapid spot diffraction patterns idendification through template matching[END_REF], thus allowing the reconstruction of crystal orientation maps (Fig. 1). The template matching was performed using the ASTAR software package from NanoMEGAS SPRL. In its native state, bone mineral is calcium phosphate close to a well-known hydroxyapatite, Ca 10 (PO 4 ) 6 (OH) 2 , a subset of the widespread geological apatite minerals [START_REF] Wopenka | A mineralogical perspective on the apatite in bone[END_REF]. Hence, our initial model for the crystal structure is a hexagonal space group (P63/m) with 44 atoms per unit cell and lattice parameters of a = 9.417 Å; c = 6.875 Å [START_REF] Hughes | Structural variations in natural F, OH, and Cl apatites[END_REF]. Every i-th acquired diffraction pattern collected at position (x i ,y i ) was compared to the full set of templates through image correlation (template matching) and the best fit gave the most probable corresponding crystallographic orientation. This first result can thus be represented in the form of a color map representing the crystalline orientation (Fig. 1f). To assess the quality of the fit, a second map of the correlation index Q i can be used, defining Q i as:

Q i = ∑ j=1 m P ( x j , y j ) T i ( x j , y j ) √ ∑ j=1 m P 2 ( x j , y j ) √ ∑ j=1 m T i 2 ( x j , y j )
where P(x j ,y j ) is the intensity of measured diffraction patterns and T i (x,y) corresponds to the intensity in every ith template. Q i compares the intensities of the reflection contained in the diffraction pattern, denoted by the function P(x,y), to the corresponding modeled intensities T i (x,y) in every i-th template in order to select the best match [START_REF] Rauch | Rapid spot diffraction patterns idendification through template matching[END_REF]. The degree of matching is represented in an 'index map' that plots the highest matching index at every location (Fig. 1g).

This parameter therefore weights the degree of correlation between the acquired and simulated diffraction patterns. If more than one phase is expected to be present, several sets of templates can simultaneously be fitted to the data in order to identify the best one. This allows constructing 'phase maps' in which each crystallographic phase is associated to a given color.

A critical aspect of the ACOM analysis is to judge the quality of the proposed crystal orientation. Indeed, it is worth emphasizing that the template matching algorithm always provides a solution, which requires evaluating the fidelity of the phase/orientation assignment, especially in the case of overlapping crystals. A reliability parameter R i was proposed to address this point. It is proportional to the ratio of the correlation indices for the two best solutions and is defined by:

R i =100(1- Q i2 Q i1 )
where Q i1 is the best solution (represented as a red circle on the stereographic projection in Fig. 1e) and Q i2 is the second best solution (shown with a green circle). Reliability values range between 0 (unsafe/black) and 100 (unique solution/white). In practice, a value above 15 is sufficient to ascertain the validity of the matching (Fig. 1h).

Results

.

Individual nanocrystal visualization.

Bone nanocrystal orientation maps were derived for the set of heated bone samples (Fig. 2) with the corresponding collective orientations on the 0001 stereographic projection. A moderate increase in crystal size was observed below 800 °C, followed by a rapid growth at higher temperatures with a crystal shape change from platelet to polyhedral.

An important advantage of ACOM is that the size and geometry measurements in bright-field TEM such as in Fig. 1a are generally performed on the whole image. In our case, the nanocrystals can exhibit a broad distribution of orientations, such that size estimation is impractical due to the platelet-shaped crystal geometry and leads to overestimated values as pointed out in earlier studies [START_REF] Ziv | Bone Crystal Sizes: A Comparison of Transmission Electron Microscopic and X-Ray Diffraction Line Width Broadening Techniques[END_REF]. The additional visualization of the nanocrystal orientation therefore allows restricting the dimensional analysis to crystals in the same orientation (the same color-code). For example, crystals displayed in red in Fig. 1f are oriented with their c-axis (longest axis) perpendicular to the scanning plane. The crystals displayed in green and blue have a 90° misorientation from this particular zone axis. While the colors are representative of the crystal orientations, the overall difference between the ratio of red vs green and blue colors at different temperatures reflect the spatial coherence in the phase index and, thus, the level of heterogeneity within a particular sample or between different samples. This fact is confirmed by additional scans with larger fields of view acquired for each sample of the temperature series (supplementary information, Fig. S4).

The crystallographic texture can be inferred from the stereographic projection along the 0001 direction showing that, below 800 °C, most regions consist of crystals mainly aligned along the c-axis, which is in good agreement with other larger scale XRD studies [START_REF] Voltolini | Hydroxylapatite lattice preferred orientation in bone: A study of macaque, human and bovine samples[END_REF]. At higher temperatures, we observed randomly oriented crystals resulting from a phase transition which nature will be discussed in a later section. Bone mineral nanocrystal size estimation.

The nanocrystal size was measured assuming two models: platelet (anisotropic) for the low temperature (LT) phase below 750 °C and polyhedral (isotropic) for the high temperature (HT) phase above 750 °C.

In the LT phase, the smallest platelet dimensions are obtained by line profiling (as displayed in Fig. 3a) of crystals having their c-axis aligned with the beam direction (represented by a red color). Crystal overlapping effects are not expected for this orientation as the long axis length of a platelet is comparable to the sample thickness. Therefore, the electron beam is predicted to pass through one crystal only, thus providing reliable size estimation.

For the HT phase, the nanocrystal size was determined using a spherical approximation (crystal diameters) based on the definition of grain boundaries, i.e. the locations where the misorientation between two pixels at the orientation map (Fig. 3c) is higher than a user selected threshold value (5° in the present case). An example of a grain boundary map for the 800 °C sample is shown in Fig. 3d. An average grain size is then estimated using a sphere diameter weighted by the grain's area in order to avoid misindexing from numerous small grains (mainly noise). For the HT phase, the nanocrystal size was determined using a spherical approximation (crystal diameters) based on the definition of grain boundaries, i.e. the locations where the misorientation between two pixels at the orientation map (Fig. 3c) is higher than a user selected threshold value (5° in the present case). An example of a grain boundary map for the 800 °C sample is shown in Fig. 3d. An average grain size is then estimated using a sphere diameter weighted by the grain's area in order to avoid misindexing from numerous small grains (mainly noise).

The grain size distribution was then obtained for the HT phase (Fig. 3e). According to the two described models for the nanocrystal sizes (summarized in Fig. 3b), we found that the smallest bone mineral particle size rises, on average, from 3.5 nm in the control state to 5.1 nm at 700 °C. Subsequently, the average particle diameter dramatically increases upon further heating: up to 70 nm at 800 °C and 94 nm at 1000 °C.

Identification of the high temperature apatite phase.

While stoichiometric hydroxyapatite has a calcium-to-phosphate ratio of 1.67, bone mineral is known to accommodate ~ 7 wt.% of carbonate and numerous other ionic substitutions [START_REF] Wopenka | A mineralogical perspective on the apatite in bone[END_REF]. To test the sensitivity of ACOM-TEM to changes in stoichiometry and therefore of space group, we used the data set containing the largest grains (1000 °C). We compared the fits obtained with the hydroxyapatite template against four template structures of different minerals occurring in nature with similar chemical composition and stoichiometry: alphatetra calcium phosphate (TCP) (α-Ca 3 (PO 4 ) 2 ) (space group -P21/a), beta-Ca pyrophosphate (β-Ca 2 P 2 O 7 ) (P4 1 ), tetra-calcium phosphate (Ca 4 (PO 4 ) 2 O) (P2 1 ), CaO (Fmm) and whitlockite (R3c), which chemical compositions and structures are shown in Fig. 4a-e. Those phases were previously encountered in synthetic apatites subjected to heat treatments and could, therefore, be potential candidates for heated bone [START_REF] Berzina-Cimdina | Research of Calcium Phosphates Using Fourier Transform Infrared Spectroscopy[END_REF]. Since these minerals have different space groups, all possible orientations are described by different fractions of the stereographic projection allowed by symmetry, as shown by the color map shapes in Fig. 4.

Two criteria can be used to conclude that the hydroxyapatite template (Fig. 4a) provides the best solution: 1) the highest index value that characterizes the quality of the solution is nearly twice larger for hydroxyapatite than for other apatite structures; in addition, 2) a given particle is expected to be fitted with the same orientation if monocrystalline, resulting in a uniform color, which is only fulfilled for hydroxyapatite. This analysis provides a first proof-of-concept that ACOM-TEM has a sufficient sensitivity to identify subtle variations of the crystal lattice that can be expected in highly disordered biological mineral structures, such as bone mineral nanocrystals. Other apatite minerals which are not expected to be found in bone tissue but have close to hydroxyapatite stoichiometry and chemical composition such as brushite (space group Ia), monetite ( P 1 ) and tuite ( R 3 m ) [START_REF] Schofield | The role of hydrogen bonding in the thermal expansion and dehydration of brushite, di-calcium phosphate dihydrate[END_REF][START_REF] Catti | Hydrogen bonding in the crystalline state. CaHPO4 (monetite), P1 or P1? A novel neutron diffraction study[END_REF][START_REF] Sugiyama | Structure and crystal chemistry of a dense polymorph of tricalcium phosphate Ca3 (PO4)2: A host to accommodate large lithophile elements in the earth's mantle[END_REF][START_REF] Calvo | The crystal structure of whitlockite from the Palermo quarry[END_REF] were also used to test the ACOM-TEM sensitivity. I.e. if ACOM-TEM had resulted in an equal probability to find these phases, it would clearly have suggested a lack of precision of the method. The analysis shows that this is not the case, i.e. these phases did not allow describing bone data as well as hydroxyapatite (see Figure S5a-d in supplementary information).

Space group: monoclinic or hexagonal ?

A common issue in the identification of a hydroxyapatite phase at different temperatures is the hypothesis of the existence of a hexagonal (P6 3 /m) to monoclinic (P2 1 /b) phase transition above T cr = 750 °C. The corresponding structures are shown in Fig. 5a. However, such a transition was mainly predicted by theoretical models [START_REF] Slepko | Hydroxyapatite: Vibrational spectra and monoclinic to hexagonal phase transition[END_REF][START_REF] Corno | Periodic ab initio study of structural and vibrational features of hexagonal hydroxyapatite Ca10(PO4)6(OH)2[END_REF] and was only observed in artificially synthesized hydroxyapatite [START_REF] Ma | Hydroxyapatite: Hexagonal or monoclinic?[END_REF][START_REF] Ikoma | Phase Transition of Monoclinic Hydroxyapatite[END_REF][START_REF] Suda | Monoclinic -Hexagonal Phase Transition in Hydroxyapatite Studied by X-ray Powder Diffraction and Differential Scanning Calorimeter Techniques[END_REF]. From the theoretical point of view, the monoclinic hydroxyapatite structure is thermodynamically more stable than the hexagonal one. Nevertheless, the hexagonal phase allows an easier exchange of OH-groups with other ions, which is necessary for bone tissue functions.

This issue was therefore addressed by matching the hydroxyapatite templates of the hexagonal [START_REF] Hughes | Structural variations in natural F, OH, and Cl apatites[END_REF] (Fig. 5b) and the monoclinic [START_REF] Elliott | Monoclinic hydroxyapatite[END_REF] (Fig. 5c) structures with the 1000 °C bone data set. The phase map in Fig. 5d represents the structure with the highest index at each scan point. Based on the index values, as well as on the uniform colorcode criterion for single crystals, one can conclude that the hexagonal space group in bone mineral HT phase is more probable than the monoclinic one. Discussion.

Our results provide a first demonstration that a structural analysis is possible at the single nanocrystal level within bone tissue using ACOM-TEM. This constitutes a valuable improvement combining the advantages of selected area electron diffraction (SAED) and TEM. While SAED produces a global diffraction pattern from all the nanocrystals probed by an electron beam defined by a micron sized aperture, ACOM-TEM allows a detailed structural analysis within a similar field of view.

The use of the artificially heated bone model shows that the phase and orientation can be unambiguously determined for temperatures above ~ 400 o C where the crystal size is > 5 nm (Fig. 2). Even below this temperature, where the situation is less clear, since the scanning resolution given by the beam size (~ 2 nm) is close to the actual nanocrystals size (~ 4-5 nm), the observation of larger fields of view reveals the presence of coherently oriented domains in the control sample with characteristic sizes of ~ 100 -200 nm which are compatible with the diameter of collagen fibrils (supplementary information, Fig. S4). Such crystallographic information is not available from standard bright-field TEM and can be used to obtain more rigorous estimates of the nanocrystals dimensions. An important limitation of a bright-field TEM assessment of the thickness of the mineral platelets is that the nanocrystals can be sectioned in different orientation, hence resulting in an artificial broadening in projection [START_REF] Ziv | Bone Crystal Sizes: A Comparison of Transmission Electron Microscopic and X-Ray Diffraction Line Width Broadening Techniques[END_REF]. These artifacts can be avoided using ACOM-TEM by only considering the platelets oriented on-edge, i.e. with the c-axis normal to the observation plane (Fig. 3). Similarly, depending on the thickness of the sample section, several nanocrystals may partly overlap, which further complicates the analysis. This case can be avoided by discarding the areas corresponding to a low correlation index and reliability parameters, i.e. to a poor structural refinement.

Our analysis reveals three stages of crystal growth upon heating: a first, moderate growth (from ~ 3.5 to 5.1 nm) between ambient temperature and ~ 700 o C, an order of magnitude increase in size (from ~ 5.1 to 70 nm) between 700-800 o C, followed by an additional growth (from 70 to 94 nm) between 800-1000 o C. This is in very good agreement with previous X-ray studies [START_REF] Rogers | An X-ray diffraction study of the effects of heat treatment on bone mineral microstructure[END_REF][START_REF] Piga | A new calibration of the XRD technique for the study of archaeological burned human remains[END_REF][START_REF] Piga | Is X-ray diffraction able to distinguish between animal and human bones?[END_REF]. However, X-ray scattering provides average information from all crystals illuminated by the beam, while ACOM-TEM allows a precise mechanistic interpretation of the heating process. In particular, the appearance of new uncorrelated orientations at temperatures > 700 o C strongly suggests a recrystallization process by fusion-recrystallization of smaller grains. Interestingly, this process can qualitatively be observed from 400 o C onwards (i.e. following total collagen degradation) in Fig. 2, as the structure qualitatively becomes more heterogeneous (larger polydispersity in crystal sizes) and disordered. However, there is a sharp transition from platelets to polyhedral crystals between 700-800 o C, clearly indicating a non-linear crystal growth.

One major difficulty in assessing the crystalline structure of bone is that the crystal chemistry is known to fluctuate, potentially giving rise to modulations of the intensity and breadth of the Bragg peaks which thus tend to overlap in XRD [START_REF] Posner | Crystal chemistry of bone mineral[END_REF][START_REF] Sakae | Historical Review of Biological Apatite Crystallography[END_REF]. For the same reason, the precise interpretation of Raman and FTIR spectra is still a matter of debate after decades of studies [START_REF] Wopenka | A mineralogical perspective on the apatite in bone[END_REF]. Since the crystal structure is used as input in the ACOM-TEM analysis, such fine deviations to an ideal crystal structure could not be assessed reliably. Nevertheless, the analysis conducted with different templates shows that this method permits a reliable distinction between different phases on the basis of intrinsic quality metrics (correlation index) and extrinsic ones, e.g. the spatial coherence (color uniformity) of the phase and orientation determination.

This allowed, in particular, testing the hypothesis proposed in previous XRD studies that diffraction patterns of heated bone could be equally well indexed with a monoclinic space group instead of a hexagonal one. The lattice parameters of the two structures were found to be very close, with a β angle close to 120 o for the monoclinic phase [START_REF] Piga | Is X-ray diffraction able to distinguish between animal and human bones?[END_REF]. The main difference is that the length of the b-axis can fluctuate significantly from the a-axis in the monoclinic case (contrary to the hexagonal case where the a-and b-axis are identical by definition). Thus, the difference between the two structures is relatively subtle but, in principle, a monoclinic structure should better account for a higher degree of crystallinity generated by heating, as demonstrated with synthetic apatites [START_REF] Ma | Hydroxyapatite: Hexagonal or monoclinic?[END_REF].

Our results conclusively show that even for samples treated at 1000 o C, the mineral is better represented by a hexagonal structure. This was further confirmed by a close manual examination of the proposed solutions for a number of representative diffraction patterns (example in supplementary information, Fig. S6). Because the monoclinic phase is more representative of stoichiometric hydroxyapatite, the fact that bone mineral is better indexed by a hexagonal group implies that there is still a substantial degree of disorder in the crystal structure even for bone heated at high temperatures.

It is important to note that the ACOM setup can be readily implemented in standard existing TEM instruments and can therefore provide a close-to routine basis for biological, medical and archeological studies. Given the resolution level, we believe that ACOM-TEM could be advantageously exploited to analyze the interface layer (typically < 1 µm) between biomaterials and bone formed at the surface of implants, a critical aspect of osseointegration [START_REF] Davies | Bone bonding at natural and biomaterial surfaces[END_REF][START_REF] Legeros | Calcium phosphate-based osteoinductive materials[END_REF]. TEM was widely used to investigate the tissue structure at this interface [START_REF] Grandfield | High-resolution three-dimensional probes of biomaterials and their interfaces[END_REF][START_REF] Palmquist | Bone--titanium oxide interface in humans revealed by transmission electron microscopy and electron tomography[END_REF], but the collective mineral nanocrystals structure and organization was never analyzed. Similarly, in the biomedical field, severe pathological perturbations of mineral nanocrystals have been reported in many bone diseases, e.g. osteoporosis [START_REF] Rubin | TEM analysis of the nanostructure of normal and osteoporotic human trabecular bone[END_REF], osteogenesis imperfecta [START_REF] Fratzl-Zelman | Unique micro-and nano-scale mineralization pattern of human osteogenesis imperfecta type VI bone[END_REF] and rickets [START_REF] Karunaratne | Significant deterioration in nanomechanical quality occurs through incomplete extrafibrillar mineralization in rachitic bone: Evidence from in-situ synchrotron X-ray scattering and backscattered electron imaging[END_REF], for which a detailed nanoscale description is still lacking [START_REF] Gourrier | Scanning small-angle X-ray scattering analysis of the size and organization of the mineral nanoparticles in fluorotic bone using a stack of cards model[END_REF]. Additionally, this method could also have a positive impact in the archaeological field, since diagenetic effects associated with long burial time of bone remains are known to affect the mineral ultrastructure in numerous ways, hence impacting the identification and conservation of bone artifacts [START_REF] Reiche | The crystallinity of ancient bone and dentine: new insights by transmission electron microscopy[END_REF]. Finally, it should be mentioned that ACOM-TEM would most likely benefit from more advanced sample preparation methods such as focused ion beam milling coupled to scanning electron microscopy (FIB-SEM) which has been shown to better preserve the tissue ultrastructure [START_REF] Jantou | Focused ion beam milling and ultramicrotomy of mineralised ivory dentine for analytical transmission electron microscopy[END_REF][START_REF] Mcnally | A Model for the Ultrastructure of Bone Based on Electron Microscopy of Ion-Milled Sections[END_REF][START_REF] Reznikov | Three-dimensional structure of human lamellar bone: the presence of two different materials and new insights into the hierarchical organization[END_REF].

Conclusion.

In the present work, we showed that both insights, the direct visualization of individual bone nanocrystals and structural information, can simultaneously be accessed using ACOM-TEM analysis. The mineral nanocrystal orientation, crystallographic phase and symmetry can be quantified, even in biological samples such as bone tissue that are known to be very heterogeneous down to the nanoscale. Our analysis of a heated bone model points to a crystal growth by fusion and recrystallization mechanisms, starting from ̴ 400 o C onwards with a sharp transition between 700 o C and 800 o C. By testing different phases corresponding to deviations from the hydroxyapatite stoichiometry as input for the structural refinement, we were able to assess the sensitivity of ACOM-TEM. We tested the hypothesis of a monoclinic space group attribution to the bone sample heated at 1000 o C and found that a hexagonal structure was more probable, suggesting the presence of crystalline defects even after heating at high temperatures. We therefore believe that ACOM-TEM could have a positive impact on applied research in biomaterials development, biomedical investigations of bone diseases and, possibly, analysis of archaeological bone remains. Bright field images with corresponding ACOM orientation maps.

Fig. 1 :

 1 Fig. 1: Generalized scheme of ACOM-TEM acquisition and data interpretation. (a) Bright-field (BF) image of bone tissue with the illustration of the scan area (dots diameter is 2 nm, enhanced for visibility), (b) recorded set of diffraction patterns, (c) example of fit using (d) the structure template for hydroxyapatite. (e) Fraction of stereographic projection with inverse pole figure color map, (f) orientation, (g) index and (h) reliability maps (high values appear brighter). Scale bar: 100 nm.

Fig. 2 :

 2 Fig. 2: Mineral crystal orientation distribution in bone tissue as a function of temperature. Orientation maps reconstructed from ACOM-TEM from control to 1000 °C with corresponding inverse pole figure color map (scale bar: 100 nm). Inset: collective crystal orientations on 0001 stereographic projection with the color bar normalized to the total number of crystals per scan.

Fig. 3 :

 3 Fig. 3: Bone mineral crystal size evolution with temperature. (a) an example of size measurement by line profiling along the c-axis of the platelet-shaped crystals for the LT regime (700 °C); (b) average crystal size vs temperature (black -smallest crystal size, red -crystal diameter). Examples of polyhedral-shaped crystal diameter measurement for the HT regime (800 °C): (c, d) indicate orientation and grain boundary maps, respectively, and (e) represents the distribution of crystal diameters.

Fig. 4 :

 4 Fig.4: Phase sensitivity. a-f, orientation maps for 1000 °C heated bone data fitted with six apatite structures with corresponding inverse pole figure color maps. Maximum index values (Imax) are given for comparison. The hydroxyapatite structure produces the best fit indicated by the highest Imax value and homogeneous colors within single crystals. Scale bar: 50 nm.

Fig. 5 :

 5 Fig. 5: Hexagonal vs monoclinic symmetry. (a) crystal structures of hexagonal and monoclinic hydroxyapatite (view along c-axis); (b, c) corresponding orientation maps for the 1000 °C sample with color code and maximum index values; (d) phase map showing mainly the presence of hexagonal phase. Scale bar: 50 nm.

Fig. S3 :

 S3 Fig. S3: Regions of scans for the set of heat-treated bone samples. Bright field images with corresponding ACOM orientation maps. Scale bar in orientation maps: 100 nm.

Fig. S4 :

 S4 Fig. S4: Larger fields of view scan regions for the set of heat-treated bone samples.Bright field images with corresponding ACOM orientation maps.

Acknowledgements.

The authors would like to acknowledge M. Morais from SIMaP for the support with the heating apparatus, D. Waroquy (ABAG, Grenoble, France) for providing the bovine samples, and the NanoBio-ICMG Platform (FR 2607, Grenoble) for granting access to the TEM sample preparation facility.

Author contributions.

M.V. ‡ , C.L.P. and J.L.P. prepared samples; M.V. ‡ , E.F.R., M.V., M.P. and A.G. performed the research; A.G., M.P. and E.F.R. provided the financial support for the project; M.V. ‡ and E.F.R. analyzed data; M.V. ‡ and A.G. wrote the paper with contributions from all authors; A.G. and M.P. designed the research.

Supplementary information


Mohamed El Demirci

Hassan Mansori

S Mezghani
email: sabeur.mezghani@chalons.ensam.fr

I Demirci

M El Mansori

H Zahouani

Energy efficiency optimization of engine by frictional reduction of functional surfaces of cylinder ring-pack system

Keywords: Honing process, Surface roughness, elastohydrodynamic friction, Cylinder engine

Energy efficiency optimization of engine by frictional reduction of functional surfaces of cylinder ring-pack system

Sabeur Mezghani, Ibrahim

Penalty term parameter z r Pressure viscosity index (Roelands), z r = p r /(ln(η 0 ) + 9.67)

Introduction

The surface features of a cylinder liner engine are the fingerprint of the successive processes the surface has undergone, and they influence the functional performance of the combustible engine [START_REF] Caciu | Parametric textured surfaces for friction reduction in combustion engine[END_REF][START_REF] Tomanik | Friction and wear bench tests of different engine liner surface finishes[END_REF][START_REF] Pawlus | A study on the functional properties of honed cylinders surface during runningin[END_REF][START_REF] Mcgeehan | A literature review of the effects of piston and ring friction and lubricating oil viscosity and fuel economy[END_REF][START_REF] Srivastava | Effect of liner surface properties on wear and friction in a non-firing engine simulator[END_REF]. Therefore, surfaces and their measurement provide a link between the manufacture of cylinder bores and their functional performances [START_REF] Whitehouse | Surfacesa link between manufacture and function[END_REF]. Hence, the quantitative characterization of surface texture can be applied to production process control and design for functionality [START_REF] De Chiffre | Quantitative Characterisation of Surface Texture[END_REF]. The optimum surface texture of an engine cylinder liner should ensure quick running-in, minimum friction during sliding, low oil consumption, and good motor engine operating parameters in terms of effective power and unitary fuel consumption. Increasingly stringent engine emissions standards and power requirements are driving an evolution in cylinder liner surface finish [START_REF] Lemke | Characteristic parameters of the abbot curve[END_REF]. Unfortunately, the full effect of different cylinder liner finishes on ring-pack performance is not well understood [START_REF] De Chiffre | Quantitative Characterisation of Surface Texture[END_REF].

In mass production of internal combustion engine cylinder liners, the final surface finish on a cylinder bore is created by an interrupted multistage abrasive finishing process, known as the plateau-honing process. In honing, abrasive stones are loaded against the bore and simultaneously rotated and oscillated. Characteristically, the resulting texture consists of a flat smooth surface with two or more bands of parallel deep valleys with stochastic angular position. Figure 1 shows a typical plateau-honed surface texture from an engine cylinder.

To guarantee efficient production at industrial level of a cylinder liner of specific shape with acceptable dimensional accuracy and surface quality, three honing stages are usually required: rough honing, finish honing, and plateau honing. The surface texture is presumably provided by the ″finish honing″ [START_REF] Sabri | Multiscale study of finish-honing process in mass production of cylinder liner[END_REF][START_REF] Sabri | Functional optimisation of production by honing engine cylinder liner[END_REF]. Thus careful control of this operation is central to the production of the structured surface so that the cylinder liner will fulfil its mechanical contact functionalities in piston ring/cylinder liner assemblies (i.e. running-in performance, wear resistance, load-carrying capacity, oil consumption, etc.) [START_REF] Sabri | Multiscale study of finish-honing process in mass production of cylinder liner[END_REF][START_REF] Pawlus | Effects of honed cylinder surface topography on the wear of piston-piston ring-cylinder assemblies under artistically increased dustiness conditions[END_REF]. ring-pack friction reduction through cylinder liner finish optimization it is necessary to be able to distinguish the effect of each process variable on the roughness of these honed surfaces [START_REF] De Chiffre | Quantitative Characterisation of Surface Texture[END_REF].

In this work, strategies for piston ring-pack friction reduction through cylinder liner finish optimization were analyzed with the goal of improving the efficiency of selection of the honing process variable. The fundamental aim was to find a relation between the honing operating variables and the hydrodynamic friction at the piston rings/cylinder interface. An additional aim was to determine how the cylinder surface micro-geometry of plateau-honed cylinders affects the predicted friction. Thus, an experimental test rig consisting of an industrial honing machine instrumented with sensors to measure spindle power, expansion pressure, and honing head displacement was developed. Honing experiments were carried out using honing stones with varying sizes of abrasive grits and varying expansion speeds, that is, the indentation pulse of the honing stone's surface against the liner wall. Furthermore, a numerical model of lubricated elastohydrodynamic contact was developed to predict the friction performances and lubricant flow of the various liner surface finishes. It uses the real topography of the liner surface as input. In fact previous studies have found that the detailed nature of the surface finish plays an important role in ring friction and oil film thickness predictions [START_REF] Jocsak | The effects of cylinder liner finish on piston ring-pack friction[END_REF]. An appreciation of the limitations of the surface roughness parameters commonly used in automotive industries in providing a link between the honing process and the generated surface performance in the hydrodynamic regime is presented.

Experimental procedure

In this work, honing experiments were carried out on a vertical honing machine with an expansible tool (NAGEL no. 28-8470) (Figure 2). The workpiece consists of four cylinder liners of a lamellar gray cast iron engine crankcase.

The steps involved in the fabrication of the cylinder liners before the finish-honing operation are boring and rough honing, respectively (Table 1 treatment by impregnation with sulfur. Another interesting variation in the feed system is the expansion mechanism in the honing head, where three expansion velocities "V e " (1.5µm/s, 4µm/s and 8µm/s) were considered. All the other working variables were kept constant (Table 1). Note that the rough and finish honing operations use a mechanical expansion system and the plateau honing uses a hydraulic system. For each combination of grit size and expansion velocity, tests were repeated five times. Thus, the sensitivity of the produced surface finish to its generation process was considered.

Negative surface replicas made of a silicon rubber material (Struers, Repliset F5) were used to assess the texture of honed surfaces after the plateau-honing stage at the mid-height of the cylinder bore specimen. Topographical features of replica surfaces were measured in three locations by a three-dimensional white light interferometer, WYKO 3300 NT (WLI). The surface was sampled at 640 × 480 points with the same step scale of 1.94 μm in the x and y directions. Form component is removed from acquired 3D data using least square method based on cubic Spline function.

We can assume that the initial roughness of the cylinder bore has no influence on the obtained surface texture in this study. It affects only the honing cycle and the stone life, that is, the wear of the abrasive grits. In fact, the thickness of the removed material after finish honing (32.17 ± 2.21 µm) is greater than the total height of the original surface, which was about 24.56 ± 6.75 µm. This means that the finish honing operation completely penetrates the original surface topography and generates a new surface texture.

Numerical model for hydrodynamic friction simulation in piston ringpack system

A numerical model was developed to estimate friction at the ring-liner-piston contact. It takes into account the real topography of the cylinder liner. The scope of this model is to predict qualitatively the friction coefficient obtained to optimize the performances when the groove characteristics of cylinder liner surfaces are variables.

Geometry definition

An incompressible viscous fluid occupying, at a given moment, a field limited by a smooth plane surface P and by a rough surface R is considered. This field is represented on figure 3 (we did not represent the profile in the x2 direction). It extends from 0 to l1, 0 to l2 and h (x1, Figure 3 The separation field between a smooth surface and a rough one

EHL Equations

To estimate the pressure distribution, film thickness, and friction coefficient, a full system approach for elastohydrodynamic lubrication (EHL) was developed. The Reynolds equations have been written in dimensionless form using the Hertzian dry contact parameters and the lubricant properties at ambient temperature. To account for the effects of non-Newtonian lubricant behaviour effectives viscosities are introduced.

With the boundary condition P = 0 and the cavitation condition (or free boundary condition) P ≥ 0 are used everywhere Special treatment is used for cavitation condition as explained below. In this equation, is equal to . and are the effective viscosities in the X and Y directions, respectively. For point contact, it is not possible to derive these effective viscosities analytically. The perturbational approach described by Ehret et al. [START_REF] Ehret | On the lubricant transport condition in elastohydrodynamic conjunctions[END_REF] is used.

This analysis is based on the assumption that the shear stresses are only partially coupled and that the mean shear stress is negligible in the y direction [START_REF] Ehret | On the lubricant transport condition in elastohydrodynamic conjunctions[END_REF][START_REF] Greenwood | Two-dimensional flow of a non-Newtonian lubricant[END_REF].

In our model, the Eyring model is used. The perturbational approach leads to the following dimensionless effective viscosities:

where the dimensionless mean shear stress is written as:

where S is the slide to roll ratio (2(u 1u 2 )/(u 1 + u 2 )) and N is given by :

The constant parameters of Equation ( 5) are given in the nomenclature.

The lubricant's viscosity and density are considered to depend on the pressure according to the Dowson and Higginson relation [START_REF] Dowson | Elastohydrodynamic lubrication. The fundamentals of roller and gear lubrication[END_REF] (Eq. 6) and Roelands equation [START_REF] Roelands | Correlational aspects of the viscosity-temperature-pressure relationships of lubricant oil[END_REF] (Eq. 7):

where 0 is the density at ambient pressure.

where η 0 is the viscosity at ambient pressure, p r is a constant equal to 1.96 x 10 8 , and z r is the pressure viscosity index (z r = 0.65).

The film thickness equation is given in dimensionless form by the following equation:

is the height of the liner surface topography at each position (X,Y). H 0 is a constant determined by the force balance condition:

The normal elastic displacement of contacting bodies is obtained by solving the linear elasticity equations in three-dimensional geometry with appropriate boundary conditions [START_REF] Habchi | A full-system approach of the elastohydrodynamic line/point contact problem[END_REF][START_REF] Habchi | A full-system approach to elastohydrodynamic lubrication problems: application to ultra-low-viscosity fluids[END_REF]. The geometry (Ω) used (figure 4) is large enough compared to contact size (Ω c ) to be considered as semi-infinite structures. The linear elastic equations consist of finding the displacement vector U in the computational domain with the following boundary conditions:

In order to simplify the model, the equivalent problem defined by [START_REF] Habchi | A full-system approach to elastohydrodynamic lubrication problems: application to ultra-low-viscosity fluids[END_REF] is used to replace the elastic deformation computation for both contacting bodies. One of the bodies is assumed to be rigid while the other accommodates the total elastic deformation. The following material properties of the bodies are used in order to have (w is the dimensionless absolute value of the Z-component of the displacement vector):

where E i and i are the Young's modulus and Poisson's coefficient, respectively, of the material for contacting bodies (i = 1, 2).

Finally, the friction coefficient is evaluated by the following formula: these negative pressures are not relevant. In such cases, the fluid will evaporate and the pressure is limited by the vapor pressure of the fluid. This process is the cavitation. This problem is usually solved by setting the negative pressure to zero. This ensures that there will be zero pressure and gradient pressure on the free boundary. In the full system approach, this solution is not possible and the penalty method is used as an alternative, as explained in [START_REF] Habchi | A full-system approach to elastohydrodynamic lubrication problems: application to ultra-low-viscosity fluids[END_REF].

Cavitation problem

This method was introduced in EHL by Wu [START_REF] Wu | A penalty formulation and numerical approximation of the Reynolds-Hertz problem of elastohydrodynamic lubrication[END_REF]. An additional penalty term was introduced in the Reynolds equation:

where is a large positive number and is the negative part of the pressure distribution. This penalty term constrains the system to P 0 and forces the negative pressure to zero.

Numerical procedure

The iterative process is repeated until the maximum relative difference between two consecutive iterations reaches 10 -6 . Table 2 summarizes the fluid properties and contact parameters used in our simulation. The difference between our model and Venner and Lubrecht work [START_REF] Venner | Multilevel Methods in Lubrication[END_REF] is less than 1%.

This test confirms the validity of the model presented in this paper. Figure 5 show an example of a pressure distribution and film thickness profiles along the central line in the X direction for rough surface like one presented in Figure 6.

Table 3 Comparison of the current model with the Venner & Lubrecht model [START_REF] Venner | Multilevel Methods in Lubrication[END_REF] Venner and Lubrecht [START_REF] Venner | Multilevel Methods in Lubrication[END_REF]

Results and discussion

The numerical model presented in Section 3 was used to predict friction in the ring-linerpiston contact and to analyze possible friction reduction strategies in the piston ring-pack.

Table 4 regroups all the experimental and numerical results. As a result of simulation of the cylinder ring-pack contact, only the average friction coefficients were compared.

Cylinder liner surface roughness distribution can differ between surfaces with the same root mean square roughness. This difference can have a significant effect on the performance and behavior of the surface within the piston ring-pack system. The surface finish created by the honing process is controlled by the size and dispersion of abrasive particles adhering to the surface of the honing sticks. Figure 6 shows the effect of varying the grit size of abrasive honing stones used in the finish honing stage on the threedimensional topographical features of the produced surface. As shown in Figure 6, coarse abrasive grits yield deeper and larger lubrication valleys and consequently rougher surfaces. working variables: V e in finish honing stage = 4 µm/s; all others parameters are kept constant and are given in Table 1.)

Influence of abrasive grit size and expansion velocity on the impregnated surface texture and its friction performance within the piston ring-pack system

As a result of these honing experiments carried out with various sizes of abrasive grits, Figure 7 presents predicted values of the coefficient of friction and mean oil film thickness as a function of abrasive grit size of honing stone used at the finish honing stage. It demonstrates clearly that the surface texture achieved with finer abrasive grits yields a lower hydrodynamic friction coefficient in a cylinder ring-pack system than that obtained by coarse abrasive grits. Since the generated honed surfaces have the same honing cross-hatch angle, the differences in predicted hydrodynamic friction observed between these different finishes are mainly a result of surface peak and valley characteristics. Hence, the increase in valleys volume may increase the oil through the valleys of the surface, yielding a decrease in the oil film thickness, which will in turn induce an increase in hydrodynamic friction.

A reduction of the expansion speed operating condition leads to lower valley depth on the surface texture impregnated during coarse honing as observed in Figure 8. This figure also

shows that the expansion velocity has no effect on the spatial morphology of the generated surface texture and hence the roughness scale, as demonstrated by multiscale surface analysis in [START_REF] Lemke | Characteristic parameters of the abbot curve[END_REF]. 1).

Relationship between liner surface friction performance and honing process efficiency

The specific energy is used as a fundamental parameter for characterizing the honing process.

It is defined as the energy expended per unit volume of material removed. The specific honing energy defines the mechanisms of removal of material from the operated workpiece. It is calculated from the following relationships:

* honing Pm t Esp Qw ( 12
)
where honing t is the effective honing time, Pm is the average power absorbed by the honing process, calculated as the difference between on-load power recorded during the finishing and average off-load power recorded before and after the test, and Qw is the volumetric removal given by the following equation:

22 c Qw H D d ( 13
)
where c H is the cylinder height, and D and d are the cylinder diameter before and after the finish honing operation, respectively. It suggests that the optimum coefficient of friction with a good honing efficiency is reached by using a grit size of 80-100 µm and an expansion velocity equal to 4 µm/s. Furthermore, significant smooth surfaces are produced by plateau-honing with fine abrasive grit sizes due to the low indentation capacity of the fine abrasive grain. This generated surface texture also yields a low predicted coefficient of friction in the piston-ring-liner interface.

However, the use of fine abrasives has the lowest efficiency. In fact, it presents lower material removal and consumes a large specific energy due to the predominance of the plowing abrasion mechanism [START_REF] Sabri | Multiscale study of finish-honing process in mass production of cylinder liner[END_REF][START_REF] Sabri | Functional optimisation of production by honing engine cylinder liner[END_REF]. This yields a lower stone life and generates undue tool wear.

Thus, to ameliorate the honing efficiency, conventional abrasives can be replaced by superabrasive crystals which do not wear or break rapidly.

Roughness characteristics of optimal plateau-honed surface texture

To give a rough estimate of the potential side effects of the surface optimization, surface roughness has been evaluated using the functional roughness parameters R k (height of the roughness core profile), R pk (reduced peak height), and R vk (reduced valley height) given by the ISO 13565-2 standard [START_REF] Sabri | Multiscale study of finish-honing process in mass production of cylinder liner[END_REF][START_REF] Jocsak | The effects of cylinder liner finish on piston ring-pack friction[END_REF]. These parameters are obtained based on the analysis of a bearing curve (the Abott-Firestone curve), which is simply a plot of the cumulative probability distribution of surface roughness height [START_REF] Sabri | Multiscale study of finish-honing process in mass production of cylinder liner[END_REF]. The peak height is an estimate of the These bubble plots show that all optimal surfaces in the hydrodynamic lubrication regime belong to the domain defined by: R pk < 1 µm R k < 3 µm R vk < 2.5 µm However, this critical domain does not guarantee the optimal behavior of the honed surface finish. For example, in Figure 12, the criterion R vk < 2.5 µm cannot exclude honed surfaces that induce a high friction coefficient. For the R pk and R k criteria, the same observation can be expressed. This suggests that the standard functional roughness parameters commonly used in the automotive industry cannot give a good classification of plateau-honed surfaces according to their functional performance. Table 5 shows linear correlation coefficients between roughness parameters and predicted coefficients of friction.

Table 5 The linear correlation coefficient between roughness parameters and the coefficient of friction

Correlation coefficient between Rpk and µ 0.666

Correlation coefficient between Rk and µ 0.664

Correlation coefficient between Rvk and µ 0.658

Hence, these standard functional parameters are not sufficient to give a precise and complete functional description of "ideal" honed surfaces. This can be attributed to the fact that bearing curve analysis is one-dimensional and provides no information about the spatial characteristics and scale of surface roughness.

Conclusion

This work focused on developing ring-pack friction reduction strategies within the limitations of current production honing processes. First, three-dimensional honed surface topographies were generated under different operating conditions using an instrumented industrial honing machine. Then, the three-dimensional surface topography of each honed cylinder bore is input into a numerical model which allows the friction performance of a cylinder ring-pack system in an EHL regime to be predicted. The strategy developed allows manufacturing to be related to the functional performance of cylinder bores through characterization. The results show that an increase in grit size will lead to an increase in surface roughness, with deeper valleys leading to an increase in hydrodynamic friction. They also show that the standard functional surface roughness parameters which are commonly used in the automotive industries do not provide a link between the honing process and the generated surface performance in the hydrodynamic regime.

Note that the analysis presented in this study does not take into consideration the effects of cylinder surface topography on its ability to maintain oil, that is, the oil consumption level. Experimental studies using a reciprocating bench tester will be carried out to evaluate the effect of honing operating conditions and cylinder surface topography on scuffing and oil consumption.

  coefficient (Gpa -1 ) Elastic deflection of the contacting bodies (m) Dimensionless elastic deflection of the contacting bodies D, d Cylinder diameter before and after the finish honing operation, respectively (mm) E eq , eq Equivalent Young's modulus (Pa) and Poisson's coefficient, respectively E i , i Young's modulus (Pa) and Poisson's coefficient of component i, respectively E r Reduced modulus of elasticity Esp Specific honing energy (J/mm 3 ) temperature zero-pressure viscosity (Pa.s) Effective viscosities in the X and Y of curvature in x direction (m) Rg height of the liner surface topography at each position (m) Lubricant density (kg.m -3 ) Dimensionless lubricant density (= / 0 ) 0 Lubricant's density under ambient condition (kg.m -3 ) S Slide to roll ratio: S = 2(u 1u 2 )/(u 1 + u 2 ) Surface velocity of body I in x-direction (m.s -1 ) u m Mean entrainment velocity (m.s -1 ) x,y,z Space coordinates (m) X,Y Dimensionless space coordinates (=x/a, y/a) w Value of Z-component of the displacement vector (m)

Figure 1

 1 Figure 1 Typical plateau-honed surface texture

Figure 2

 2 Figure 2 (a) Vertical honing machine with expansible tool; (b) Schematic representation of the honing head in continuous balanced movement.

  x2) respectively according to x1, x2 and x3. h (x1, x2) represents the fluid thickness. The smooth body is animated by a movement at the constant velocity "u1" along Ox1 axis whereas rough surface is static.

Figure 4

 4 Figure 4 Scheme of geometric model for computations of the elastic deformation (1) and of the Reynolds equation (2).

  pressure appears in the resolution of the Reynolds equation. Physically,

  Reynolds equation, linear elastic equations, and load balance equation are simultaneously solved using a Newton-Raphson resolution procedure. The dimensionless viscosity , density , and film thickness H in the Reynolds equation are replaced by the expression given above. Except for the load balance equation, a standard Galerkin formulation is used. For the load balance equation, an ordinary integral equation is added directly with the introduction of an unknown H 0 . Unstructured variable tetrahedral meshing is used for both Reynolds and linear elastic equations. A total number of 100000 degrees of freedom are used in the simulation. An

Figure 5

 5 Figure 5 Pressure (P) and film thickness (H) profiles along the central line in the X direction.

Figure 6

 6 Figure 6 Three-dimensional topographies of plateau-honed surfaces produced using different abrasive grit sizes in the finish honing stage: (a) 40 µm, (b) 110 µm, and (c) 180 µm. (Process a b c

Figure 7

 7 Figure 7 Evolution of friction coefficient of the cylinder ring-pack as a function of the abrasive grit size in the honing finish stage.

Figure 8

 8 Figure 8 Three-dimensional topographies of plateau-honed surfaces produced using three different expansion velocities in the finish honing stage: (a) 1.5 µm/s, (b) 4 µm/s, and (c) 8 µm/s. (Process working variables: Abrasive grit size in the finish honing stage equal to 110µm; all other parameters are kept constant and are given in Table1).

Figure 9 ,

 9 Figure 9, which presents a plot of the friction coefficient versus specific energy, highlights the link between honing process operating conditions and the functional behavior of plateau-honed surfaces in the hydrodynamic lubrication regime.

Figure 9

 9 Figure 9Predicted coefficient of friction of surface texture generated with different honing grit sizes and indentation pulse as a function of consumed specific energy during finish honing. (The size of circles is proportional to the size of the honing abrasives which varies from 30 µm to 180 µm.)

  Figures 10,[START_REF] Pawlus | Effects of honed cylinder surface topography on the wear of piston-piston ring-cylinder assemblies under artistically increased dustiness conditions[END_REF], and 12 display, for different abrasive grit sizes and at various expansion velocities, the existing correlation between the predicted friction coefficient within the cylinder ring-pack system and the functional roughness parameters of the plateau-honed surfaces of the cylinder bore, R pk , R k , and R vk , respectively.

Figure 10 Figure 11 Figure 12

 101112 Figure 10 Predicted coefficient of friction vs. functional roughness parameter R pk . (The size of circles is proportional to the size of the honing abrasives which varies from 30 µm to 180 µm.)

Table 1

 1 Honing working conditions

	Honing process variables	Rough honing Finish honing Plateau honing
	V a : Axial speed (m/min)	28	28	28
	V r : Rotation speed (rpm)	230	230	230
	Honing time (sec)	20	15	2
	Expansion type	Mechanical	Mechanical	Hydraulic
	V e : Expansion velocity (µm/s) 5	1.5, 4, and 8
	Number of stones	6	6	6
	Abrasive grit type	Diamond	Silicon carbide Silicon carbide
	Grain size (µm)	125	30-180	30
	Bond type	Metal	Vitrified	Vitrified
	Abrasive stone dimensions	2 × 5 × 70	6 × 6 × 70	6 × 6 × 70
	(mm × mm × mm)

Table 2

 2 Parameters for our simulation used with rough surfaces

	Parameter Value Parameter Value
	F N (N)	500	(GPa -1 ) 22.00
	u m (m.s -1 ) 10.0 R x (m)	0.04
	η 0 (Pa.s)	0.04 E i (GPa)	210
	i	0.3	τ 0 (MPa)	0.5

Table 3

 3

gives dimensionless central and minimum oil film thickness for the following dimensionless Moes and Venner parameters M=200 and L=10.

Table 4

 4 Honing process variables, roughness parameters of honed surfaces and theirs predicted coefficient of friction

	Ve	Grit size	Rpk	Rk	Rvk	µ (%)
	1.5	180	0.663	1.859	1.960	2.495
	1.5	145	0.777	1.798	2.429	2.466
	1.5	110	0.679	1.813	2.025	2.474
	1.5	90	0.564	1.954	1.825	2.459
	1.5	80	0.566	1.628	1.678	2.426
	1.5	50	0.253	0.625	0.598	2.441
	1.5	40	0.217	0.659	0.426	2.446
	4	180	1.091	2.680	3.650	2.543
	4	145	1.091	2.418	3.455	2.468
	4	110	0.979	2.521	2.969	2.465
	4	90	0.913	2.045	2.720	2.447
	4	80	0.838	2.016	2.800	2.434
	4	50	0.353	0.622	0.748	2.459
	4	40	0.287	0.931	0.530	2.450
	8	180	1.292	3.403	4.276	2.462
	8	145	1.308	2.922	4.057	2.487
	8	110	1.135	2.828	3.469	2.462
	8	90	1.096	2.534	3.213	2.488
	8	50	0.521	1.007	1.263	2.461
	8	40	0.359	0.510	0.644	2.442