Fodil Meraghni

Qiang Chen

George Chatzigeorgiou
email: georges.chatzigeorgiou@ensam.eu

Extended mean-field homogenization of unidirectional piezoelectric nanocomposites with generalized Gurtin-Murdoch interfaces

Keywords: Extended Mori-Tanaka Homogenization, Surface Piezoelectricity, Nanocomposites, Gurtin-Murdoch Interface, Electromechanical Loading

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

Introduction

The role of interface/surface in affecting the stiffness and strength of heterogeneous microstructures is well recognized. The interface can have an important contribution to not only the homogenized moduli but also the local stress field distribution, thereby impacting the plastic flow and failure modes [START_REF] Barenblatt | The formation of equilibrium cracks during brittle fracture. General ideas and hypotheses. Axially-symmetric cracks[END_REF][START_REF] Needleman | A continuum model for void nucleation by inclusion debonding[END_REF][START_REF] Tu | Cohesive zone-based damage evolution in periodic materials via finitevolume homogenization[END_REF][START_REF] Benveniste | Imperfect soft and stiff interfaces in two-dimensional elasticity[END_REF][START_REF] Hashin | Thin interphase/imperfect interface in elasticity with application to coated fiber composites[END_REF][START_REF] Jiang | Physics-informed deep neural network enabled discovery of size-dependent deformation mechanisms in nanostructures[END_REF][START_REF] Wang | Interface cracks with surface elasticity in anisotropic bimaterials[END_REF]. In many cases, the assumption of the perfect interface in composite materials considering the inclusion and the matrix are well-bonded, cannot adequately describe the physical nature of the interface region. Therefore, different interface models have been developed according to their kinematic and kinetic characteristics [START_REF] Chatzigeorgiou | Generalized interfacial energy and size effects in composites[END_REF][START_REF] Firooz | Understanding the role of general interfaces in the overall behavior of composites and size effects[END_REF][START_REF] Firooz | Homogenization accounting for size effects in particulate composites due to general interfaces[END_REF][START_REF] Benveniste | A general interface model for a three-dimensional curved thin anisotropic interphase between two anisotropic media[END_REF][START_REF] Tu | Evolution of interfacial debonding of a unidirectional graphite/polyimide composite under off-axis loading[END_REF][START_REF] Buryachenko | Local and Nonlocal Micromechanics of Heterogeneous Materials[END_REF], aiming at better identifying the underlying stress transfer mechanism across the interface region and eventually capturing the macroscopic behavior of the microstructures.

Surface stresses in solid bodies, which are the focus of the present work, are caused by the differences in configuration and in coordination numbers between atoms at the surface and in the bulk materials [START_REF] Gibbs | The collected works of J. Willard Gibbs[END_REF]. They are generally non-hydrostatic and non-homogeneous Moreover, in material processing, the surface of a material is frequently exposed to, e.g., oxidation, aging, grit blasting, plasma jet treatment, etc., thus possessing different properties from its adjacent bulk material [START_REF] Javili | A finite element framework for continua with boundary energies. Part I: The two-dimensional case[END_REF]. In nanocomposite microstructures studied in this manuscript, the molecular dynamics simulations [START_REF] Chen | Finite-volume homogenization and localization of nanoporous materials with cylindrical voids. Part 2: New results[END_REF] reveal that the large surface-to-volume ratio could cause significant interfacial stress and energy along the surface that stretches between the inclusion and the matrix.

A related feature of the nanostructures with surface stress and energy is the fact that the mechanical behaviors are size-dependent. More specifically, the local and homogenized responses of nanocomposites depend on the properties and curvature of the interface [START_REF] Chen | Derivation of the generalized Young-Laplace equation of curved interfaces in nanoscaled solids[END_REF][START_REF] Kushch | Elastic fields and effective moduli of particulate nanocomposites with the Gurtin-Murdoch model of interfaces[END_REF]. The importance of surface stress and energy in nanocomposite structures has attracted considerable research attention over the past two decades through extensive atomistic simulations in the literature [START_REF] Kushch | Atomistic vs. continuum models of nanoporous elastic solid: Stress fields, sizedependent effective stiffness and surface constants[END_REF][START_REF] Davydov | On molecular statics and surface-enhanced continuum modeling of nano-structures[END_REF]. From a continuum mechanics approach's point of view, surface stresses can be phenomenologically interpreted as an elastic membrane equipped with its own potential energy.

The surface elasticity model developed by Gurtin and Murdoch [START_REF] Gurtin | A continuum theory of elastic material surfaces[END_REF][START_REF] Gurtin | Surface stress in solids[END_REF] is widely accepted in the literature. The Gurtin-Murdoch interface model views the interface as an infinitesimally thin sheet between the inclusion and the matrix, characterized by its own elastic moduli. The interface is kinematically coherent, indicating the tangential strain components are equal to strains on the adjacent sides of the bulk materials.

The effect of surface stress and energy on the effective moduli and elastic fields of nanocomposites with the Gurtin-Murdoch interface has been extensively quantified through various micromechanics and homogenization approaches, defined by the representative volume element (RVE) concept for statistically homogeneous microstructures and the repeating unit cell (RUC) concept for periodic microstructures [START_REF] Pindera | Micromechanics of spatially uniform heterogeneous media: A critical review and emerging approaches[END_REF][START_REF] Chen | Extended mean-field homogenization of viscoelastic-viscoplastic polymer composites undergoing hybrid progressive degradation induced by interface debonding and matrix ductile damage[END_REF]. Miller and Shenoy [START_REF] Miller | Size-dependent elastic properties of nanosized structural elements[END_REF] obtained the surface elastic constants of single-crystal metals with two different crystallographic orientations using molecular dynamic simulations. Their contributions have spurred extensive use of the Gurtin-Murdoch interface model in the literature. For instance, employing the modified Mori-Tanaka approach, Duan and coworkers [START_REF] Duan | Thermo-elastic properties of heterogeneous materials with imperfect interfaces: Generalized Levin's formula and Hill's connections[END_REF][START_REF] Duan | Nanoporous materials can be made stiffer than non-porous counterparts by surface modification[END_REF] were perhaps the first to investigate systematically the homogenized elastic moduli of the nanoporous aluminum with the Gurtin-Murdoch interface model. Following this work, Chen et al. [START_REF] Chen | Size-dependent elastic properties of unidirectional nanocomposites with interface stresses[END_REF] developed the composite cylinder assemblage and the generalized self-consistent method with energetic surfaces. More recent contributions are the elasticity-based homogenization approaches such as the equivalent inhomogeneity method, Mogilevskaya et al. [START_REF] Mogilevskaya | Equivalent inhomogeneity method for evaluating the effective elastic properties of unidirectional multi-phase composites with surface/interface effects[END_REF], the locally-exact homogenization theory, Wang et al. [START_REF] Wang | Homogenized moduli and local stress fields of unidirectional nano-composites[END_REF]. The finiteelement, Javili et al. [START_REF] Javili | Relationships between the admissible range of surface material parameters and stability of linearly elastic bodies[END_REF], Chatzigeorgiou et al. [START_REF] Chatzigeorgiou | Generalized interfacial energy and size effects in composites[END_REF], Chen and Pindera [START_REF] Chen | Homogenization and localization of elastic-plastic nanoporous materials with Gurtin-Murdoch interfaces: An assessment of computational approaches[END_REF], Kushch [START_REF] Kushch | Atomistic vs. continuum models of nanoporous elastic solid: Stress fields, sizedependent effective stiffness and surface constants[END_REF], and finitevolume-based, Chen et al. [START_REF] Chen | Homogenization and localization of elastic-plastic nanoporous materials with Gurtin-Murdoch interfaces: An assessment of computational approaches[END_REF][START_REF] Chen | Homogenization and localization of nanoporous composites-A critical review and new developments[END_REF], numerical homogenization theories are also favorable techniques due to their flexibility to handle more complicated microstructures and phase constitutive laws. All the above studies indicate that depending on whether the surface energy is positive or negative, using the Gurtin-Murdoch interface model, nanoporous composites may produce either enhancement or reduction of the homogenized moduli when the pore radii are sufficiently small (typically below 20 nm), relative to materials with the same porosity volume fraction but much larger pore sizes. The readers are referred to the recent articles for more complete references in this area, cf. Wang et al. [START_REF] Wang | Surface stress effect in mechanics of nanostructured materials[END_REF], Chen et al. [START_REF] Chen | Homogenization and localization of nanoporous composites-A critical review and new developments[END_REF], Firooz et al. [START_REF] Firooz | Homogenization of Composites With Extended General Interfaces: Comprehensive Review and Unified Modeling[END_REF], Mogilevskaya et al. [START_REF] Mogilevskaya | Fiber-and Particle-Reinforced Composite Materials With the Gurtin-Murdoch and Steigmann-Ogden Surface Energy Endowed Interfaces[END_REF].

Apart from the elastic behavior, there is a growing interest in developing micromechanics models for understanding the multiphysics response of piezoelectric/dielectric nanostructures due to their promising applications in nanoactuators and nanogenerators [START_REF] Yu | Monitoring on triboelectric nanogenerator and deep learning method[END_REF][START_REF] Yang | Particle-Laden Droplet-Driven Triboelectric Nanogenerator for Real-Time Sediment Monitoring Using a Deep Learning Method[END_REF][START_REF] Isaeva | Lead-free 0-3-type composites: From piezoelectric sensitivity to modified figures of merit[END_REF][START_REF] Thakor | Dielectric spectroscopy and structural characterization of nano-filler-loaded epoxy resin[END_REF][START_REF] Thabet | Assessment of dielectric strength and partial discharges patterns in nanocomposites insulation of single-core power cables[END_REF]. Such work has the potential to facilitate the efficient design and analysis of nano-sized functional components/devices.

In piezoelectric nanostructures, in addition to surface stress, the surface energy can also depend on the surface electric displacement [START_REF] Hoang | First-principles based multiscale model of piezoelectric nanowires with surface effects[END_REF][START_REF] Dai | Surface piezoelectricity: Size effects in nanostructures and the emergence of piezoelectricity in non-piezoelectric materials[END_REF]. Thus far, relatively few works have been devoted to incorporating surface piezoelectricity into the homogenization techniques albeit many micromechanics models are already available in the literature with multiphysics capabilities. An exception is a recent work by Chen et al. [START_REF] Chen | Homogenization of size-dependent multiphysics behavior of nanostructured piezoelectric composites with energetic surfaces[END_REF] who investigated in a systematical manner a full set of homogenized properties and local stress/electric fields of piezoelectric nanoporous composites with incorporated surface piezoelectricity for various pore radii and volume fractions using the multiphysics CCA/GSCM and finite-element as the gold solution. Prior to this study, Chen [START_REF] Chen | Exact size-dependent connections between effective moduli of fibrous piezoelectric nanocomposites with interface effects[END_REF] developed the exact size-dependent connections between effective moduli of piezoelectric nanocomposites with interface effects using the generalized Hill's connections and Levin's formula. Xiao et al. [START_REF] Xiao | Size-dependent effective electroelastic moduli of piezoelectric nanocomposites with interface effect[END_REF] studied the effective electroelastic moduli of nanocomposites using the generalized self-consistent approach. Gu and Qin [START_REF] Gu | Variational principles and size-dependent bounds for piezoelectric inhomogeneous materials with piezoelectric coherent imperfect interfaces[END_REF] established the bounds for the effective properties of an inhomogeneous material with coherent imperfect interfaces in the setting of piezoelectricity using the classical minimum potential principles of linear piezoelectricity.

In this contribution, a generalized Mori-Tanaka method enhanced with surface piezoelectricity has been developed for the first time for identifying the homogenized moduli of unidirectional piezoelectric nanocomposites. The bulk material is transversely isotropic while the interface is modelled using the extended Gurtin-Murdoch model considering the in-plane surface stress and electric displacement. Kinematically coherent interface is considered wherein the mechanical displacement and electric potential continuity conditions hold between the fiber and matrix phase while the interfacial stress and electric displacement jump are governed by the surface equilibrium and electrostatic conditions, also known as the Young-Laplace equations [START_REF] Chen | Derivation of the generalized Young-Laplace equation of curved interfaces in nanoscaled solids[END_REF].

Eshelby's inhomogeneity problems are solved analytically for four systems of equations, generated from the known macroscopic loading conditions, thereby permitting the exact computation of the dilute concentration piezoelectric tensors. The accuracy of the extended theory is validated extensively vis-à-vis the CCA/GSCM and finite-element predictions of Chen et al. [START_REF] Chen | Homogenization of size-dependent multiphysics behavior of nanostructured piezoelectric composites with energetic surfaces[END_REF].

The rest of the manuscript is organized as follows: Section 2 recalls the constitutive relation of the piezoelectric materials and introduces the generalized Gurtin-Murdoch coherent interface model. Sections 3 and 4 describe the generalized Eshelby's problem and Mori-Tanaka homogenization for inhomogeneity enhanced with energetic surfaces, respectively. Section 5 derives the analytical solutions for the determination of dilute concentration tensors for electromechanical fields under known macroscopic loading. Numerical results are generated in Section 6 to demonstrate the capabilities of the extended theory for accurate identification of the homogenized moduli under various loading conditions. Section 7 draws pertinent conclusions.

Useful equations for the analytical formalisms are provided in the Appendices.

Preliminaries

For a two-phase piezoelectric composite in which nano-inhomogeneities are randomly dispersed in the matrix phase, as shown in Figure 1, assuming small deformations/small rotations, electrostatic conditions with no electric charge, the strain tensor and the electric field vector

  q ε
in the phase in the Cartesian coordinates can be expressed with the help of the

  q e th q   1 2 3
, , x x x displacement vector and the electric scalar potential respectively through the vector-type

  q u   q  equations: (1)
            T 1 , 2        ε u u q q q q q
e where indicates the gradient operator. The superscripts denote the matrix and  0,1 q  fiber/porosity, respectively. , , , , ,

q q        σ      T 1 2 3 , , q q  d d d d
electrostatic equations:

(2)

    0, 0       σ q q d
where represents the divergence operator.   A linear piezoelectric, transversely isotropic material obeys the following constitutive law, written in matrix-type notation:

(3)

                    T q q q q q q q q q q         σ L ε e e ε κ e e d
where , and denote the 4 th order elastic, 2 nd order dielectric, and 3 rd order piezoelectric

  q L   q κ   q e
tensors for the phase, respectively, whose matrix-type forms are given below for a transversely th q isotropic material: , ,

    0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 q q k m k m l k m k m l l l n m m p                        L    
                 e     11 11 33 0 0 0 0 0 0 q q               κ
where , , , and are Hill's moduli [START_REF] Hill | Theory of mechanical properties of fiber-strengthened materials: I. elastic behavior[END_REF] for plane strain bulk modulus, cross modulus, k l n m p unit uniaxial straining modulus, transverse shear modulus, and axial shear modulus, respectively.

      q q q   E   (4)
where is the symmetric matrix. and are the

T         L e e κ =  9 9    T ,  σ d    T ,   E ε e 12 1 
vectors. Also, the generalized vector (5)

     T 1 2 3
, , ,

q q u u u   U is introduced.
In addition to the bulk materials, it is assumed that an infinitesimally thin energetic interface between the fibers and the matrix can be modelled by the generalized Gurtin-Murdoch model [START_REF] Chen | Exact size-dependent connections between effective moduli of fibrous piezoelectric nanocomposites with interface effects[END_REF]. For an energetic interface with a cylindrical shape of a constant radius , the interface a equations are better represented in cylindrical coordinates . If the normal vector of the

 

, , r z  interface is parallel to the local axis , the constitutive law of the surface reads:

r (6)     ( ) s s s   E  
where the superscript denotes an energetic surface and s (7)

         
                                                                Ε  e e d d



The interface is assumed to be kinematically coherent, indicating no atomic bonds are broken along the surface separating the two phases. The coherency of the surfaces can be interpreted as the tangential and axial components of the strain and electric field, in the cylindrical coordinates , to be continuous across the interface as follows:

  , , r z  (8)                           0 1 , , , , , , , 2 , , 2 , , 2 , , , , , , , , s
zz zz zz z z z z z z a a z a z a z a z a z a z a z a z a z                                                                                     e e e
e e e

The interface equilibrium and electrostatic equations are written as:

(9)                       0 1 0 1 0 1 0 s rr s s z r s s z zz rz s s z r a a z a z a z                                          d d d
where

. The displacement and electric potential continuity

              0 1                 
conditions are expressed through the relations (10)

                1 0 1 0 , , , , , , , , , , , i
i u a z u a z i r z a z a z          

Eshelby's problem for inhomogeneity with energetic interface

Identifying the overall properties of a composite sensitive to electromechanical loading conditions is the task of homogenization strategies and techniques. In this manuscript, the mean-field approach of the Mori-Tanaka method is followed. This method relies on two aspects: i) solving a basic (Eshelby-type) problem of an inhomogeneity embedded in an infinite medium, allowing to compute the so-called ''dilute concentration'' tensors, and ii) using this information at the level of the composite for identifying the ''concentration'' tensors, whose knowledge leads to link the microscopic and macroscopic fields and to obtain the macroscopic matrix of the overall

 =  composite.
Consider a cylindrical inhomogeneity, embedded in an infinite medium, as shown in Figure 2. The medium and the inhomogeneity are characterized by constant generalized moduli and

  0 
respectively. The inhomogeneity occupies the space with volume and is bounded by

  1  1  1 V
the energetic interface with the unit normal vector . The medium occupies the space ,  n 0  which is extended to infinity (boundary surface ). At far distance, a linear field

  (11)
0 0 0 0 11 1 12 2 13 3 1 0 0 0 0 12 1 22 2 23 3 2 0 0 0 0 13 1 23 2 33 3 3 0 0 0 0 1 1 2 2 3 3 ext x x x u x x x u x x x u x x x                                                U e e e
with , constant values applied.

0 ij  0 i e   , 1, 2,3 i j 
For this problem, which is a generalized version of the famous Eshelby inhomogeneity problem [START_REF] Eshelby | The determination of the elastic field of an ellipsoidal inclusion, and related problems[END_REF] accounting for multiphysics phenomena, the constitutive law is position-dependent and reads

(12)           0 0 1 1 , ,             E E    x x x x x
The boundary conditions correspond to uniform at a far distance. The medium is subjected to linear displacement and electric potential at a far distance

For the pure mechanical problem, Gu et al. [START_REF] Gu | Size-dependent effective elastic moduli of particulate composites with interfacial displacement and traction discontinuities[END_REF], based on the methodology of Duan et al. [START_REF] Duan | A unified scheme for prediction of effective moduli of multiphase composites with interface effects. Part I: Theoretical framework[END_REF], proposed to substitute the inhomogeneity/interface system with an equivalent particle, using energy principles. Instead of that, we propose to identify dilute concentration tensors according to the formulas:

(13) 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0

    1 1 0 0 ,         E T E H E  with (14) 1 T 1 1 dS       E U   (15) 1 1 ( ) 1 1 1 1 d d s V S               and (16)
n n n n n n n n n n n n             



where the superscripts " " and " " signs represent the fiber and matrix side of the interface,   respectively. We note that while represents the classical dilute strain concentration tensor, T H differs from the definition of a classical dilute stress concentration tensor [START_REF] Böhm | A short introduction to basic aspects of continuum micromechanics[END_REF]. It is recalled that for an energetic interface, the continuity of displacement across the interface holds.

Mori-Tanaka method for composites with energetic interfaces

Consider now an RVE that occupies the space with volume .

 0 V 1 0      1 0 V V V  
inhomogeneity has a volume fraction . According to the Mori-Tanaka theory, the

1 f v V V 
          0 0 0 0 0 ( ) 1 d 1 , 1 1 d d 1 f f s f f V v v V V S v v V V                        E E E T E E H E       
From these formulas one has (18)

      1 0 0 1 0 1 0 0 1 + , , 1 + f f f f v v v v                                   E T E A E E T A E A E H A E E    
where denotes the fourth-order identity tensor. The dilute concentration tensors are obtained  through the CCA method.

Analytical solutions for four special boundary-value problems

This approach is motivated by the studied boundary value problems by Hashin and Rosen [START_REF] Hashin | The Elastic Moduli of Fiber-Reinforced Materials[END_REF] in their famous composite cylinders assemblage theory. The necessary modifications on these problems, discussed in the sequel, provide the exact solution for the dilute concentration tensors.

Similar technique has been utilized in various articles for studying the mechanical and piezoelectric response of long fiber composites [START_REF] Benveniste | Stress fields in composites with coated inclusions[END_REF][START_REF] Wang | Annular Coated Inclusion model and applications for polymer nanocomposites-Part II: Cylindrical inclusions[END_REF][START_REF] Chatzigeorgiou | Micromechanical method for effective piezoelectric properties and electromechanical fields in multi-coated long fiber composites[END_REF][START_REF] Chatzigeorgiou | Elastic and inelastic local strain fields in composites with coated fibers or particles: theory and validation[END_REF]. Before passing to the actual boundary value problems, it is essential to express all the fields and the conservation laws in cylindrical coordinates.

Expressing Eshelby's inhomogeneity problem in cylindrical coordinates

Inside the representative volume element, the various fields generated at every phase ( ) 0,1 q  depend on the spatial position, i.e ( 19)

            , , , q q q q    U E  x x x x
Due to the geometry of the inhomogeneities, the problem can be transformed into cylindrical coordinates, using a system of concentric cylinders for the inhomogeneity and the infinite matrix, according to the equations in Appendix A. In cylindrical coordinates, the axes are transformed to , according to the relations

  1 2 3 , , x x x   , , r z  (20) 1 2 3 cos , sin , x r x r x z     
The vectors , and in cylindrical coordinates are expressed as All the fields are expressed in cylindrical coordinates

U E  (21)             2 , , 2 2 q q rr rr q zz zz r r r q q q rz rz z z z r r z z u u u                                                                                                 U E e e e
In the cylindrical coordinate system, the equilibrium and electrostatic equations are written as ( 22)

                                  1 0, 2 1
0, 1 0, 1 0 q q q q q r rr rr rz q q q q r r z q q q q z rz rz zz

q q q q r r z r r r z r r r z r r r z r r r z                                                                   d d d d
while the components at each phase are given by the expressions E (23)

                                      1 1 1 1 q r q q r q z q q q r q q q z r q q z q q q u r r u u r u z u r r u u r u r u z r u u z r r z                                                                        

E

The inhomogeneity is considered to have a constant radius . r a 

Due to the transverse isotropy of all phases, the generalized stiffness matrix retains the  same form in cylindrical coordinates. In the following subsections, the boundary value problems are presented and the analytical form of the solutions (in terms of the generalized vector

) is

  q U obtained.

Axial shear/in-plane electric field

Under this loading mode, displacement and electric potential boundary conditions at the far field ( ) read:

        1 2 , 1 1 2 , 1 0, cos , cos i i q q r q z q i i q e q i i u u r u r a r r a                               
with and . The unknown coefficients and are determined from the

1 1   2 1    , q i  , e q i
 boundary conditions, the interface relations Eqs. ( 9) and ( 10) and the condition that the displacements and the electric potential are finite at . 0 r 

Transverse shear strain

Under this loading mode, displacement boundary conditions at the far field ( ) read:

      1 4 , , 1 1 4
, 1 sin 2 , cos 2 ,

0 i i q r q i q i i q q i i q z r u r a r u r a u                              with ( 28
)               1 , 1 2 ,2 3 ,3 4 ,4 3, 2 1, 1, 3, 1 1
,

q q q q q q q q q q q k m k m k m m                      
The known coefficients are determined from the boundary conditions, the interface relations , q i  Eqs. ( 9) and [START_REF] Firooz | Homogenization accounting for size effects in particulate composites due to general interfaces[END_REF] and the condition that the displacements and the electric potential are finite at . 0 r 

Plane strain/axial electric field

Under this loading mode, displacement and electric potential boundary conditions at the far field (

) read:

        1 1 4 2 , , 1 1 1 4
, , 1 cos 2 , sin 2 , 0,

i i i q r q i q i i i q q i q i i q z q e r r u r r Z a a r u r a u z                                             
In the above expressions, and for are given by the relations described in the

, q i  i  1, 2,3, 4 i 
previous boundary value problem. Moreover, and . The unknown coefficients for

1 1   2 1
   the phase and are determined from the boundary conditions, the interface relations th q , q i  , q i Z Eqs. ( 9) and [START_REF] Firooz | Homogenization accounting for size effects in particulate composites due to general interfaces[END_REF] and the condition that the displacements and the electric potential are finite at . 0 r 

Hydrostatic strain

For this type of loading, displacement conditions at far field ( , ) correspond to ext r   ext z   hydrostatic strain:

(31) , 0,

ext ext ext ext ext r z u r u u z      
The analytical solution for the displacement field at every , and takes the form

r  z (32)       1 2 , 1 , 0, i q r q i i q q z r u r Z a u u z               
with and . The unknown coefficients are determined from the boundary

1 1   2 1    , q i
Z conditions, the interface relations Eqs. ( 9) and ( 10) and the condition that the displacements and the electric potential are finite at . 0 r 

Computing the dilute concentration tensors

Solving the previously described boundary value problems, the average fields can be

  q E
computed in both the extended fiber (including the interface) and the matrix. In these studies, the constants , and are arbitrary, thus one can choose them properly in order to construct the   e  dilute concentration tensors (for instance, setting one equal to 1 and zeros to the rest). The four discussed boundary value problems are sufficient for identifying and . The computations are T H included in Appendix B for the reader's convenience. Similar computational procedure is described in Chatzigeorgiou and Meraghni [START_REF] Chatzigeorgiou | Elastic and inelastic local strain fields in composites with coated fibers or particles: theory and validation[END_REF]; Chatzigeorgiou et al. [START_REF] Chatzigeorgiou | Micromechanical method for effective piezoelectric properties and electromechanical fields in multi-coated long fiber composites[END_REF]; Chen et al. [START_REF] Chen | Extended mean-field homogenization of viscoelastic-viscoplastic polymer composites undergoing hybrid progressive degradation induced by interface debonding and matrix ductile damage[END_REF] .

Numerical Results

In this section, the new Mori-Tanaka scheme considering energetic surfaces is utilized to obtain a full set of homogenized properties of nanoporous composites where the surface piezoelectricity effect prevails. In the following examples, the matrix phase is considered to be PZT-7A. The piezoelectric properties of the considered materials can be found in Chen et al. [START_REF] Chen | Parametric multiphysics finite-volume theory for periodic composites with thermo-electro-elastic phases[END_REF], which are also listed in Table 1. The surface parameters are typically identified from extensive molecular dynamics simulations of atomic lattices of finite thickness by progressively reducing the number of through-thickness atomic planes and extrapolating the limit to zero [START_REF] Miller | Size-dependent elastic properties of nanosized structural elements[END_REF]. Most publications in the literature used two sets of surface elastic constants, obtained by the molecular dynamics simulations for two different crystallographic orientations [START_REF] Miller | Size-dependent elastic properties of nanosized structural elements[END_REF], which yield stiffening and softening behavior of such materials at the nanoscale with decreasing pore sizes due to the positive and negative surface energy, respectively. Since a full set of surface parameters for PZT-7A are not available in the open literature, the surface parameters are taken as known quantities. As customary, to vigorously test the accuracy of the extended Mori-Tanaka approach, two types of surface parameters labeled as surface A and surface B, listed in Table 2, are considered, which yield negative and positive surface energies respectively. It should be noted that a justification of the applicability of the surface parameters with negative surface energy density remains an open question due to the lack of reliable experimental methods. Nevertheless, they serve as good candidates for critically assessing the predictive capability of the continuum mechanics-based simulation of elastic nanostructured materials with the extended Gurtin-Murdoch interfaces. It is worth mentioning that thermodynamically there is no restriction about the positiveness of the surface properties (contrarily to a bulk material, whose stiffness tensor needs to be positive definite). However, one needs to verify that the chosen surface parameters do not lead to an unstable overall response [START_REF] Javili | Relationships between the admissible range of surface material parameters and stability of linearly elastic bodies[END_REF][START_REF] Javili | Thermomechanics of Solids With Lower-Dimensional Energetics: On the Importance of Surface, Interface, and Curve Structures at the Nanoscale. A Unifying Review[END_REF][START_REF] Wang | The surface effects on solid and hollow nanowires under diametral loading[END_REF]. The extended theory is validated by comparison with the results using the CCA/GSCM technique and the finite-element approach developed in our recent work [START_REF] Chen | Homogenization of size-dependent multiphysics behavior of nanostructured piezoelectric composites with energetic surfaces[END_REF]. This validation is performed for nanoporous PZT-7A with different pore radii and pore volume fractions in large parameter space. It should be noted that, while real porous composites were modelled in Chen et al. [START_REF] Chen | Homogenization of size-dependent multiphysics behavior of nanostructured piezoelectric composites with energetic surfaces[END_REF] by excluding directly porosity from the bulk materials, we mimic porosities by assigning the fiber phase with very low moduli. The accuracy of the extended Mori-Tanaka approach to capture the size-dependent multiphysics behavior of nanoporous composites is further demonstrated by comparison with the finite-element results. It should be emphasized that a fundamental difference between the two approaches lies in the fact that while the former provides microstructural detail-free estimates for statistically homogeneous media, the latter explicitly takes into account realistic microstructures under the periodic material assumption. Chen et al. [START_REF] Chen | Homogenization of size-dependent multiphysics behavior of nanostructured piezoelectric composites with energetic surfaces[END_REF] were perhaps the first to investigate in a systematic manner the multiphysics moduli of this class of materials with incorporated surface piezoelectricity using an extended multiphysics finite-element asymptotic homogenization technique. The latter follows the classical homogenization framework described by Cavalcante et al. [START_REF] Cavalcante | Homogenization of elastic-plastic periodic materials by FVDAM and FEM approaches -An assessment[END_REF] and Chen and Wang [START_REF] Chen | Computationally-efficient homogenization and localization of unidirectional piezoelectric composites with partially cracked interface[END_REF]. Namely, the displacement field and electric potential are decomposed into homogenized and fluctuating components considering the scale separation and the unit cell problem is solved using the variational principle, subjected to the periodic boundary condition. In addition, the surface piezoelectricity effect is accounted for through additional surface energies associated with the elements separating different phases. The Mori-Tanaka method with surface piezoelectricity, however, is based on the strong form solution of the interface equilibrium and electrostatic conditions. Therefore, the solution techniques employed in the extended finite-element and Mori-Tanaka approaches with surface piezoelectricity are totally different, lending to the extended Mori-Tanaka method's rigorous validation.

The unit cell for the finite-element analysis is discretized into 48×12 quadratic quadrilateral elements as shown in Piezoelectric ceramics like PZT have very high piezoelectric properties but are difficult to bond onto curved surfaces due to their high stiffness and brittleness. This issue can be circumvented by developing piezoceramic-piezopolymer composites which possess flexibility from the polymer matrix and high piezoelectric properties from piezoelectric ceramics [START_REF] Jain | PVDF-PZT Composite Films for Transducer Applications[END_REF].

Therefore, in the last case of the example, we focus on a nanocomposite with unidirectional PZT fiber embedded in a less-active PVDF polymer matrix and study the surface piezoelectricity effect on the effective properties of this material system. The material properties of the PVDF found from [START_REF] Chen | Parametric multiphysics finite-volume theory for periodic composites with thermo-electro-elastic phases[END_REF] were listed in Table 3. The PZT fiber has a fixed volume fraction of 50% but with varying sizes. Figure 11 presents the comparison of the two mechanical constants, , , and  of surface piezoelectricity, as a function of the fiber radius, generated by the extended Mori-Tanaka with incorporated Gurtin-Murdoch interface. It is observed that for this materials system, the effective moduli associated with the fiber axial direction ( and ) exhibit remarkable size-n  * 33  dependent effects for both types of surfaces. However, due to the presence of the soft PVDF matrix, the in-plane moduli ( and ) are less sensitive to the PZT fiber size. k  * 11  Table 3 Material parameters for the PVDF matrix An extended Mori-Tanaka homogenization approach has been developed for the first time for nanocomposites with transversely isotropic piezoelectric phases by incorporating the surface piezoelectricity effect based on the generalized Gurtin-Murdoch coherent interface model. This extension enables efficient parametric studies of the not-well-characterized unidirectional piezoelectric composites in the nanoparameter range where surface piezoelectric effects have the potential to substantially alter the homogenized moduli under electromechanical loadings.

n  k  two dielectric constants
  GPa k   GPa l   GPa n   GPa p   GPa m
The main novelty of this work lies in the analytical nature of computing the dilute concentration tenors for both nanofibers and energetic surfaces, derived from the exact solution of Eshleby's inhomogeneity problems under specific macroscopic electromechanical loading conditions. The significant results generated in this work demonstrate that the extended theory is capable of reproducing the size effect of piezoelectric nanocomposites when the characteristic size of the inhomogeneity is sufficiently small. Except for the transverse shear loading, the new meanfield framework yields indistinguishable homogenized moduli relative to the CCA and FE reference results even at high volume fractions and small pore sizes for unidirectional piezoelectric nanocomposites under axisymmetric, axial shear, and electric field loading. The proposed analytical approach not only provides a useful tool for rapid identification of the size-dependent homogenized moduli of nanocomposites under electromechanical loading, but also paves the way for the further extension of the model to account for the inelastic mechanisms through the transformation field analysis technique against which other micromechanics approaches may be compared.

Equation Section 1 Appendix A: Transformation between coordinate systems and useful integrals

The base unit vectors in cylindrical coordinates are written as (A34)

cos sin 0 sin , cos , 0 0 0 1 r z                                        n n n
Transformation between Cartesian and cylindrical coordinate systems for the fields is obtained according to the formulas [START_REF] Chatzigeorgiou | 1 -Mathematical Concepts[END_REF][START_REF] Chatzigeorgiou | Chapter 1 -Tensors[END_REF] (A35)

cyl Car cyl Car E cyl Car        U R U E W E W   where 3 1 6 3 6 3 1 3 3 6 3 6 , , 1 v E v v                               R 0 Q 0 Q 0 R W W 0 0 R 0 R and , cos sin 0 sin cos 0 0 0 1 v                 R , 2 2
2 2 1 cos sin 0 sin 2 0 0 2 1 sin cos 0 sin 2 0 0 2 0 0 1 0 0 0 sin 2 sin 2 0 cos 2 0 0 0 0 0 0 cos sin 0 0 0 0 sin cos                                           Q 2 2 2 2 cos sin 0 sin 2 0 0 sin cos 0 sin 2 0 0 0 0 1 0 0 0 1 1 sin 2 sin 2 0 cos 2 0 0 2 2 0 0 0 0 cos sin 0 0 0 0 sin cos                                        

Q

The average in the inhomogeneities/interfaces system is expressed in cylindrical

  1 E coordinates as (A36)     2 2 1 T T 2 2 0 0 0 1 1 d d d d 2 2 L a r z z L z L L r a a z r r z L a L a                       E U U U   with , ,   1 T r z u u u                  U R T cos
0 0 0 0 r                                       T 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 z                             



For the interface, we introduce the extended stress ( )

s  (A37) T ( ) ( ) ( ) ( ) ( ) ( )
0, , , 0, 0, , 0, ,

s s s s s s zz z z             d d
Thus the average in the inhomogeneities/interfaces system is expressed in cylindrical

  1  coordinates as (A38)   ( ) 2 2 1 T T 2 2 0 0 0 1 1 d d d d d 2 2 s L a L E E L L r r z a r z L a L a                 W W   
Equation Section (Next)

Appendix B: Analytical computations for the four boundary value problems

The unknown constants , and are identified using i) the boundary conditions at , q i  , e q i  , q i Z , ii) the consistency condition that the fields should be finite at , and iii) the interface ext r r  0 r  conditions Eqs. ( 9) and [START_REF] Firooz | Homogenization accounting for size effects in particulate composites due to general interfaces[END_REF] between the material phases.

B.1. Axial shear/in-plane electric field

For this boundary value problem, the boundary conditions and the fact that all fields should be finite at yield

0 r  (B40) 1,2 1,2 0,1 0,1 0, , e e e          
The rest of the unknown constants are given by the solution of the system (B41)

e e xz xz xz      K Ξ F F where , , , 1,1 0,2 1,1 0,2 e e
                 Ξ     0 0 15 1 0 xz p e                F     0 15 0 11 1 0 e xz e                  F                         1 0 1 0 15 15 66
                                 

K

The solution of this system can be written in the form

(B42) mec e elc     Ξ Ξ Ξ
The average and in the fiber are given by E  (B43)

    1 1 0 0 0 0 0 0 0 0 0 0 0 0 1 , 0 0 0 0 0 0 0 0 0 s xz xz xz e e se xz xz xz M M a M M                                                                                          E  with                 1,1 1,1 1 1 1,1 15 1,1 1 1 15 1,1 11 1,1 66 1,1 26 1,1 26 1,1 22 1,1 , , , , , xz

B.2. Transverse shear

For this boundary value problem, the boundary conditions and the fact that all fields should be finite at yield

0 r  (B45) 1,3 1,4 0,1 0,2 0,         
The rest of the unknown constants are given by the solution of the system (B46)

xy xy    K Ξ F with     1,1 0 1,2 0,3 0 0,4 1 2 , 1 2 xy m m                                 Ξ F                 1,1 0,4 1 0 0 22 1 0 0 41 22 1 1 0 0 0 0 0 2 6 4 0 0 1 1 1 1 1 0 0 0 0 2 6 2 2 2 0 0 s s xy s s m m k L a k m m k L                                           K with         1 1 ( ) 41 1,1 22 1 1 6 , 2 2 s s k m k L k m          
The solution of this system can be written in the form

(B47) mec   Ξ Ξ
The average and in the fiber are given by E  (B48)

    1 1 0 0 0 0 0 0 0 0 0 2 2 , 2 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 s xy xy xy M M a                                                                                         E  with     1,1 1,1 1,2 1 22 1,1 1,1 1,2
M m L M                       
All the above constants can be written in the form

(B50) 1,3 1,4 0,1 1,2 0,2 0,1 0, , Z Z            
The rest of the unknown constants are given by the solution of two systems:

1. The first is written as

(B51) e e xx xx xx      K Z F F with       1,1 ( ) 1 0 0 0,2 32 31 31 0 1 0 1 , , 2 e xx xx s Z Z e a e e k                              Z F F       1 0 22 0 0 1 1 1 0 2 2 xx s a L k m                K
The solution is written in the form (B52) The rest of the unknown constants are given by the solution of the system:

B54)                   (
                                                                                     E    33 0 0 0 0 2 s se e xx M                                  with (B55)             1,1 1 1,1 1 1,1 1 31 1,1 22 1,1 23 1,1 32 1,1 , 2 , 2 , 2 , , 2 , 2
(B59) xx h    K Z F with         1,1 0 0 1 0,2 23 0 1 1 , 2 h s Z Z a L k l l                         Z F
and is given by Eq. (B51).

                                                                                        E  with                       1,1 1 1 1,1 1 1 1,1 1 31 1,1 22

T T T T T T T T T T T T T T T T T

                             T (B64)
                           
where the elements of each matrix are obtained by solving the boundary volume problems detailed in Appendices B.1-B.4.

Figure 1 A

 1 Figure 1 A piezoelectric solid with randomly-dispersed nano-fibers/porosities with a constant radius . The Gurtin-a

33 

 33 constants. In a compact form, Eq. (3) can be written as:

Figure 2

 2 Figure 2 Cross section of cylindrical fiber with the energetic interface (highlighted in blue) inside an infinite medium.

  macroscopic and tensors are related to the corresponding fields in the inhomogeneity and E  the interface according to the formulas (17)

Figure 3

 3 Figure 3 Cross section of cylindrical fiber with the energetic interface (highlighted in blue) inside an infinite medium.

  arbitrary known constant. The analytical solution for the displacement field at every 

Figure 4

 4 Figure 4 illustrates the differences of the axisymmetric moduli , and the axial shear * k l 

Figure 4 Figure 5 Figure 6

 456 Figure 4 Normalized Hill's axisymmetric and shear moduli of a transversely isotropic porous PZT nanocomposite with 0.30 porosity volume fraction as a function of the porosity radius. Comparison with the CCA and three-phase model predictions of Chen et al.[START_REF] Chen | Homogenization of size-dependent multiphysics behavior of nanostructured piezoelectric composites with energetic surfaces[END_REF]

Figure 7 .

 7 The hexagonal repeating unit cell retains the most suitable microstructure for direct comparison with the extended Mori-Tanaka approach as it yields transversely isotropic material properties. The plots of homogenized moduli generated by the finite-element and extended Mori-Tanaka approaches as a function of porosity volume fraction at two different radii a=5nm and a=20nm are presented in Figures 8-10. In the case of transverse shear moduli, the GSCM predictions are also enclosed in the figure for comparison. It is observed the extended Mori-Tanaka and finite-element predictions show remarkable accordance for the axisymmetric and axial shear moduli, piezoelectric and dielectric constants, even at very highvolume fractions. In comparison with the finite element results, the discrepancies of the transverse shear moduli generated by the extended Mori-Tanaka are smaller relative to the GSCM technique.

Figure 7 Figure 8 Figure 9 Figure 10

 78910 Figure 7 Repeating unit cell discretization using the quadratic quadrilateral element for the finite-element simulation

Figure 11

 11 Figure 11 Normalized effective constants of a transversely isotropic PZT/PVDF nanocomposite with 0.50 porosity volume fraction as a function of the fiber radius for two types of surfaces.

3 .

 3 Plane strain/axial electric fieldFor this boundary value problem, the boundary conditions and the fact that all fields should be finite at yield 0 r 

  as those of the expressions in Eq. (B48).For this boundary value problem, the boundary conditions and the fact that all fields should be

Table 1

 1

				Material parameters for the bulk PZT-7A ceramics
	k	  GPa	l	  GPa	n	  GPa	p	  GPa	m	  GPa
		112.1		74.24		131		25.30		35.9

Data Availablity

The raw/processed data required to reproduce these findings cannot be shared at this time as the data also forms part of an ongoing study

Declaration of interests

☒ The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

☐ The authors declare the following financial interests/personal relationships which may be considered as potential competing interests:


Adrien Girard

Lucien Saviot

Silvia Pedetti

Mickaël D Tessier

Jérémie Margueritat

Hélène Gehan

Benoit Mahler

Benoit Dubertret

Alain Mermet
email: alain.mermet@univ-lyon1.fr

Mass load effect on the resonant acoustic frequencies of colloidal semiconductor nanoplatelets

Resonant acoustic modes of ultrathin

CdS and CdSe colloidal nanoplatelets (NPLs) with varying thicknesses were probed using low frequency Raman scattering. The spectra are dominated by an intense band ascribed to the thickness breathing mode of the 2D nanostructures. The measured Raman frequencies show strong deviations with respect to values expected for simple bare plates, all the more so as the thickness is reduced. The deviation is shown to arise from the additional mass of the organic ligands that are bound to the free surfaces of the nanoplatelets. Calculating the eigenfrequencies of vibrating platelets weighted down by the mass of the organic ligands provides very good agreement with the observed experimental behaviours. This finding opens up a new possibility of nanomechanical sensing such as nanobalances.

Resonant acoustic modes are typical manifestations of size effects that occur as the wavelength of the vibration approaches the size of the system that sustains it. Macroscopic examples of resonant acoustic modes are 1D standing waves of a string, 2D flexural modes of plates or 3D breathing-like oscillations of spheres. Nowadays, the nanoscale equivalent of these types of modes raises strong interest particularly for nanoelectromechanical system (NEMS) applications. The most common way of exploiting nanomechanical modes for sensing purposes is to resonantly excite them through electrical actuators and track resonance frequency changes due to external perturbations (mass loading, pressure, deformation...). The typical frequencies of the so-solicited mechanical resonances lie in the sub-GHz or GHz

range and the mode shapes can cover a wide variety of displacement fields (bending motion of a cantilever or of a long nanotube 1 , flexural modes of nanosheets 2 ,...). A key challenge in NEMS applications is to develop smaller and lighter resonators operating at higher frequencies with good quality factors. Since resonance frequencies scale with inverse sizes, nano-objects with very small dimensions (of a few nanometers) potentially fulfill these requirements; yet, detecting their very high resonance frequencies approaching the THz range require different experimental setups than the classical electromechanical devices.

Inelastic Raman light scattering at low frequencies has proved as an efficient tool to probe resonant acoustic mode frequencies of nano-objects with very small sizes. Initially developed from quasi-spherical nanocrystals nucleated in glasses 3 , this technique has been successfully employed to characterize resonant acoustic modes (also called Lamb modes) from a variety of spherical nanostructures like quantum dots 4,5 , metallic nanospheres [6][7][8] or dielectric nanoparticles [9][10][11] .

Compared to spherical or quasi-spherical nanostructures, investigations of resonant acoustic modes from 2D nano-objects with nanometer thicknesses are more scarce. Recent low frequency Raman works on a few layer graphene sheets 12 identified interlayer shear modes below 50 cm -1 (1.5 THz) thus enabling to probe the interlayer interactions. Similar measurements on nanosheets of transition metal dichalcogenides 13 also led to the detection of interlayer shear modes together with a new type of mode ascribed to the breathing thickness vibration of the nanosheets. This latter mode is the 2D equivalent of the spherically symetrical breathing vibration of a nanosphere (also referred to as the � = 0 Lamb mode) that is routinely detected in low frequency Raman spectra from spherical nanoparticles [START_REF] Saviot | Nanoparticles and Quantum Dots[END_REF] . Like for spherical nanoparticles, the experimental Raman breathing mode frequencies of the 2D nanosheets were found to well agree with the values derived from the continuum elasticity approach whereby the nanosheet is modeled as a freely vibrating plate.

The present paper is dedicated to the low frequency Raman

1-5 | 1
characterization of the thickness breathing mode from atomically flat colloidal semi-conductor nanoplatelets. At variance with comparable measurements on the exfoliated nanosheet systems, strong deviations from freely vibrating plate frequencies are observed. We show that these deviations are caused by the additional mass load of the organic ligands lying at the surface of the nanoplatelets, differently from known acoustic impedance effects inherent of a continuous surrounding medium.

CdS and CdSe NPLs were synthesized through colloidal chemistry routes in the zinc-blende (ZB) structure [START_REF] Ithurria | [END_REF]16 . NPLs with thickness ranging between 3 and 5 monolayers (MLs) were obtained using hot injection reactions previously described 17 . Sulfur and selenium precursors, respectively, were injected at high temperatures (between 180 and 250 • C) into a reaction mixture of cadmium carboxylate complexes in 1-octadecene. CdS and CdSe NPLs possessing thicknesses larger than 5 MLs were prepared through colloidal atomic layer deposition reactions (c-ALD 16 ). This technique allows the deposition of several CdS and CdSe monolayers on thinner CdS and CdSe NPLs, respectively. The reaction is performed at room temperature using inorganic ions precursors dispersed in organic polar solvent such as N-methylformamide (NMF) 16 . The monolayers growth occurs by the alternate reaction of anionic and then cationic precursors at the large top and bottom lateral surfaces of thinner NPLs cores. For example, CdSe NPLs with 11 MLs were obtained growing four additional CdSe layers on CdSe NPLs with 3 MLs.

While the lateral dimensions of the NPLs (of the order of tens to a hundred nanometers) are difficult to control, their thicknesses can be precisely adjusted between 3 and 14 ML, i.e. between approximately 1 nm and 4 nm respectively 18 . The large bottom and top faces of the nanoplatelets are (100) faces bearing Cd + cations that are bound to oleate (OA) ligands. Detailed morphological characterizations and thickness evaluations through transmission electron microscopy and high-angle annular dark field TEM are given elsewhere 16,18 .

UV-visible absorption spectra of the NPLs dissolved in their native solvents (hexane or toluene) are displayed in Figure 1(A),(C). These spectra show remarkably narrow absorption bands associated with the strong confinement of the exciton across the NPL thicknesses and with a negligible dependence on their lateral dimensions 19 . The sharpness of the excitonic transitions, relative to the heavy hole-and light hole-to electron transitions (hh-e and lh-e), reflects the high thickness monodispersity of the NPL assemblies which is a definite asset for the low frequency Raman characterizations. While the heavy hole absorption edge in CdSe NPLs shifts from 450 nm up to 650 nm as thicknesses increase from 3 ML up to 11 ML respectively, in CdS it only shifts from 405 nm up to 500 nm over the 4 -14 ML thickness range.

The Raman measurements were performed from evaporated drops of solutions on copper substrates at room temperature. The Raman frequencies were determined from measurements with different spectrometers operating either at fixed (532 nm and 633 nm) or adjustable (660 nm, 560 nm, 514.5 nm and 496.5 nm) excitation wavelengths in order to work in nearly resonant conditions with the hh exciton (see Supplementary Information). citation wavelengths did not reveal any significant resonance dependence, i.e. frequency values agreed within ±5% whatever the excitation wavelength. Figure 1(B) and (D) show the low frequency Raman spectra of the CdS and CdSe NPLs. For both systems, one observes a well defined Raman peak shifting from about 10 cm -1 up to about 35 cm -1 as the thickness of the NPLs decreases. This low frequency Raman signal can be safely associated with the NPLs as the larger range Raman spectra (insets of Figure 1(c),(d)) attest the detection of the longitudinal optical modes near 300 cm -1 and 200 cm -1 for CdS and CdSe NPLs respectively. The detailed study of the dependence of the LO modes upon the NPL thicknesses was reported previously 20,21 ; it lies outside the scope of this paper.

Figure 2 reports the dependence of the low frequency mode versus the inverse thickness (1/h) of the nanoplatelets. At first order, the dependence is observed to be linear. Along similar inelastic light scattering experiments from thin films 22,23 , the most natural explanation for this behavior is a 'thickness mode' which consists of a longitudinal standing wave across the thickness of the platelet. In the case of a freely vibrating plate the fundamental frequency of such mode depends on the elastic constant C 11 and on the density ρ of the plate according to

ν = 1 2h � C 11 ρ (1)
where � C 11 ρ = V L is the longitudinal sound velocity along the (100) thickness direction of the NPL material. This expression is typically used to evaluate the resonance frequencies of piezoelectric microplates [START_REF] Fujishima | CRC handbook of electrical filters[END_REF] (notwithstanding the change of C 11 for C 33 whenever the piezoelectric material structure is non-cubic). Note that the inverse thickness dependence of the low frequency Raman mode implies that its spectral width decreases with increasing thickness (considering a nearly identical thickness distribution), as is observed in Figure 1(c) and(d). Comparing the experimental Raman frequencies (symbols in Figure 2) with those expected from Equation 1 using the macroscopic values of C ZB 11 (larger slope lines in Figure 2) reveals that the resonance frequencies strongly deviate from the nominal values, all the more so as the thickness of the NPLs is reduced. This behaviour strikingly constrasts with that observed from exfoliated MoS 2 and WSe 2 nanosheets 13 (bottom graph of Figure 2) for which the Raman frequencies conform with those expected from the corresponding macroscopic values C 33 and ρ downto thicknesses of 1.3 nm (see Supplementary Information).

The observed deviations with respect to the nominal frequencies derived from the bulk sound velocities � C 11 ρ show two important characteristics: the reduction is not constant with thickness and it is quite significant (it reaches about -50% for the thinnest platelets). This reduction cannot be explained by a possible relaxation of lattice parameters which was shown to amount to a few percents 26 . It can neither be explained by different boundary conditions of the thickness mode, namely the case of a fully supported nanoplatelet with the resting face fixed and the opposite one free. Indeed, in these latter conditions, the prefactor 1/2e in Equation 1 should be replaced by 1/4e, leading to a constant frequency reduction of -50% whatever the NPL thickness. Finally, assuming the NPLs as fully embedded in an organic matrix of acoustic impedance Z m (Z m = ρ m V m L ) would result essentially in a damping, i.e. a broadening of the Raman peaks without significant frequency shifts, as observed with embedded CdS spherical nanoparticles 27 .

The fact that the observed resonance frequencies are all the more downshifted as the thickness of the NPLs is reduced points to a typical lumped mass effect of a resonator under an additional mass load. In the present case of colloidal NPLs, the mass of the organic ligands covering the free surfaces of the platelets, which is not taken into account in Equation 1, is expected to play the role of an additional mass.

In order to verify this assumption, we can model a ligand capped NPL as a freely vibrating * slab of thickness h centered at z = 0 with two identical surface layers of surface mass densities σ at z = ± h 2 . Assuming a purely longitudinal motion, the equation of motion for an element of mass δ M = σ δ S writes as:

δ SC 11 ∂ u ∂ z � z = h 2 ,t � = -δ M ∂ 2 u ∂t 2 � z = h 2 ,t � (2)
where the left-hand expression is the pressure force acting on δ M. Searching for symmetric standing wave solutions (because the anti-symmetric ones are not Raman-active):

u(z,t) = A sin � ω V L z � [B cos(ωt) +C sin(ωt)] (3)
and applying the previous boundary condition leads to the transcendental equation:

C 11 ω V L cos � ω V L h 2 � = σ ω 2 sin � ω V L h 2 � (4)
which, since ρV 2 L = C 11 , can be written as:

cos � ω V L h 2 � = σ ω ρV L sin � ω V L h 2 � ( 5
)
where ω = 2πν. The fundamental thickness breathing mode frequency is given by the smallest solution of this equation, which can only be obtained numerically on the basis of the input parameters C ZB 11 , ρ ZB and σ OA . σ OA can be estimated from the number density of OA ligands. Indeed, charge neutrality of the NPLs imposes that each surface Cd + cation binds to one oleic acid (M OA = 282.46 g.mol -1 ) so that according to the nominal surface density of 2 Cd + per a 2 in (100) planes of the ZB structure † , σ OA = 2 a 2 M OA N A (N A Avogadro number), i.e. 5.9 OA/nm 2 and 5.4 OA/nm 2 for CdS and CdSe respectively. Figure 3 compares the results of the calculations with the experimental frequencies. Considering the absence of any fit parameter, the agreement is remarkable both qualitatively and quantitatively thereby confirming the assumed mass effect. The agreement between the experimental data points and the calculations may be improved by adjusting the input parameters ρ ZB , C ZB 11 and σ OA . Since there is no analytical expression for the lumped mass frequencies, we have tested the separate influences of the three parameters for the CdS NPLs (see Supplementary Info). It comes out that a better agreement can be obtained by either increasing ρ ZB or σ OA by 30% or reducing C ZB 11 by 10% or by varying simultaneously the three parameters in lower ratios. A sole effect of a 30% larger mass density seems however unlikely as the relaxation of the lattice parameters in NPLs is not expected to exceed a few percents 26 . Conversely, a 10% reduction of the elastic constant C ZB 11 with respect to the assumed value of 88 GPa (estimated from ab-initio calculations [START_REF] Deligoz | [END_REF] ) appears more plausible considering the expected effects of size reduction on the elastic properties of nanoplates [30][31][32] . The fact that both CdS and CdSe NPLs have the thanks to slightly interdigitized ligands, the most natural boundary conditions to apply for an average NPL is free motions at both ends † where a is the lattice parameter of the ZB structures, i.e. a = 0.582 nm for CdS and a = 0.608 nm for CdSe tendency to curl at very small thicknesses 29,33 is indeed an indication of perturbed elasticity. Finally, the range of acceptable adjustability for σ OA remains to be better characterized. Comparing the assumed OA number density of 5.4 OA/nm 2 for CdSe NPLs with that of 4.6 OA/nm 2 measured through NMR spectrometry from CdSe quantum dots capped with OA ligands 34 points to an already large value. Yet, the smaller value observed in QDs can be explained by the fact that the different facets of QDs are not all (100) planes resulting in a lower average Cd + density compared to NPLs. It therefore appears inappropriate to increase σ OA to obtain a better agreement between the experimental and the calculated frequencies. Note that the purification steps during the synthesis of the NPLs may lead to localized removals of the ligands or to extra mass loads from adsorbed reaction by-products. This might explain the relative data dispersion for CdSe NPLs; more systematic studies using fresh samples and different ligands can help clarify this situation.

From the above considerations, it occurs that the colloidal shell of a nano-object can significantly affect its vibration resonance frequency. Such finding contrasts with the widely studied cases of ligand capped metallic nanospheres (e.g. thiol capped Au 8 or Ag 6 nanoparticles) for which Lamb mode frequencies evaluated on the basis of freely vibrating bare nanospheres were found to agree well with the Raman experimental values. Considering that the surface-to-volume ratios of a nanosphere of radius r and of a nanoplatelet of thickness h are comparable (3/r vs 2/h respectively), what appears as a decisive factor in the observation of the mass effect is the ratio of the soft shell surface mass density σ lig to the hard core mass ρ core L (where L is the confining dimension of the nano-object, i.e. h or r): the lower σ lig /ρ core L, the weaker the effect on the resonance frequency. As an example, this ratio decreases by about a factor 5 going from an oleic acid capped CdS NPL of 1 nm thickness to a dodecanethiol capped Au nanosphere of 3 nm diameter (assuming equal ligand surface number densities for NPLs and nanospheres). Our finding suggests that previous interpretations of significantly downshifted Lamb mode frequencies from very small (r < 2 nm) CdSe nanoparticles capped with comparable organic ligands 35 may have omitted this effect.

The evidencing of a lumped mass effect in the resonance frequencies of NPLs brings new insight into the effect of a surrounding on the vibration of nano-objects. So far, the effect of a surrounding on the Lamb modes of nanoparticles has been essentially described by acoustic mismatch effects that have proved well suited for nanoparticles embedded in a matrix 36,37 . The impact of an embedding medium on the resonance frequencies of nanospheres can be estimated thanks to the Complex Frequency Method 27 ; it essentially consists in a broadening of the Raman peaks. The presently evidenced mass effect differs from this approach in the sense that the close environment of the nano-object behaves more as an additional inert load rather than as a continuous medium with an acoustic impedance. This implies that at the probed frequencies the dynamics of the ligands is closer to single particle dynamics rather than to phonon-like collective dynamics. This novel ingredient of core-shell nano-object vibrations deserves better characterization.

As a conclusion, the thickness breathing modes of CdS and CdSe colloidal nanoplatelets with thicknesses varying between 1 and 4 nm were probed thanks to low frequency Raman scattering. Unlike similarly investigated systems (nanospheres or exfoliated nanosheets), the resonance frequencies were found to be strongly downshifted with respect to the expected frequencies of freely vibrating plates. Taking into account a lumped mass effect due to the presence of organic ligands at the surface of the nanoplatelets proves to well describe the deviation. The strong reduction (reaching nearly 50%) of the resonance frequencies observed for the thinnest nanoplatelets allows to envisage their use as molecular mass sensors. Compared to other nano-objects like binanopyramids 38 , the nanobalance effect using nanoplatelets is expected to be much more sensitive thanks to their extended lateral faces combined very small thicknesses.

These works were supported by the ANR Nanovip project, grant ANR-13-JS10-0002 of the French Agence Nationale de la Recherche.

Mass load effect on the the resonant acoustic frequencies of colloidal semiconductor nanoplatelets: Supplementary Information

Adrien Girard, a Lucien Saviot, b Silvia Pedetti, c Mickaël D. Tessier, c Jérémie Margueritat, a Hélène Gehan, a Benoit Mahler, a Benoit Dubertret, c and Alain Mermet * a

Raman measurements from CdS and CdSe NPLs

In order to enhance the low frequency Raman signal, nearly resonant excitation conditions were chosen with Raman excitations lying slightly above the hh exciton absorption edge. This choice is a compromise to guarantee an enhanced interaction with light through strong absorption with a somewhat minimum photoluminescence which tends to dominate the Raman spectra and thus hide the weaker low frequency Raman signals. Because of the operation constraints of the Raman spectrometers, several setups were used.

A multichannel micro-Raman HORIBA Labram HR spectrometer equipped with Notch filters operating at the fixed excitation wavelength of 532 nm was used for all CdS samples. Differently, the strong shift of the NPL excitonic edge for CdSe samples required to operate with a wavelength versatile macro setup built of a 5-grating monochromator equipped with single channel detection. The Raman spectra of the 3 ML, 4 ML, 5 ML and 11 ML CdSe NPLs were recorded with the 496.5 nm, 514.5 nm, 560 nm and 660 nm excitation lines respectively. The 6 ML CdSe NPL sample was investigated with a Renishaw Invia micro-spectrometer operating at 633 nm.

For several samples, it was checked that the Raman frequencies did not depend significantly upon wavelength excitation (Figure 1). For CdSe NPLs, the frequencies were deduced from the average values between 532 nm 633 nm. ) CdSe @ 532 nm CdSe @ 633 nm Since bulk elastic constants of zinc-blende phase CdS and CdSe have, to our knowledge, not been characterized (available elastic data pertain essentially to the macroscopically stable wurtzite phase), we experimentally checked the value of C ZB 11 = 98 GPa which was derived from ab-initio calculations 1 for bulk ZB-CdS.

To this end, bare (i.e. ligand free) CdS quantum dots in the ZB phase were synthesized in order to measure, through low frequency Raman scattering, one known Lamb mode whose frequency depend on the value of C ZB 11 .

Synthesis of the bare ZB-CdS

Raman measurements

The Raman measurements were performed in identical conditions as those described in the paper. Measurements were performed at 532 nm (using the HORIBA Labram HR micro-spectrometer) and at 633 nm (using the Renishaw Invia micro-spectrometer). Figure 3 displays the corresponding low frequency Raman spectra; once again, no dependence upon excitation wavelength is observed.

In accord with comparable measurements performed on CdS QDs embedded in glasses 3 , the prominent peaks are ascribed to the fundamental breathing Lamb modes (described by the angular momentum � = 0) of the quasi-spherical QDs. For bare spherical nanoparticles, the Lamb � = 0 mode frequency is given by:

ν �=0 � 0.9 V L D = 0.9 D � C 11 ρ (1)
A linear fit of the experimental frequencies (Figure 3-bottom) along this equation yields C ZB-CdS

11

= 98 GPa (ρ ZB-CdS = 4870 kg.m -3 ). This value is in excellent agreement with that derived from ab-initio calculations 1 . 3 Analysis of published Raman data from MoS 2 and WSe 2 nanosheets 4 Zhao et al 4 recently reported low frequency Raman measurements from MoS 2 and WSe 2 nanosheets with thicknesses lying between 1 and 6 nm. Compared to our colloidal NPLs, these systems have no ligands. Among the detected low frequency modes, two of them (named B1 and B2 in Ref. 4 ) show frequencies decreasing with increasing thickness. DFT calculations 4 helped assign these modes to the fundamental and third overtone of so-called interlayer breathing modes. These modes are in fact thickness breathing modes, along the descriptions in Figure 4 of Ref. 4 . Accordingly, plotting the fundamental B1 frequencies of MoS 2 and WSe 2 nanosheets as a function of the corresponding reciprocal nanosheet thicknesses * (bottom graph of Figure 2 in the paper), one observes linear trends as expected from Equation 1 of the paper. Linear fits of these data yield, using macroscopic densities † , C 33 values of 55±2 GPa and 51±1 GPa for MoS 2 and WSe 2 respectively (note that differently from cubic structures, the speed of sound in such hexagonal close packed structures along the thickness direction is given by � C 33 ρ ). These values perfectly match the C 33 values derived from the mass-spring model developed in Ref. 4 .

4 Evaluation of the influences of the parameters C 11 , ρ and σ on the NPL frequencies mented in the main paper.

Fig. 2

 2 Fig. 2 Top graph: Evolution of the low frequency Raman shift with the inverse nanoplatelet thickness for CdS and CdSe NPLs (symbols). If not visible, error bars are within symbol sizes. Low slope lines are best linear fits to the data. Larger slope lines are calculated expected frequencies with corresponding bulk values of C ZB 11 and ρ ZB of CdS and CdSe. Bottom graph: Experimental Raman shifts (symbols) and calculated expected frequencies with bulk C 33 values for MoS 2 and WSe 2 exfoliated nanosheets from Ref. 13 .

Fig. 3

 3 Fig. 3 Calculated thickness breathing mode frequencies (solid lines) of CdS and CdSe NPLs derived from a mass sandwich nanoplate model assuming a stoechiometric Cd + /OA -1:1 ratio for the ligand surface mass densities σ OA . Symbols are experimental frequencies and dashed lines are expected frequencies for bare NPLs.

a

  Institut Lumière Matière, Université de Lyon, Université Claude Bernard Lyon 1, UMR CNRS 5306, 69622 Villeurbanne, France b Laboratoire Interdisciplinaire Carnot de Bourgogne, UMR 6303 CNRS-Université de Bourgogne Franche Comté, 9 Av. A. Savary, BP 47 870, F-21078 Dijon Cedex, France c Laboratoire de Physique et d'Etudes des Matériaux, CNRS, Université Pierre et Marie Curie, ESPCI, 10 rue Vauquelin, 75005 Paris, France * alain.mermet@univ-lyon1.fr

Fig. 1

 1 Fig. 1 Thickness mode frequency vs 1/h measured with 532 nm (full circles) and 633 nm (empty circles) excitation lines for selected CdS (top) and CdSe (bottom) NPL samples. No significant dependence upon wavelength excitation is observed. Lines are expected values from bare NPLs.

Fig. 2

 2 Fig. 2 Optical absorption spectra of bare ZB-CdS quantum dots with different diameters.

1 )Fig. 3

 13 Fig. 3 Top : Low frequency Raman spectra of bare CdS quantum dots recorded with two different excitation wavelengths, 532 nm (green lines) and 633 nm (red lines). Indicated sizes are diameters. Bottom : Frequency dependence of reciprocal diameter 1/D. The line is a zero intercept fit along Equation 1, yielding C CdS 11 = 98 ± 2 GPa.

Figure 4 11 =Fig. 4

 4114 Figure 4 displays the separate dependencies of the lumped mass frequencies for CdS NPLs on the input parameters C 11 , ρ and σ , i.e. changing one of the three parameters while fixing the two other ones at their nominal values (C ZB-CdS 11 = 98 GPa, ρ ZB-CdS = 4870 kg.m -3 and σ CdS OA = 5.9 OA/nm 2 ). The results are com-

  Double checking the Raman frequency values with different ex-

			A -Optical absorption CdS		B -Raman CdS			CdS 12ML LO
	Absorbance (a.u)		8ML 14ML 12ML 10ML			Intensity (a.u)	8ML 10ML 12ML 14ML			532 nm	0	100	200 cm -1	300
			4ML				4ML
	400	500	600	700	800	-60	-40	-20	0	20	40	60
		Wavelength (nm)						Raman shift (cm -1 )
			C -Optical absorption CdSe		D -Raman CdSe			CdSe 4ML LO
	Absorbance (a.u)			6ML 5ML 4ML	11ML	Intensity (a.u.)	6ML 11ML 5ML 4ML			660 nm 633 nm 560 nm 514.5 nm	0	100	cm -1	200
			3ML				3ML			496.5 nm
	400	500	600	700	800	-60	-40	-20	0	20	40	60
		Wavelength (nm)						Raman shift (cm -1 )
	Fig. 1 Top graphs: Optical absorption spectra of CdS (A) and CdSe (C)
	NPLs nanoplatelets with varying thicknesses (ML=monolayer). Dashed
	vertical lines indicate the Raman excitation wavelengths. Bottom
	graphs: corresponding low frequency Raman spectra of CdS (B) and
	CdSe (D) NPLs. For CdSe NPLs the different Raman excitation lines
	are indicated.

Nanoscale