domain: phys spi

text:
Cécile Baron
email: cecile.baron@univ-amu.fr

Vu-Hieu Nguyen

Salah Naili

Carine Guivier-Curien

Interaction of ultrasound waves with bone remodelling: a multiscale computational study

Keywords: Bone remodelling, ultrasound stimulation, osteocyte, poroelasticity, wall shear stress

Ultrasound stimulation is thought to influence bone remodeling process. But recently, the efficiency of ultrasound therapy for bone healing has been questioned. Despite an extensive literature describing the positive effect of ultrasound on bone regenerationcell cultures, animal models, clinical studies -there are more and more reviews denouncing the inefficiency of clinical devices based on low intensity pulsed ultrasound stimulation (LIPUS) of the bone healing. One of the reasons to cause controversy comes from the persistent misunderstanding of the underlying physical and biological mechanisms of ultrasound stimulation of bone repair. As ultrasonic waves are mechanical waves, the process to be studied is the one of the mechanotransduction. Previous studies on the bone mechanotransduction have demonstrated the key-role of the osteocytes in bone mechanosensing. Osteocytes are bone cells ubiquitous inside the bone matrix, they are immersed in the interstitial fluid (IF) inside the lacuno-canalicular network (LCN). They are considered as particularly sensitive to a particular type of mechanical stress: wall shear stress on osteocytes due to the IF flow in the LCN. Inspired from these findings and observations, the

Introduction

Bone is a living tissue that is constantly being remodelled, adapting to its mechanical environment and capable of repair. Bone is one of few tissues that can heal without forming a fibrous scar. Trauma is the most expensive medical condition after heart conditions, costing the United States 56 billion dollars every year. Of that, 21 billion is used for the treatment of fractures. For these reasons, the efficacious and expedient treatment of fractures is of paramount importance to the patient, physician, and healthcare system as a whole [START_REF] Buza | Bone healing in 2016[END_REF]. Unfortunately, sometimes the healing process fails and non-unions or delayed unions occur. Delayed fracture healing and nonunion are observed in up to 5-10% of all fractures, and can present a challenging clinical scenario for the treating physician. The beneficial effect of ultrasound on bone remodelling was discovered in the 1950s [START_REF] Corradi | Effect of ultrasonics on the development of osseous callus in fractures[END_REF][START_REF] Yasuda | The classic: Fundamental aspects of fracture treatment by Iwao Yasuda, reprinted from[END_REF]. Since the 1980s, several authors have studied this phenomenon in different processes of adaptation of the bone to its environment: growth [START_REF] Duarte | The stimulation of bone growth by ultrasound[END_REF], targeted remodeling [START_REF] Chan | Low intensity pulsed ultrasound accelerated bone remodeling during consolidation stage of distraction osteogenesis[END_REF] and healing [START_REF] Schortinghuis | Ultrasound stimulation of maxillofacial bone healing[END_REF].

However, within the recent years, several systematic reviews and meta-analyses [START_REF] Busse | Trial to re-evaluate ultrasound in the treatment of tibial fractures (TRUST): a multicenter randomized pilot study[END_REF][START_REF] Schandelmaier | Low intensity pulsed ultrasound for bone healing: systematic review of randomized controlled trials[END_REF][START_REF] Griffin | Low intensity pulsed ultrasound for fractures of the tibial shaft[END_REF] provide moderate to high-certainty evidences against the contribution of USS in bone healing. The controversy on the subject continues and many questions remain unanswered [START_REF] Aspenberg | Comment to a BMJ EditorialIs LIPUS the baby in the bathwater?[END_REF][START_REF] Mortazavi | Mode & mechanism of low intensity pulsed ultrasound (LIPUS) in fracture repair[END_REF]. Actually, the real issue is that the underlying physico-biological mechanisms of the of ultrasound on bone remodeling remain unclear, as reported in several reviews [START_REF] Claes | The enhancement of bone regeneration by ultrasound[END_REF][START_REF] Martinez De Albornoz | The evidence of lowintensity pulsed ultrasound for in vitro, animal and human fracture healing[END_REF][START_REF] Padilla | Stimulation of bone repair with ultrasound: a review of the possible mechanic effects[END_REF]. One of the objectives of this work is to gain insight into the mechanical impact of ultrasound stimulation on bone from tissue to cell and explore its capacity to trigger bone remodeling. On one hand, clinical studies, animal models and cell cultures put on evidence an effect of LIPUS on bone repair, but the underlying physico-biological mechanisms remain unclear and feed an important literature as reported in several paper reviews [START_REF] Claes | The enhancement of bone regeneration by ultrasound[END_REF][START_REF] Martinez De Albornoz | The evidence of lowintensity pulsed ultrasound for in vitro, animal and human fracture healing[END_REF][START_REF] Padilla | Stimulation of bone repair with ultrasound: a review of the possible mechanic effects[END_REF]. On the other hand, the influence of a mechanical loading on bone remodelling has been widely explored by [START_REF] Cowin | Candidates for the mechanosensory system in bone[END_REF]; [START_REF] Weinbaum | A model for the excitation of osteocytes by mechanical loadinginduced bone fluid shear stresses[END_REF]; [START_REF] Claes | Magnitudes of local stress and strain along bony surfaces predict the course and type of fracture healing[END_REF]; [START_REF] Mullender | Mechanotransduction of bone cells in vitro: mechanobiology of bone tissue[END_REF]; [START_REF] Klein-Nulend | Mechanosensation and transduction in osteocytes[END_REF], for instance. Some of these works have revealed the osteocytes as the cornerstone of the mechanosensing and the pilots of the bone mechanotransduction. The osteocytes are ubiquituous bone cells sited in the lacuno-canalicular network (LCN), immerged in the interstitial fluid (IF) which fullfills the voids between osteocytes and the extra-cellular matrix (ECM). Authors consider that when mechanical loading is applied to bone at macroscale, it induces pressure gradients inside the bone, resulting in IF flow around the osteocytes which are in turn submitted to several mechanical loadings such as radiation forces, hydrostatic pressure or wall shear stress (WSS). The latter one has been introduced by [START_REF] Weinbaum | A model for the excitation of osteocytes by mechanical loadinginduced bone fluid shear stresses[END_REF] and experimentally demonstrated by [START_REF] Klein-Nulend | Sensitivity of osteocytes to biomechanical stress in vitro[END_REF] as the dominant factor, in bone mechanotransduction more than streaming for example [START_REF] Klein-Nulend | Mechanobiology of bone tissue[END_REF], 2013). Therefore, in this preliminary study we chose to focus on the WSS. The main assumption is that the bone remodelling is partly triggered by the level of fluid-induced wall shear stress especially on the process of the osteocytes [START_REF] Bakker | The production of nitric oxide and prostaglandin E2 by primary bone cells is shear stress dependent[END_REF][START_REF] Anderson | Nano-microscale models of periosteocytic flow show differences in stresses imparted to cell body and processes[END_REF].

A lot of numerical studies showed the impact of the physiological loading on the bone to interpret the mechanobiology of bone healing [START_REF] Claes | Magnitudes of local stress and strain along bony surfaces predict the course and type of fracture healing[END_REF][START_REF] Lacroix | A mechanoregulation model for tissue differentiation during fracture healing: analysis of gap size and loading[END_REF][START_REF] Isaksson | Comparison of biophysical stimuli for mechano-regulation of tissue differentiation during fracture healing[END_REF][START_REF] González-Torres | Influence of the frequency of the external mechanical stimulus on bone healing: a computational study[END_REF][START_REF] Gómez-Benito | Influence of high-frequency cyclical stimulation on the bone fracture-healing process: mathematical and experimental models[END_REF][START_REF] Nguyen | Numerical study of deformation-induced fluid flows in periodic osteonal matrix under harmonic axial loading[END_REF]Nguyen et al., , 2010a[START_REF] Nguyen | Influence of interstitial bone microcracks on strain-induced fluid flow[END_REF]. These numerical models usually considered a cyclic compressive loading applied to bone at frequency corresponding to daily activity (walk, run, jump) and related the mechanical stimulus induced at the cell scale to a biological response in terms of proliferation and differentiation of bone cells following the model of Prendergast and colleagues [START_REF] Prendergast | Biophysical stimuli on cells during tissue differentiation at implant interfaces[END_REF]. Ultrasound stimulation is a type of mechanical loading and as such may trigger reactions similar to physiological loading. However the characteristics of ultrasonic loading are very different from those of physiological loading and this raises a number of questions such as: do the mechanotransduction processes identified for physiological loading also act on osteocyte at US frequencies? Are there any others phenomena which can be implied in ultrasound stimulation? To understand how ultrasound could affect bone remodelling, a computational model is developed in order to estimate the level of mechanical stimulus applied on the osteocytes by the mesoscopic ultrasound stimulation. To our knowledge, it is the first numerical study on the multiscale interaction of ultrasound with the bone remodelling process. This paper does the link between biological observations and physical interpretation integrating both bone levels of porosity: the vascular network and the LCN network. However, the aim is not to infer the biological mechanisms acting inside the cell, but to concentrate on the characterization of the mechanical forces applied to the cell via the study of the WSS. To appropriately tackle the question, this mechanical model integrates two levels of porosity: the vascular network (Havers and Volkman canals -100 µm) reconstructed from CT scans and the lacuno-canalicular network containing the osteocytes (1 µm) taken into account through the poroelasticity theory. On the basis of previous works, this paper proposes to estimate the osteocyte wall shear stress identified as the relevant stimulus in bone mechanotransduction. To do this properly, one of the key parameters is the permeability of the LCN, the value of which remains controversial [START_REF] Cardoso | Advances in assessment of bone porosity, permeability and interstitial fluid flow[END_REF].

As a consequence, the goal of this paper is focused on the mechanical aspect of the bone mechanotransduction under ultrasound stimulation. The relevancy of the proposed mechanical model is tested looking at the influence of the LCN permeability value on the wall shear stress induced on the osteocytes by the ultrasound stimulation. The wall shear stress levels assessed for two permeability values are compared and discussed in relation to the results of the literature.

After this introduction, the paper is organized as follows. The materials and methods are presented in section 2. We give in subsection 2.1 the in vivo configuration, the geometry description is presented in subsection 2.2, the governing equations are given in subsection 2.3, the weak formulation is given in subsection 2.4 and the computational model is given in subsection 2.5. We show and discuss the main results on the estimation of wall shear stress applied to the osteocyte by ultrasonic stimulation in section 3. Section 4 finishes the manuscript with the conclusion.

Materials and Methods

In vivo configuration

Ultrasound stimulation is mainly dedicated to fracture healing in clinical applications, however, the goal of this preliminary study is to investigate the effect of ultrasound as a trigger for bone remodelling process and not its direct influence on the healing callus. Consequently, we choose to represent a slice of intact cortical bone which can be above or below the fracture site and still in interaction with ultrasonic waves since the US transducer has a diameter of 1 cm, greater than the fracture gap size (few mm). Considering that the remodelling would be managed by the osteocytes, the bone healing should be initiated from the intact parts of bone surrounding the fracture site. Moreover, the ultrasound stimulation of bone remodelling could be used for skeletal events other than fracture (localized bone loss, bone metastases, osteolytic tumors) and could also be represented by the model developed in this study. This model is therefore a 2D model considering a bone cross-section perpendicular to the bone axis. Cortical bone is considered as a biphasic medium (fluid-saturated poroelastic solid) taking into account two levels of porosity. The first one is the vascular porosity level (pores ∅ ∼ 100 µm). The vascular pores geometry is extracted from µCT images (pixel size = 22.5 µm). The vascular pores are supposed to be full of fluid. The second one is the LCN where the osteocytes are sited surrounded by the interstitial fluid. This level of porosity is constituted of an ubiquitous pore network of lacunae (pores ∅ ∼ 15 µm) and canaliculi (pores ∅ ∼ 0.5 µm). The size of these voids are too small to be extracted from the above-mentioned µCT images, consequently the bone matrix around the vascular pores is modeled as a poroelastic medium using Biot's theory to describe the wave propagation in an equivalent medium [START_REF] Biot | Theory of elasticity and consolidation for a porous anisotropic solid[END_REF][START_REF] Biot | The elastic coefficients of the theory of consolidation[END_REF]. To be close to the in-vivo configuration, the model takes into account the soft tissues around the bone considered as water in first approximation.

Geometry description

For this problem, we consider a two dimensional (2D) plane strain model as shown in Fig. 1. In this model, the poroelastic bone domain is denoted Ω b . The soft tissue, the marrow and the fluid-filled lacuna are all considered as fluid and denoted by Ω f 0 , Ω f j (j = 1, .., N -1) and Ω f N , respectively. The interfaces between the bone Ω b and the fluid domains Ω f j are denoted by Γ j (j = 0, .., N ). The external surface of the soft tissue domain Ω f 0 consists of 2 parts denoted by Γ F 0 (free surface) and Γ P 0 (imposed by a pressure), respectively.

Fig. 1: Geometry description

In what follows, we denote the coordinates of a point x by (x 1 , x 2 ) and the time by t. A superposed dot stands for time derivative, ∇ and ∇• stand for the gradient and divergence operators in 2D space, respectively. The symbol "•" denotes the scalar product and the symbol ":" between tensors of any order denotes their double contraction.

Governing equations

Equations in the fluid domains Ω f j

The fluid occupying the domain Ω f j is an acoustic fluid whose mass density and bulk modulus at equilibrium are denoted by ρ f j and K f j , respectively. By using the linear acoustic theory and by neglecting the body forces (other than inertia), the wave equation for the fluid in the domain Ω j may be expressed only in terms of the pressure field p f j (x, t) as follows:

1

K f j pf j - 1 ρ f j ∇ 2 p f j = 0, ∀x ∈ Ω f j (j = 0, ..., N ).
(1)

Note that the index j is used for indicating the number of the fluid domain, then there is no summation over j in the above equation as well as in the following.

The velocity vector v j (x, t) is relied to the gradient of pressure field p f j (x, t) following the Euler equation:

ρ f j vj + ∇p f j = 0, ∀ x ∈ Ω f j , (j = 0, ..., N ).
(2)

2.3.2 Equations in the anisotropic poroelastic bone (Ω b )
The cortical bone material is assumed to consist of a solid skeleton (with mass density ρ s ) and a connected pore network saturated by fluid (with mass density ρ f ). The Biot's anisotropic poroelastic model [START_REF] Biot | The elastic coefficients of the theory of consolidation[END_REF]Coussy, 2005;[START_REF] Thompson | A reformation of the equations of anisotropic poroelasticity[END_REF]) could be used to describe the wave propagation problem in bone [START_REF] Cowin | A recasting of anisotropic poroelasticity in matrices of tensor components[END_REF][START_REF] Nguyen | Simulation of ultrasonic wave propagation in anisotropic poroelastic bone plate using hybrid spectral/finite element method[END_REF].

Neglecting the body forces (other than inertia), the equations describing the wave propagation within the bone domain (Ω b ) read:

∇ • σ = ρ üs + ρ f ẅ , (3)
-∇p = ρ f üs + κ -1 ẇ + b ẅ , (4)
where ρ = φ ρ f + (1 -φ) ρ s is the mixture density, φ is the porosity, σ(x, t) is the total stress tensor and p(x, t) is the interstitial fluid pressure in the pores; the vectors of displacement of the solid skeleton and of the fluid are denoted by u s (x, t) and u f (x, t), respectively; the vector of relative displacement between the fluid and the solid frame weighted by the porosity is denoted by w = φ(u fu s ); κ is the permeability tensor defined by κ = k/µ where k is the intrinsic permeability tensor and µ is the fluid's dynamic viscosity; the tensor b is defined as b = (ρ f /φ) a where a is the tortuosity tensor. Note that both tensors k and b are symmetric second-order tensors.

The constitutive equations for the anisotropic linear poroelastic material are given by:

σ = C : ǫ -α p , (5)
- 1 M p = ∇ • w + α : ǫ , ( 6
)
where C is the elasticity fourth-order tensor of drained porous material; α, which is a symmetric second-order tensor, is the Biot effective tensor; M is the Biot scalar modulus; ǫ(x, t) is the infinitesimal strain tensor which is defined as the symmetric part of ∇u s .

Boundary and interface conditions

The surface of the soft tissue consists of 2 parts: Γ 0 = Γ P 0 ∪Γ F 0 which correspond to the zone excited by pulsed pressure and the free surface, respectively:

p f 0 (x, t) = P 0 g(t), ∀x ∈ Γ P 0 , (7)
p f 0 (x, t) = 0, ∀x ∈ Γ F 0 , (8)
where P 0 is the amplitude and g(t) is pulsed time function with a central frequency f c : g(t) = sin(2πf c t). At interfaces Γ j ((j = 0, ..., N )), the continuity of pressure and stress fields between the poroealstic medium and the fluid domains imposes:

p f = p f j , σn j = -p f j n j , ∀x ∈ Γ j (j = 0, ..., N ), (9)
where n j is the normal unit vector to Γ j pointing from the porous domain Ω b towards the fluid domain Ω f j (see Fig. 1). The periosteal interface (Γ 0 ) is assumed to be impervious [START_REF] Goulet | Influence of vascular porosity on fluid flow and nutrient transport in loaded cortical bone[END_REF], requiring:

w = 0, ∀x ∈ Γ 0 . ( 10
)
At the other interfaces Γ j (j = 1, ..., N ), open-pore condition is assumed, requiring:

ẇ • n j = (v j -us ) • n j , ∀x ∈ Γ j (j = 1, ..., N ). ( 11
)
In view of the Euler equation ( 2), the interface condition (11) may be rewritten as: The pressure field in the fluid occupying the domain Ω f j is described by Eq. ( 1). The weak form of this equation may be obtained by multiplying it by a scalar valued test function pf j (with support in Ω f j ) and then integrating over the domain Ω f j . By applying the Green-Gauss theorem and taking into account the boundary conditions ( 12), the weak formulation of the acoustic wave problem in Ω f j reads:

1 ρ f j ∇p f j + ẅ + üs • n j = 0, ∀x ∈ Γ j (j = 1, ..., N ). ( 12
Ω f j 1 K f j pf j pf j + Ω f j 1 ρ f ∇p f j • ∇p f j + Γj pf j ( ẅ + üs ) • n j = 0, ∀p f j ∈ C ad , (13)
where C ad is the admissible function space of the problem constituted by the sufficiently differentiable realvalued functions.

Note that the n j is pointing toward the interior of the fluid domain Ω f j .

Weak formulation in the bone domain Ω b

In order to resolve the time-dependent problem which is a high frequency problem, the finite element analysis were based on the (u sw) formulation. The balance equations of linear momentum (3)-( 4) are already written in terms of u s and w. Accordingly, constitutive equations ( 5)-( 6) may be restated as:

σ = C u ǫ + M α ∇ • w, (14)
p = -M (α : ǫ + ∇ • w) , (15)
where C u = C + M α ⊗ α is the undrained elasticity tensor. (The symbol ⊗ designates the tensorial product between two tensors.)

Let ũs and w be two vector valued test functions with support in Ω b . The weak form of balance equations ( 3)-( 4) may be obtained by multiplying them by ũs and w, respectively, integrating over Ω b and applying the Green-Gauss theorem (Nguyen et al., 2010b):

Ω b ũs • ρü s + Ω b ũs • ρ f ẅ + Ω b ∇ũ s : σ + j Γj ũs • p f j n j = 0, ∀ũ s ∈ C ad , ( 16
)

Ω b w • ρ f üs + Ω b w • (b ẅ) + Ω b w • κ -1 ẇ + Ω b ∇ w : (p I) + j Γj w • p f j n j = 0, ∀ w ∈ C ad , ( 17
)
where I is the second-order identity tensor.

Computational model

Ultrasound signal. The ultrasound signal emitted from Γ P 0 (see Fig. 1) is a pulsed signal with a central frequency f c =1 MHz, a duty cycle of 20% and a repetition frequency f r =1 kHz. The spatial average-time average acoustic intensity delivered, I SATA , is of 30 mW/cm 2 corresponding to an acoustic pressure amplitude P 0 =67 kPa (assuming pressure amplitude is constant over the transducer face). The transducer diameter is equal to 2 cm.

Material properties.

-Fluid phase: soft tissue, marrow and fluid inside the vascular pores are supposed to be perfect fluids mechanically equivalent to water: speed of sound c f =1500 m/s; mass density ρ f =1000 kg/m 3 ; bulk modulus K=2.25 GPa.

-Poroelastic phase: bone matrix around the vascular pores is considered as a transverse isotropic poroelastic medium. The LCN porosity is taken from the literature φ=0.05 [START_REF] Smit | Estimation of the poroelastic parameters of cortical bone[END_REF][START_REF] Lemaire | On the paradoxical determinations of the lacuno-canalicular permeability of bone[END_REF], the fluid inside the pores is considered as water with a dynamic viscosity µ=10 -3 Pa.s -1 .

Around the LCN, the extra-cellular matrix (ECM) has a mass density ρ ECM equal to 2 g/cm 3 . The model is constructed in two dimensional space in a plane perpendicular to the bone axis where the elastic properties of the ECM are assumed to be isotropic. The elasticity tensor C in this plane is expressed from the Voigt notation as [START_REF] Scheiner | Poromicromechanics reveals that physiological bone strains induce osteocyte-stimulating lacunar pressure[END_REF]:

C =   C 1111 C 1122 0.0 C 1122 C 1111 0.0 0.0 0.0 C 1212   , (18)
with C 1212 = 1 2 (C 1111 -C 1122 ). The values used in the simulations are C 1111 = 22.88 GPa and C 1122 = 10.14 GPa which lead to the Young modulus E = 16.56 GPa and Poisson ratio ν = 0.308. The Biot's parameters are calculated from this elasticity tensor. The component values of the Biot effective tensor and the Biot scalar modulus are given by α 11 = α 22 = 0.15 and M = 35.6 GPa (for more details see [START_REF] Rosi | Numerical investigations of ultrasound wave propagating in long bones using a poroelastic model[END_REF]). Concerning the permeability of this porous medium, two values are investigated. Inspired from the permeability estimation of the vascular porous network in cortical bone [START_REF] Zhang | Estimates of the peak pressures in bone pore water[END_REF][START_REF] Swan | Micromechanically based poroelastic modeling of fluid flow in haversian bone[END_REF], a first value of LCN permeability is calculated from the Kozeny-Carman relation, assuming a Poiseuille flow in a network of cylindrical capillaries aligned in the direction e i of radius R i = 50 × 10 -9 m and representing a porosity φ=5%: k KC = φR 2 i /8 = 1.56 × 10 -17 m 2 . A second one is taken from literature [START_REF] Smit | Estimation of the poroelastic parameters of cortical bone[END_REF], estimation based on the best fit between finite element predictions and experimental data on canine bone:

k S = 2.2 × 10 -22 m 2 .
Computation parameters. The convergence of the numerical results was tested to choose a mesh size (∆x) smaller than the tenth of the wavelength λ in each material. To ensure a relevant sampling of the timedependent phenomena, the solver time stepping ∆t was fixed to T/40 where T=1/f c (f c is the central frequency of pulsed time function and ∆t = 2.5 × 10 -8 s). Note that this time step value verifies the calculation stability condition ∆t < ∆x √ 2Vmax , where V max is a maximum velocity estimated. Numerical simulations were run with the commercial finite-element software Comsol Multiphysics.

Results and Discussion

As it has been widely reported in litterature, the osteocytes are the cornerstone of bone remodelling. One of their main characteristics is that they are mechanosensitive cells (see for instance [START_REF] Bonewald | The amazing osteocyte[END_REF]). It has been demonstrated that the dendritic processes of the osteocytes are more sensitive to mechanical loading than the cell body [START_REF] Han | Mechanotransduction and strain amplification in osteocyte cell processes[END_REF][START_REF] Burra | Dendritic processes of osteocytes are mechanotransducers that induce the opening of hemichannels[END_REF] and in particular to fluid shear stress [START_REF] Bakker | The production of nitric oxide and prostaglandin E2 by primary bone cells is shear stress dependent[END_REF][START_REF] Anderson | Nano-microscale models of periosteocytic flow show differences in stresses imparted to cell body and processes[END_REF]. Inspired from these assumptions, this study focuses on the IF wall shear stress (WSS IF ) induced by the ultrasound stimulation into the canaliculi. However, the present poroelastic approach does not allow to directly assess WSS IF acting on individual osteocytes, but we can still calculate the WSS IF from the velocity of the interstitial fluid relative to the ECM obtained from the Biot's theory ( ẇ). We do this using two different models of fluid flow associated with two different permeability values of the LCN. KC-model. In line with the first permeability value studied (k KC ), the wall shear stress τ KC can be firstly assessed using Kozeny-Carman theory, through equation ( 19):

τ KC = R i µ| ẇ| 2k KC , ( 19
)
with R i is the typical radius of canaliculi in the LCN (50×10 -9 m), µ is the IF dynamic viscosity (10 -3 Pa.s -1 ), k KC is the LCN permeability (1.56×10 -17 m 2 ) and ẇ is the velocity of the interstitial fluid relative to the ECM calculated from the Biot's theory (see §2.3.2). Equation ( 19) assumes that the LCN is modelled as parallel straight tubes full of fluid embedded in a solid matrix.

WT-model. Wall shear stress levels can also be defined from [START_REF] Wang | Modeling interstitial flow in an artery wall allows estimation of wall shear stress on smooth muscle cells[END_REF] as:

τ WT = µ| ẇ| √ k S , (20)
where k S is the LCN permeability equals to 2.2×10 -22 m 2 [START_REF] Smit | Estimation of the poroelastic parameters of cortical bone[END_REF]. This equation given in [START_REF] Goulet | Influence of vascular porosity on fluid flow and nutrient transport in loaded cortical bone[END_REF] considers the fluid flow in an annular space through a net of transverse fibers taking into consideration the presence of the osteocyte process inside the canaliculus and the glycosaminoglycans (GAG) fibers which are transverse fibrils which anchor and center the cell process in its canaliculus (see Fig. 2).

KC-model vs WT-model

We calculate the moving average of the WSS IF over a signal period (1 µs). Data are smoothed with the "loess" method (LOcally Estimated Scatterplot Smoothing), which is a method using linear regression of smooth data. The method is altered to assign a zero weight to data outside six mean absolute standard deviations. We compare the time evolution of WSS IF (τ KC or τ WT ) over 3 cycles for both LCN permeability values k KC and k S (see Fig. 3). In what follows, the quantities τ KC and τ WT are evaluated in points shown on the bone geometry.

Firstly, the computational model developed in this study provides results on WSS IF induced by ultrasound stimulation at cell scale which can be considered as reasonable. Comparing both models (KC and WT) used to represent the canalicular flow, a very significant difference can be observed: the WSS IF levels τ KC , range from 0 to 600 Pa i.e. 3 orders of magnitude higher than those of τ WT between 0 and 0.6 Pa. It is therefore clear that the way to calculate WSS IF in the porous medium has a great influence on the results. In comparison with literature, only the WT-model provides WSS IF values in agreement with a potential mechanical signal to activate osteocyte: 0.8<τ <3 Pa [START_REF] Weinbaum | A model for the excitation of osteocytes by mechanical loadinginduced bone fluid shear stresses[END_REF], even if they are close to the lower limit of the interval. These results confirm that the KC-model, bundle of aligned cylindrical pores full of fluid, neglecting the presence of the process, is inappropriate to correctly represent the interstitial fluid behavior in the LCN, as it has already been proven by [START_REF] Weinbaum | A model for the excitation of osteocytes by mechanical loadinginduced bone fluid shear stresses[END_REF]. They showed that one has to consider the fibers surrounding the cell process. The WT-model, firstly developed to model the flow through muscle cells, has been already used for modeling the fluid shear stress inside the LCN of the cortical bone [START_REF] Goulet | Influence of vascular porosity on fluid flow and nutrient transport in loaded cortical bone[END_REF] using the analogy of the system osteocyte process and GAG fibers with muscle cell and muscle fibers but also in scaffolds designed for bone regeneration [START_REF] Ouyang | Hydromechanical mechanism behind the effect of pore size of porous titanium scaffolds on osteoblast response and bone ingrowth[END_REF][START_REF] Ali | Permeability and fluid flow-induced wall shear stress of bone tissue scaffolds: computational fluid dynamic analysis using newtonian and non-newtonian blood flow models[END_REF] including under LIPUS [START_REF] Wu | Enhancement of cell ingrowth, proliferation, and early differentiation in a three-dimensional silicon carbide scaffold using lowintensity pulsed ultrasound[END_REF].

According to these observations, only results from the WT-model are analysed further.

Focus on WT-model

WSS IF distribution over 10 cycles. The stimulation duration is extended to 10 cycles (10 ms) for the WTmodel. On Fig. 4, the WSS IF distribution is represented for three different times: 1 ms, 5 ms and 10 ms. The most stimulated areas are located around the endosteum and the vascular pores interfaces. Globally the WSS IF levels increase with time. To further explore evolution in time of the WSS IF potentially experienced by the osteocytes, τ WT is evaluated on 12 points distributed on the bone surface during 10 cycles (see Fig. 5). The curves are obtained in the same procedure as described in § 3.2. Some points, located near the endosteum interface and nearby the vascular pores, reach higher τ WT values than the six points investigated in Fig. 3: 0 < τ WT < 1.5 Pa, still in the range defined by [START_REF] Weinbaum | A model for the excitation of osteocytes by mechanical loadinginduced bone fluid shear stresses[END_REF]. This WSS IF distribution may raise questions about the possible role of lining cells into the bone remodelling stimulated by ultrasound. Lining cells are bone cells covering over 90% of the surfaces of adult bone, and which are also recognized as mecanosensory cells [START_REF] Mullender | Osteocytes and bone lining cells: which are the best candidates for mechano-sensors in cancellous bone?[END_REF][START_REF] Robling | Mechanical Signaling for Bone Modeling and Remodeling[END_REF]. They are liable to play a role in bone remodelling triggering, especially since ultrasonic stimulation transmits compression waves perpendicular to the bone axis and thus a form of mechanical stress different from physiological loading, preferentially oriented parallel to the bone axis.

The obtained values (see Fig. 5) and the spatial heterogeneity of the WSS IF distribution (see Fig. 4) are in agreement with the results reported in [START_REF] Fan | A multiscale 3d finite element analysis of fluid/solute transport in mechanically loaded bone[END_REF] on mouse bone with a similar canalicular flow model as the WT-model, even if the loading conditions are not the same, in particular, the frequency loading (around 0.5 Hz) which is lower than the pulsed ultrasound frequencies.

Another factor that can influence the distribution of the mechanical stimulus is absorption that could lead to wave attenuation. The Biot's poroelastic model used in this work to represent the interaction of the LCN with ultrasonic waves assumes that the solid phase is elastic and the fluid is viscous (µ = 10 -3 Pa/s). This absorption is taken into account via the fluid's dynamic viscosity in the permeability tensor (α). On the other hand, absorption related to the solid phase of the porous medium is neglected in this first approach.

Frequency loading dependency. Indeed, the frequency of cyclic loading is regarded as a decisive factor for bone formation [START_REF] Weinbaum | A model for the excitation of osteocytes by mechanical loadinginduced bone fluid shear stresses[END_REF][START_REF] Bacabac | Nitric oxide production by bone cells is fluid shear stress rate dependent[END_REF][START_REF] Klein-Nulend | Mechanobiology of bone tissue[END_REF][START_REF] González-Torres | Influence of the frequency of the external mechanical stimulus on bone healing: a computational study[END_REF][START_REF] Gómez-Benito | Influence of high-frequency cyclical stimulation on the bone fracture-healing process: mathematical and experimental models[END_REF]. It is noteworthy that in our case, two frequencies can be considered and can play a role in the bone cell stimulation: the frequency of the ultrasonic signal, 1 MHz; and the pulse repeti-tion frequency, 1 kHz. Both are much higher than the physiological loading frequencies (1-20 Hz). However, it has been demonstrated that a same load applied on bone at mesoscopic scale, at high frequency, could induce sufficient stress on the osteocytes to activate a biological response when it would have had no effect at low frequencies [START_REF] Weinbaum | A model for the excitation of osteocytes by mechanical loadinginduced bone fluid shear stresses[END_REF][START_REF] Han | Mechanotransduction and strain amplification in osteocyte cell processes[END_REF]. However the frequency range (0-100 Hz) of these observations remains below the two frequencies considered here. Following [START_REF] Bacabac | Nitric oxide production by bone cells is fluid shear stress rate dependent[END_REF], one of the determinant parameter for osteocytes response is the rate of fluid shear stress (the product of fluid shear stress amplitude in Pa and frequency loading in Hz). In their work on cell culture in a flow chamber, they assumed that low-magnitude, high-frequency mechanical stimuli may be as stimulatory as high-amplitude, low-frequency stimuli. In the present work, assuming that the trend remains similar for frequencies beyond 9 Hz (maximum frequency studied in Bacabac et al. ( 2004)), it would mean that considering the signal frequency of 1 MHz, a fluid shear stress magnitude of 50 µPa could be sufficient to activate an osteogenic response, and considering the pulse repetition frequency of 1 kHz, the fluid shear stress magnitude reliable to induce osteocyte response should be around 0.05 Pa. In both cases, τ WT is higher than these levels of fluid shear stress.

Stimulation duration. One of the limitations of this study is that the maximum stimulation duration considered is 10 ms which is much shorter than the daily treatment time (20 min) defined in the fracture repair protocols [START_REF] Poolman | Low intensity pulsed ultrasound (LIPUS) for bone healing: A clinical practice guideline[END_REF]. Consequently, one of the points to check is that there are no cumulative effects during stimulation cycles. On Fig. 5, it can be seen that values obtained on the chosen points are quite stable over time. This trend is confirmed on Fig. 6 which shows the value of τ W T averaged per cycle of 1 ms for each considered point. Nevertheless, further works are needed to resolve technical issues (simulation time and file size) and to look at longer-time phenomena.

Conclusion

This work provides a numerical model of the interaction of low pulsed intensity ultrasound and cortical bone as a two-level porous medium. The main objective of this study is to gain insight into the process of mechanotransduction induced by ultrasound stimulation. Considering the osteocytes as the pilot of the bone remodelling under mechanical stimuli, the modelling of the LCN was essential. The reconstruction of this porosity network from images not being available over the entire surface considered, the model uses the theory of poroelasticity to assess the physical quantities relative to this porosity scale and estimate the effect of ultrasonic waves at the cell scale. The results obtained tend to prove the potential of the model detailed in this study to represent the interaction of ultrasound stimulation and cortical bone taking into account two levels of porosity. Furthermore, they highlight how the value of the LCN permeability could affect the WSS IF levels.

According to literature, it tends to demonstrate that the KC-model is not relevant to represent the fluid flow into the LCN and that the WT-model seems to be more appropriate to evaluate the WSS IF levels potentially experienced by the osteocytes. Nevertheless, most of the assumptions have been interpolated from previous works on bone mechanotransduction under physiological loading, more data are needed on osteocyte activation under higher frequency stimulation to identify the physical and biological processes at work in the context of pulsed ultrasound stimulation.

Further studies are on going to investigate the influence of US parameters on the mechanical signal received by osteocytes and explore the effect of the multiscale bone anisotropy through a 3D-model including the absorption of the solid phase of the bone matrix. Investigations are under progress to experimentally validate the model.

  formulation in the fluid domains Ω f j

Fig. 2 :

 2 Fig. 2: Schematic methodology to calculate WSS IF from the output parameter ẇ given by the poroelastic theory implemented in the numerical model. Two assumptions are investigated: the KC-model of the fluid flow inside the canaliculus, ignoring the existence of the osteocyte process, surrounded by the ECM following the Kozeny-Carman relation (left) and the WT-model including the process and the GAG fibers immersed in interstitial fluid to match up to relation given by Wang and Tarbell (1995) (right). The osteocytes is in light gray, the interstitial fluid in dark gray and the GAG fibers in black; the ECM is hatched.

Fig. 3 :Fig. 4 :

 34 Fig. 3: WSS IF estimated in 6 points over 3 cycles (3 ms) for (a) KC-model and (b) WT-model.

Fig. 5 :

 5 Fig. 5: WSS IF estimated for WT-model in 12 points over 10 cycles (10 ms).

Fig. 6 :

 6 Fig. 6: WSS IF averaged over each of the ten first cycles.