Effective anisotropy of periodic acoustic and elastic composites
Vincent Laude, Julio Andres Iglesias Martinez, Yan-Feng Wang, Muamer Kadic
EEffective anisotropy of periodic acoustic and elastic composites
Vincent Laude, a) Julio Andres Iglesias Martinez, Yan-Feng Wang, and Muamer Kadic Institut FEMTO-ST, UMR CNRS 6174, Univ. Bourgogne Franche-Comt´e, 25030 Besan¸con,France Department of Mechanics, School of Mechanical Engineering, Tianjin University, 300350 Tianjin,China
The propagation of acoustic or elastic waves in artificial crystals, including the case of phononic and soniccrystals, is inherently anisotropic. As is known from the theory of periodic composites, anisotropy is directlydictated by the space group of the unit cell of the crystal and the rank of the elastic tensor. Here, we examineeffective velocities in the long wavelength limit of periodic acoustic and elastic composites as a function of thedirection of propagation. We derive explicit and efficient formulas for estimating the effective velocity surfaces,based on second-order perturbation theory, generalizing the Christofell equation for elastic waves in solids.We identify strongly anisotropic sonic crystals for scalar acoustic waves and strongly anisotropic phononiccrystals for vector elastic waves. Furthermore, we observe that under specific conditions, quasi-longitudinalwaves can be made much slower than shear waves propagating in the same direction.
I. INTRODUCTION
Artificial crystals, when considered in the long wave-length limit, can be considered a sub-class of compositematerials , to which they add the property of spatial pe-riodicity and the existence of a space group describingthe symmetries of their unit-cell. Composite materialscan be assigned effective properties obtained by a limitingprocess, in the frame of homogenization theory. Homoge-nization has a long history and has been considered fromvarious physical and mathematical viewpoints . Com-posite structural mechanics often relies on the representa-tive volume element (RVE) approach, relating the inter-nal strain and stress fields to certain assumed boundaryconditions . Two-scale homogenization has a solidmathematical foundation and has been applied success-fully in various physical fields. As a framework, it isvalid for a general partial differential equation (PDE)and ultimately gives the limiting or homogenized PDE,and hence directly the effective material constants.In the case of periodic composites, a direct approachis to consider the dispersion relation, i.e. the bandstructure. Indeed, when both the frequency ω and thewavenumber k tend to zero, propagation becomes nondispersive and the function ω ( k ) = c eff k is linear. Start-ing from the Γ point of the first Brillouin zone, there isone non dispersive band for sonic crystals and three nondispersive bands for phononic crystals. Then one can fitthe dispersion relation to the form of the elastic tensor de-duced from the symmetries described by the space groupof the crystal. This is the approach of choice for elasticcomposites . A related empirical approach is to observeFabry-Perot oscillations in the transmission through a fi-nite crystal to estimate the effective velocity Elaborating upon the plane wave expansion (PWE)method that is used to compute the band structure of a) Electronic mail: [email protected] phononic crystals, Alevi et al. obtained the long wave-length limit for periodic elastic composites . In thecase of periodic acoustic composites, or sonic crystals,Krokhin et al. similarly obtained a PWE formula thatthey used to discuss the dependence of the effective veloc-ity with the filling fraction . For periodic elastic com-posites, Nemat-Nasser et al. proposed a more generalvariational approach where an appropriate functional ba-sis satisfying Bloch boundary conditions is considered .All these works did not consider explicitly anisotropy, asthe direction of propagation does not appear in the de-rived expressions. Moreover, an issue is that there is afull matrix to be inverted for each direction, which doesnot make the formulas obtained more efficient than adirect dispersion relation computation. The PWE ho-mogenization method was tentatively extended by vari-ous authors to the phononic crystal case, or of periodicelastic composites . A firm mathematical formula-tion, however, was not obtained before Torrent et al. .An appealing approach was provided by Kutsenko et al.who obtained a generalized Christofell equation for shearelastic waves in phononic crystals . Again, they didnot consider explicitly anisotropy.Our approach to the effective anisotropy of artificialcrystals is based on a variational formulation, as inthe case of two-scale homogenization, thus replacing inthe end the PWE implementation with a finite elementmethod. Similar to Krokhin’s and Kutsenko’s ap-proaches, we work directly with a second-order perturba-tion theory of the dispersion relation in periodic media.We obtain explicit formulas generalizing the Christofellequation for plane waves in homogeneous solids, that de-pend explicitly on the direction of propagation. Theformulas can be fitted against the form of the elastictensor that results from considering the space group ofthe crystal. We apply the theory to laminate, two-dimensional, and three-dimensional crystals of variousstructures. We identify strongly anisotropic sonic crys-tals for scalar acoustic waves and phononic crystals forvector elastic waves in which quasi-longitudinal waves are1 a r X i v : . [ phy s i c s . c l a ss - ph ] F e b uch slower than shear waves. II. EFFECTIVE VELOCITY FOR PERIODICACOUSTIC CONPOSITES
Bloch waves are the eigenfunctions of sonic crystalsand in general of periodic fluid composites. They havethe form p ( r ) exp( ı ( k · r − ωt )), with ω the angular fre-quency, k the wavevector, and p ( r ) the periodic part ofthe pressure field. They can be obtained by solving thetime-harmonic acoustic wave equation −∇ · (cid:18) ρ ∇ ( p exp( − ı k · r )) (cid:19) = ω B p exp( − ı k · r ) (1)under periodic boundary conditions. The mass density ρ ( r ) and the elastic modulus B ( r ) are inhomogeneousfunctions of space coordinates.In the finite element method, the eigenproblem defin-ing the band structure is solved in weak form as (cid:104) ( ∇ − ı k ) q, ρ − ( ∇ − ı k ) p (cid:105) = ω (cid:104) q, B − p (cid:105) , ∀ q. (2)In this equation q ( r ) exp( ı k · r ) is a test function definedin the same functional space as the solution ( q ( r ) is pe-riodic) and the symbol ∀ q means ’for all test functions’.The scalar product is defined for two scalar functions as (cid:104) a, b (cid:105) = (cid:82) Ω a ∗ b , with ∗ the complex conjugation opera-tion, and for two vector functions as (cid:104) a , b (cid:105) = (cid:82) Ω a ∗ · b .The left-hand side of Eq. (2) is thus (cid:90) Ω ( ∇ + ı k ) q ∗ · ρ − ( ∇ − ı k ) p. (3)The phononic band structure depicted in Fig. 1(a) isthe functional relation ω ( k ), obtained from Eq. (2). Fora sonic crystal, there is a single band starting from the Γpoint of the first Brillouin zone. For small frequency andwavenumber, that band is non dispersive but anisotropic:its slope, the effective velocity, depends on the directionof propagation. Plotting the effective velocity as a func-tion of the unit vector ˆ k defines the effective velocitysurface depicted in Fig. 1(b). Numerically, it is sufficientin order to obtain it to consider a small value for k andsolve Eq. (2) as a function of ˆ k , keeping only the lowesteigenvalue. A closed form expression, giving more physi-cal insight into the origin of anisotropy, can be obtainedas follows.We wish to consider an expansion for small wavenum-ber k = | k | and small frequency ω . From the point ofview of perturbation theory, the first-order solution for ω is zero, implying that we consider only the lowest bandstarting at the Γ point at the center of the first Brillouinzone, so we need a second-order solution in k and ω . Forthe Bloch wave itself, the first-order solution is enough.We consider the following ansatz for the periodic pressurefield to first-order p ( r ) ≈ p + ıkp ( r ; ˆ k ) (4) M G X k w G k ^ k ^ k ^(a) (b) fq Figure 1. Definition of the effective velocity surface for a pe-riodic acoustic composite, or sonic crystal. (a) The phononicband structure plotted along high symmetry directions in thefirst Brillouin zone (figured here by points X, M and Γ) hasone band starting at the Γ point in any direction. The slopeof that band is the effective velocity v eff (ˆ k ), a function of thedirection of propagation for acoustic waves given by unit vec-tor ˆ k . (b) The effective velocity surface is the locus of v eff (ˆ k ),a closed surface in three-dimensional space. with ˆ k = k /k a unit vector in the direction of propa-gation. p is a constant field since (cid:82) Ω ∇ q ∗ · ρ − ∇ p = 0for all test functions q implies ∇ p = 0 uniformly. As aresult ∇ p ≈ ık ∇ p and for instance( ∇ − ı k ) p ≈ − ık ˆ k p + ık ∇ p (5)to first order. As a result, the gradient of pressure is alinear function of the wavenumber that also depends onthe direction of propagation. Note that we do not needto consider an explicit dependence with frequency, sinceclose to the Γ point ω depends linearly on k – and alsodepends on the direction of propagation. As a result,the effective phase velocity v eff = ωk depends only on thedirection of propagation.Developing (2) we have (cid:104)∇ q, ρ − ∇ p (cid:105) − ık (cid:104)∇ q, ρ − ˆ k p (cid:105) + ık (cid:104) ˆ k q, ρ − ∇ p (cid:105) + k (cid:104) q, ρ − p (cid:105) = ω (cid:104) q, B − p (cid:105) , ∀ q (6)Then inserting the first-order approximation for the so-lution and keeping terms up to second order ık (cid:104)∇ q, ρ − ∇ p (cid:105) − ık (cid:104)∇ q, ρ − ˆ k p (cid:105) + k (cid:104)∇ q, ρ − ˆ k p (cid:105) − k (cid:104) ˆ k q, ρ − ∇ p (cid:105) + k (cid:104) q, ρ − p (cid:105) = ω (cid:104) q, B − p (cid:105) , ∀ q (7)The first two terms are of first order and the remainingterms of second order. They must be zero independently,since the equation is continuously valid for all k and ω .The two conditions are thus (cid:104)∇ q, ρ − ∇ p (cid:105) = (cid:104)∇ q, ρ − ˆ k p (cid:105) , ∀ q ; (8) v (cid:104) q, B − p (cid:105) = (cid:104) q, ρ − p (cid:105) + (cid:104)∇ q, ρ − ˆ k p (cid:105) − (cid:104) ˆ k q, ρ − ∇ p (cid:105) , ∀ q. (9)2quation (8) defines the first order correction p in theweak sense. Setting q = p it further follows (cid:104)∇ p , ρ − ∇ p (cid:105) = (cid:104)∇ p , ρ − ˆ k p (cid:105) = (cid:104) ˆ k p , ρ − ∇ p (cid:105) . (10)The last expression holds only if ρ is a real-valued func-tion. Finally, setting q = p in Eq. (9) we obtain anestimator for the square of the effective phase velocity v (ˆ k ) = (cid:104) p , ρ − p (cid:105) − (cid:104) ˆ k p , ρ − ∇ p (cid:105)(cid:104) p , B − p (cid:105) . (11)Equation (11) gives explicitly the effective velocity sur-face for acoustic pressure waves in the long wavelengthlimit. It is equivalent to Krokhin’s PWE formula ,but it avoids refering to the inversion of a full matrix.Actually, the matrix inversion is replaced by the solu-tion of the sparse linear problem defined by Eq. (8).Anisotropy is exclusively contained in the correction term (cid:104) ˆ k p , ρ − ∇ p (cid:105) that represents the part of the elastic po-tential energy of the Bloch wave that is stored in themicrostructure; i.e. that term vanishes only for an ho-mogeneous unit cell. If both B and ρ are real-valuedfunctions, including the case of lossless media, the lat-ter term is positive per Eq. (10) and we have the upperbound v ≤ (cid:104) p , ρ − p (cid:105)(cid:104) p , B − p (cid:105) , (12)i.e. the effective velocity is always smaller than the ratioof the averaged inverses of the mass density and the mod-ulus. As a consequence, the velocity surface is containedwithin a sphere whose radius is the square root of (12). III. EFFECTIVE TENSORS FOR PERIODICACOUSTIC COMPOSITES
In the case of fluid composites, Eq. (11) leads to ascalar effective value of the elastic modulus that can bedefined as B eff = (cid:104) p , p (cid:105)(cid:104) p , B − p (cid:105) . (13)That value is independent of the direction of propaga-tion. The numerator of Eq. (11) can be checked to bea quadratic form with respect to the direction vector ˆ k ,hence it defines a rank-2 effective tensor for the inverseof mass density, i.e.ˆ k · (cid:18) ρ (cid:19) eff ˆ k = (cid:104) ˆ k p , ρ − (ˆ k p − ∇ p ) (cid:105)(cid:104) p , p (cid:105) . (14)Thus it is the effective mass density that is anisotropicin the case of fluid composites. The effective tensor canbe checked to be symmetric and has the general form (cid:18) ρ (cid:19) eff = r r r . r r . . r (15) (b)(c)(a) a XY d XYZ -1500-1000-500 0 500 1000 1500-1500 -1000 -500 0 500 1000 1500 0 500 1000 1500 v e l o c i t y ( m / s ) velocity (m/s)(XY) plane(XZ) plane (YZ) plane Figure 2. A 2D square-lattice sonic crystal composed of trian-gular steel rods in water. (a) The triangular rods are rotatedby 45 ◦ with respect to the X axis. The ratio of the length ofthe equilateral triangle to the lattice constant is d/a = 0 . When the tensor is written in its principal axes, it be-comes diagonal and positive (cid:18) ρ (cid:19) eff = r r
00 0 r (16)There is a single longitudinal wave whatever the directionof propagation, satisfying the relation v (ˆ k ) = B eff ( r α + r β + r γ ) (17)with ( α, β, γ ) = (cos θ cos φ, cos θ sin φ, sin θ ) the compo-nents of ˆ k along the principal axes. When under thisform, fitting the effective velocity surface is very easy,since only the value of the phase velocity in three differ-ent directions is required.As a first example, we consider the 2D sonic crystal ofsteel rods in water whose unit cell is depicted in Fig. 2.For simplicity, steel is in this section considered as anequivalent fluid supporting only longitudinal waves. The3 able I. Effective constants for periodic acoustic composites.Effective constant B eff r r r Units GPa m / kg m / kg m / kgFig. 2 3.034 6 .
07 10 − .
03 10 − .
58 10 − Fig. 3 2.9 6 .
46 10 − .
46 10 − .
88 10 − Fig. 4 2 .
84 10 − . − .
416 0 . .
17 10 − .
85 10 − . − Fig. 6 2.2 1 .
94 10 − .
94 10 − .
865 10 −
0 1600 v e l o c i t y ( m / s )
0 1600 v e l o c i t y ( m / s ) v e l o c i t y ( m / s ) eff (m/s) (b)(c)(a) a XY d XYZ -1500-1000-500 0 500 1000 1500-1500 -1000 -500 0 500 1000 1500 0 500 1000 1500 v e l o c i t y ( m / s ) velocity (m/s)(XY) plane(XZ) plane (YZ) plane Figure 3. (a) A 2D hexagonal-lattice sonic crystal composedof triangular steel rods in water. The triangular rods arerotated by 10 ◦ with respect to the X axis. The ratio of thelength of the equilateral triangle to the lattice constant is d/a = 0 .
7. The crystal has a C symmetry and is transverseisotropic. (b) Effective velocity surface. (c) Cross-sectionsthrough the symmetry planes of the crystal. steel inclusions have a triangular shape and are organizedaccording to a square lattice. The structure is invariantalong the Z axis and has a vertical symmetry plane pass-ing along the diagonal of the square. Hence the crystalis orthotropic, with the first two principal axes rotatedby 45 ◦ in the ( XY ) plane. The material constants usedare ρ = 1000 kg / m and B = 2 . ρ = 7780 kg / m and B = 264 GPa for steel. The ve-
0 400 v e l o c i t y ( m / s )
0 400 v e l o c i t y ( m / s ) v e l o c i t y ( m / s )
0 50 100 150 200 250 300 350v eff (m/s) (b)(c)(a) a XY XYZ -400-300-200-100 0 100 200 300 400-400 -300 -200 -100 0 100 200 300 400 0 100 200 300 400 v e l o c i t y ( m / s ) velocity (m/s)(XY) plane(XZ) plane (YZ) plane Figure 4. (a) A laminar 1D sonic crystal composed of al-ternated layers of water and air with equal thickness. Thestructure is invariant along axes Y and Z . The crystal isorthotropic with two independent tensor elements. (b) Effec-tive velocity surface. (c) Cross-sections through the symmetryplanes of the crystal. locity surface has an almost circular cross-section in the( XY ) plane and an almost elliptical cross-section in allplanes containing the Z axis. The fitted effective con-stants in Table I confirm that r and r are almostequal, whereas r has a slightly larger value. We checkedthat the results are similar for other lattices and inclusionshapes: anisotropy remains quite limited for sonic crys-tals with an inclusion fully immersed in the surround-ing matrix. In the case of the hexagonal lattice and thesame inclusion but rotated by 10 ◦ , see Fig. 3, there isa C symmetry in addition to the invariance axis (the Z axis is a rotation center of order 3). The C symme-try imposes strictly r = r , a property that is verifiednumerically in Table I.The simplest acoustic composite with very stronganisotropy is a simple alternation of two very differentmaterials, for instance water and air; see Fig. 4. Thematerial constants used for air are ρ = 1 . / m and B = 142 kPa. X is an axis of revolution and the crystalis transverse isotropic. Of course, such a theoretical sonic4
0 1500 v e l o c i t y ( m / s )
0 1500 v e l o c i t y ( m / s ) v e l o c i t y ( m / s )
700 800 900 1000 1100 1200 1300 1400 1500v eff (m/s) (b)(c)(a) a XY XYZ -1500-1000-500 0 500 1000 1500-1500 -1000 -500 0 500 1000 1500 0 500 1000 1500 v e l o c i t y ( m / s ) velocity (m/s)(XY) plane(XZ) plane (YZ) plane a . a Figure 5. (a) A 2D sonic crystal composed of a periodic arrayof waveguides containing water. The crystal is orthotropicwith three independent tensor elements. (b) Effective velocitysurface. (c) Cross-sections through the symmetry planes ofthe crystal. crystal of air and water is not easily accessible to exper-iment. For the laminar case, the effective tensor (cid:16) ρ (cid:17) eff is known analytically . We checked that the formulas r = (cid:104) ρ (cid:105) − and r = r = (cid:104) ρ − (cid:105) match with the fittedresult in Table I for Fig. 4, where (cid:104) . (cid:105) denotes the spatialaverage.A feasible solution to obtain strongly anisotropic soniccrystals is to consider a single phase material, for instancewater, contained in a periodic array of solid tubes act-ing as acoustic waveguides without a frequency cut-off.We neglect here the generation of elastic waves in thesolid waveguides containing the fluid supporting acousticwaves. For instance, the square-lattice crystal of Fig. 5defines an orthotropic crystal with three different prin-cipal velocities. The phase velocity in the Z directionis faster than the phase velocity in the Y direction, be-cause acoustic waves have to propagate for a longer dis-tance from one side of the unit cell to another, and evenfaster than the phase velocity in the X direction. Thesituation is typical of labyrinthine sonic crystals or meta-
0 950 v e l o c i t y ( m / s )
0 950 v e l o c i t y ( m / s ) v e l o c i t y ( m / s )
600 650 700 750 800 850 900 950v eff (m/s) (b)(c)(a) a XYZ XYZ -1000 0 1000-1000 0 1000 0 400 800 v e l o c i t y ( m / s ) velocity (m/s)(XY) plane(XZ) plane (YZ) plane . a . a Figure 6. (a) A 3D sonic crystal composed of a periodic arrayof waveguides containing water. The crystal is orthotropicwith two independent tensor elements. (b) Effective velocitysurface. (c) Cross-sections through the symmetry planes ofthe crystal. materials used for sound absorption. Figure 6 shows a 3Dlabyrinthine sonic crystal containing water. That crystalis orthotropic with two independent tensor elements.
IV. EFFECTIVE VELOCITIES FORPERIODIC ELASTIC COMPOSITES
The derivation of the effective velocity formula for elas-tic composites, or phononic crystals, follows the samepath as for sonic crystals in the previous section, withthe added difficulty that the displacement field u i ( r ) isa vector field with three components. The vector elas-todynamic equation, here written in component form,replaces the scalar acoustic equation − ( c ijkl ( u k exp( − ı k · r )) ,l ) ,j = ω ρu i exp( − ı k · r ) . (18)The weak form of Eq. (18), valid for Bloch waves of theform u i ( r ) exp( ı ( − k · r + ωt )), is (cid:104) ( ∇ − ı k ) q , c : ( ∇ − ı k ) u (cid:105) = ω (cid:104) q , ρ u (cid:105) , ∀ q . (19)5he notation c : means contraction of the last two indicesthe rank-4 tensor c . The left-hand-side of Eq. (19) is incomponent form (cid:90) Ω ( ∂ j + ık j ) q ∗ i c ijkl ( ∂ l − ık l ) u k . (20)One difficulty in the vector (elastic) case is that thereis not a single value for the zero-th order constant fieldat zero frequency. Instead, for elasticity we have threepossible values, for each of the three different possiblepolarizations. In the phononic band structure, there arenow three different propagating bands starting from theΓ point. Therefore the ansatz for the displacement fieldup to the first order is taken as u ≈ ξ α ( u α + ik u α (ˆ k )) (21)where summation on α = 1 , , ω, k ) = 0 is of dimension 3. The three coefficients ξ α inthe linear combination are unknown. Instead of Eq. (8),the first-order corrections are obtained as the solution ofthe linear problems (cid:104)∇ q , c : ∇ u α (cid:105) = (cid:104)∇ q , c : ˆ ku α (cid:105) , ∀ q . (22)for α = 1 , ,
3. For the second-order terms, we now haveinstead of Eq. (9) v (cid:104) q , ρ u β (cid:105) ξ β = ξ β (cid:104) (cid:104) ˆ kq , c : ˆ ku β (cid:105) (23)+ (cid:104)∇ q , c : ˆ ku β (cid:105) − (cid:104) ˆ kq , c : ∇ u β (cid:105) (cid:105) , ∀ q for β = 1 , ,
3. As before, we select the three test function q = u ( α )0 to obtain a generalization of the Christofell’sequation for elastic waves in anisotropic homogeneousmedia v (cid:104) u α , ρ u β (cid:105) ξ β = (cid:104) (cid:104) ˆ ku α , c : ˆ ku β (cid:105) −(cid:104) ˆ ku α , c : ∇ u β (cid:105) (cid:105) ξ β . (24)This expression defines a 3 × to vector elastic waves, and con-tains the full anisotropy of wave propagation in the longwavelength limit. The implementation under a varia-tional form is much more efficient than PWE formulas ,because there is no matrix that needs to be inverted, onlytwo 3 × u α contain structural anisotropy and arisebecause of discontinuities at the inclusions or at internalboundaries. The formula has an explicit dependence onthe direction of propagation: it gives the three effectivevelocity surfaces directly. Each of the velocity surfacescan be assigned to either the longitudinal wave or oneof the two shear waves that exist in the long wavelengthlimit. V. EFFECTIVE TENSORS FOR PERIODICELASTIC COMPOSITES
In the case of elastic composites, Eq. (24) leads to ascalar effective value of the mass density if the vectors u α are chosen orthogonal. Then ρ eff = (cid:104) u , ρ u (cid:105)(cid:104) u , u (cid:105) , (25)where u equals any of the three u ( α )0 .The effective elastic tensor defined by Eq. (24) is sym-metric and of rank 4. Its general form in contractednotation is then( c ) eff = c c c c c c . c c c c c . . c c c c . . . c c c . . . . c c . . . . . c (26)With phononic crystals in the long wavelength limit, thesymmetry is given by the space group describing thesymmetries of the unit-cell considered a continuous dis-tribution of matter . This is in contrast to the pointgroup for crystal lattices composed of atoms assumed tobe punctual . We will not consider all possible spacegroups in the following, but only combinations of sym-metry planes. In case there is one symmetry plane, e.g.( x , x ), then( c ) eff = c c c c . c c c . . c c . . . c c . . . . c . . . . . c . (27)In case there are two orthogonal planes of symmetry, thecrystal is orthotropic and( c ) eff = c c c . c c . . c . . . c . . . . c . . . . . c . (28)If the crystal is transversely isotropic with respect to axis x then( c ) eff = c c c . c c . . c . . . c . . . . c . . . . . ( c − c ) . (29)Transverse isotropy is a sub-case of orthotropy.6 able II. Effective tensors for periodic elastic composites.Fig. 7 Fig. 8 Fig. 10¯ ρ (kg / m ) 4461 7780 7780 c (GPa) 14.66 22.49 83.80 c (GPa) 120.38 1.12 83.80 c (GPa) 120.68 31.31 1.29 c (GPa) 42.74 2.88 0.085 c (GPa) 2.91 7.23 0.085 c (GPa) 2.91 0.24 0.071 c (GPa) 7.12 0.54 0.98 c (GPa) 7.12 8.70 0.27 c (GPa) 34.90 0.63 0.27 -6000-3000 0 3000 6000-6000 -3000 0 3000 6000 0 2500 5000 v e l o c i t y ( m / s ) velocity (m/s)(XY): LS1 S2 (b)(a) a XY -6000-3000 0 3000 6000-6000 -3000 0 3000 6000 0 2500 5000 v e l o c i t y ( m / s ) velocity (m/s)(XZ): LS1 S2 -6000-3000 0 3000 6000-6000 -3000 0 3000 6000 0 2500 5000 v e l o c i t y ( m / s ) velocity (m/s)(YZ): LS1 S2 (d)(c) Figure 7. (a) A laminar 1D phononic crystal composed ofalternated layers of steel and epoxy with equal thickness. Thestructure is invariant along axes Y and Z . The crystal isorthotropic. (b-d) Cross-sections of the three effective velocitysurfaces through the symmetry planes of the crystal. The fits in the following figures are for curves withthe following expressions, valid for the ( XY ) plane oforthotropic crystals:¯ ρV L,S ( φ ) = 12 [( c + c ) α + ( c + c ) β ± (cid:0) [( c − c ) α − ( c − c ) β ] +4( c + c ) α β (cid:1) − ] , (30)¯ ρV SH ( φ ) = c α + c β . (31)Fitting of the velocity curves then provides an esti-mator for effective parameters ( c , c , c , c , c , c ).All effective parameters can be obtained by fitting ve-locity curves in the two additional planes ( XZ ) and ( Y Z ). Equations (30) and (31) indeed remain validwith a replacement of the former set of parameters with( c , c , c , c , c , c ) and ( c , c , c , c , c , c ),respectively. Redundancy in the effective parameters inthe fitting process is not a problem and instead helpsfinding more accurate estimates for the effective elas-tic tensor. Table II gathers the effective parameters ofthe periodic elastic composites considered next. Twoisotropic solid materials are considered in examples, steeland epoxy. Independent material constants for steel are c = 264 GPa, c = 84 GPa, and ρ = 7780 kg / m ;for epoxy they are c = 7 .
54 GPa, c = 1 .
48 GPa, and ρ = 1142 kg / m .The case of phononic crystals with a solid matrix leadsto some anisotropy for square-lattice crystals but trans-verse isotropy for hexagonal-lattice crystals . An alter-nation of epoxy and steel layers in a 1D phononic crys-tals, see Fig. 7, leads as in the case of the sonic crystalof Fig. 4 to strong anisotropy with orthotropic symme-try. Propagation in the plane ( Y Z ) is further isotropic.Overall, the longitudinal velocity remains always fasterthan the two shear waves. As a note, the laminar casecan be treated analytically, resulting in explicit formulasfor the effective elastic tensor : c ∗ = (cid:104) / ( λ + 2 µ ) (cid:105) − , c ∗ = c ∗ = (cid:104) /µ (cid:105) − ,c ∗ = (cid:104) µ (cid:105) , c ∗ = c ∗ = (cid:104) λ/ ( λ + 2 µ ) (cid:105)(cid:104) / ( λ + 2 µ ) (cid:105) − ,c ∗ = (cid:104) µλ/ ( λ + 2 µ ) (cid:105) + (cid:104) λ/ ( λ + 2 µ ) (cid:105) c ∗ , (32) c ∗ = c ∗ = (cid:104) µ ( λ + µ ) / ( λ + 2 µ ) (cid:105) + (cid:104) λ/ ( λ + 2 µ ) (cid:105) c ∗ , where λ and µ are Lam´e’s constants for isotropic mate-rials ( λ + 2 µ = c , µ = c ) and (cid:104) . (cid:105) denotes the spatialaverage. We checked that the fitted values appearing inTable II for Fig. 7 are consistent with the analyticalresult.The 2D phononic crystal in Fig. 8 uses the same meshas the sonic crystal of waveguides in Fig. 5. The longbeams now play the role of elastic waveguides, however.The structure becomes quite soft for longitudinal wavespropagating in the Y direction compared to the otherprincipal axes, i.e. c is much smaller than c and c ,as Table II indicates. Remarkably, c > c , so thatpure shear waves (polarized along the Z axis) in the Y direction are significantly faster than longitudinal waves.The in-plane shear wave is coupled with the longitudinalwave by the structure and remains always slower thanthat longitudinal wave. This property is consistent with c < c in Table II.The interplay of symmetry and anisotropy in 2D struc-tures is further illustrated in Fig. 9. The hexagonal-lattice crystal is made of a single phase of steel. Theinitial configuration in Fig. 9(a) is composed of threeidentical diamonds connected at the center and at threevertices of the boundary of the hexagonal unit-cell. Ithas three symmetry planes (and C symmetry). As aresult, elastic wave propagation in the plane ( XY ) isisotropic. The in-plane shear wave is very slow, whereasthe pure shear wave is just slightly slower than the longi-7 v e l o c i t y ( m / s ) velocity (m/s)(YZ): LS1 S2 -5000-2500 0 2500 5000-5000 -2500 0 2500 5000 0 2000 4000 v e l o c i t y ( m / s ) velocity (m/s)(XZ): LS1 S2 -5000-2500 0 2500 5000-5000 -2500 0 2500 5000 0 2000 4000 v e l o c i t y ( m / s ) velocity (m/s)(XY): LS1 S2 (b)(a) (d)(c) a XY . a a Figure 8. (a) A 2D phononic crystal composed of a periodicarray of steel bars. The crystal is orthotropic. (b-d) Cross-sections of the three effective velocity surfaces through thesymmetry planes of the crystal. In plane ( XY ), the longitu-dinal velocity becomes smaller than the S2 shear velocity ina certain angular range. tudinal wave. Then the central connection point is grad-ually shifted downward in Figs. 9(b-d), leaving only onevertical symmetry plane and the structure becomes or-thotropic. The in-plane shear wave always remains veryslow and the longitudinal wave in the Y direction be-comes slower and slower, and in any case slower than thepure-shear wave. This example illustrates how structurecontrols wave anisotropy.Considering again the 3D structure of beams with cu-bic lattice of Fig. 6 leads to the velocity surfaces shown inFig. 10. The structure is again orthotropic. Anisotropyis however in the case of elastic waves quite different tothe case of acoustic waves, due to the vector characterof the polarization. As a note, there is no decoupling ofin-plane and out-of-plane elastic waves in the 3D case,in contrast to the 2D case. There is a very slow shearwave for all directions of propagation. The longitudinaland the other shear waves are strongly anisotropic, butthe longitudinal wave always remains faster. In case thewaves are coupled by the structure, the velocity surfacesrepulse and do not cross. As a result of this topolog-ical property, that must be fulfilled for all propagationdirections defined on the unit sphere that forms a closedsurface in 3D space, longitudinal and shear velocity sur-faces are strictly imbricated in the case considered. -4000-2000 0 2000 4000-4000 -2000 0 2000 4000 0 2000 4000 v e l o c i t y ( m / s ) velocity (m/s)(XY): LS1 S2-4000-2000 0 2000 4000-4000 -2000 0 2000 4000 0 2000 4000 v e l o c i t y ( m / s ) velocity (m/s)(XY): LS1 S2 -3000-2000-1000 0 1000 2000 3000-3000-2000-1000 0 1000 2000 3000 0 1000 2000 3000 v e l o c i t y ( m / s ) velocity (m/s)(XY): LS1 S2 -2500-1250 0 1250 2500-2500 -1250 0 1250 2500 0 1250 2500 v e l o c i t y ( m / s ) velocity (m/s)(XY): LS1 S2 (b)(a) (d)(c) a XY a XY (f)(e) (h)(g) a XY a XY Figure 9. A 2D hexagonal-lattice phononic crystal composedof a periodic array of steel bars. (a) In the initial configura-tion, the C symmetry implies transverse symmetry. (b) Theeffective velocity surfaces are then transversely isotropic in the( XY ) plane. (c,e,g) The central connection point of the barsis brought down in steps of 0 . a , breaking the C symmetrybut leaving the symmetry plane ( Y Z ) intact, hence makingthe crystal orthotropic. (d,f,h) Corresponding cross-sectionsof the three effective velocity surfaces through the symmetryplanes of the crystal. In plane ( XY ), the longitudinal veloc-ity gradually becomes smaller than the S2 shear velocity in acertain angular range. v e l o c i t y ( m / s ) velocity (m/s)(XY): LS1 S2
0 3500 v e l o c i t y ( m / s )
0 3500 v e l o c i t y ( m / s ) v e l o c i t y ( m / s )
0 500 1000 1500 2000 2500 3000 3500v eff (m/s) (b)(a)
0 2000 v e l o c i t y ( m / s )
0 2000 v e l o c i t y ( m / s ) v e l o c i t y ( m / s )
0 500 1000 1500 2000 2500v eff (m/s)
0 450 v e l o c i t y ( m / s )
0 450 v e l o c i t y ( m / s ) v e l o c i t y ( m / s )
100 150 200 250 300 350 400 450v eff (m/s) -3500 0 3500-3500 0 3500 0 1500 3000 v e l o c i t y ( m / s ) velocity (m/s)(XZ): LS1 S2 -3500 0 3500-3500 0 3500 0 1500 3000 v e l o c i t y ( m / s ) velocity (m/s)(YZ): LS1 S2 (d)(c) (f)(e) XYZ XYZ XYZ
L waveS1 waveL waveS2 wave
Figure 10. A 3D phononic crystal composed of a periodicarray of steel bars, with the same mesh as in Figure 6. Thecrystal is orthotropic with three symmetry planes. (a,c,e)Effective velocity surfaces for the three elastic waves. (b,d,f)Cross-sections through the symmetry planes of the crystal.
VI. CONCLUSION
The main results of this work are the formulas (11)and (24) for the effective velocities of acoustic and elas-tic waves in periodic composites. Those formulas havea variational form similar to those produced by two-scale homogenization theory, but they were directly ob-tained from a second-order perturbation analysis of thephononic band structure of the physics of waves in peri-odic media. The influence of the microstructure, that isthe details of the internals of the crystal, is encompassedin a first order perturbation obtained as the solution ofan auxiliary problem on the unit-cell. The effective ten-sors are obtained from volume averages over the unit cellinvolving the zeroth order perturbation, here either aconstant pressure field or a constant displacement vec- tor field. Effective velocities depend continuously on thedirection of propagation and form effective velocity sur-faces characteristic of the crystal anisotropy in the longwavelength limit.Periodic acoustic composites, though sustaining scalarpressure waves in a fluid medium that is isotropic at themicroscopic level, can be made quite strongly anisotropicby a proper design of the structure of the unit-cell. Weparticularly point at possible realizations with periodicarrays of hollow waveguides forming labyrinths for thefundamental acoustic guided mode, which is dispersion-less and without frequency cut-off.In periodic elastic composites, the vector character ofwave polarization plays a determinant part. For 2D elas-tic composites for which in-plane and out-of-plane (pureshear) waves are decoupled, the longitudinal wave canbe made slower than the pure shear wave over a givenangular range by structural design with a single-phasematerial. For 3D elastic composites, the coupling of allthree components of the displacement field leads to im-bricated velocity surfaces.
ACKNOWLEDGMENTS
We acknowledge support by the EIPHI GraduateSchool (contract “ANR-17-EURE-0002”).
DATA AVAILABILITY STATEMENT
The data that support the findings of this study areavailable from the corresponding author upon reasonablerequest.
REFERENCES G. W. Milton,
The theory of composites (CambridgeMonographs on Applied and Computational Mathe-matics, Cambridge, UK, 2002). R. V. Craster and S. Guenneau, eds.,
Acoustic meta-materials: Negative refraction, imaging, lensing andcloaking , Springer Series in Materials Science, Vol. 166(Springer Science & Business Media, 2012). S. J. Hollister and N. Kikuchi, “A comparison of ho-mogenization and standard mechanics analyses for pe-riodic porous composites,” Comput. Mech. , 73–95(1992). A. Sridhar, V. G. Kouznetsova, and M. G. Geers, “Ho-mogenization of locally resonant acoustic metamateri-als towards an emergent enriched continuum,” Com-put. Mech. , 423–435 (2016). G. Allaire, “Homogenization and two-scale conver-gence,” SIAM Journal on Mathematical Analysis ,1482–1518 (1992).9 R. V. Craster, J. Kaplunov, and A. V. Pichugin, “High-frequency homogenization for periodic media,” Proc.R. Soc. A , 2341–2362 (2010). Y. Chen, M. Kadic, S. Guenneau, and M. Wegener,“Isotropic chiral acoustic phonons in 3D quasicrys-talline metamaterials,” Phys. Rev. Lett. , 235502(2020). X. Chen, J. Moughames, Q. Ji, J. A. I. Mart´ınez,H. Tan, S. Adrar, N. Laforge, J.-M. Cote, S. Euphrasie,G. Ulliac, et al. , “Optimal isotropic, reusable truss lat-tice material with near-zero Poisson’s ratio,” ExtremeMech. Lett. , 101048 (2020). Y. Chen, T. Frenzel, S. Guenneau, M. Kadic, andM. Wegener, “Mapping acoustical activity in 3D chiralmechanical metamaterials onto micropolar continuumelasticity,” J. Mech. Phys. Solids , 103877 (2020). F. Cervera, L. Sanchis, J. S´anchez-P´erez, V. R.Mart´ınez-Sala, C. Rubio, and F. Meseguer, “Refrac-tive acoustic devices for airborne sound,” Phys. Rev.Lett. , 023902 (2002). Z. Hou, F. Wu, X. Fu, and Y. Liu, “Effective elasticparameters of the two-dimensional phononic crystal,”Phys. Rev. E , 037604 (2005). P. Halevi, A. Krokhin, and J. Arriaga, “Photonic crys-tal optics and homogenization of 2D periodic compos-ites,” Phys. Rev. Lett. , 719 (1999). A. A. Krokhin, J. Arriaga, and L. N. Gumen, “Speed ofsound in periodic elastic composites,” Phys. Rev. Lett. , 264302 (2003). S. Nemat-Nasser, J. R. Willis, A. Srivastava, and A. V.Amirkhizi, “Homogenization of periodic elastic com-posites and locally resonant sonic materials,” Phys.Rev. B , 104103 (2011). Q. Ni and J. Cheng, “Anisotropy of effective velocity forelastic wave propagation in two-dimensional phononic crystals at low frequencies,” Phys. Rev. B , 014305(2005). Q. Ni and J. Cheng, “Long wavelength propagation ofelastic waves in three-dimensional periodic solid-solidmedia,” J. Appl. Phys. , 073515 (2007). J. Liu, Y. Wu, F. Li, P. Zhang, Y. Liu, andJ. Wu, “Anisotropy of homogenized phononic crystalswith anisotropic material,” Europhys. Lett. , 36001(2012). D. Torrent, Y. Pennec, and B. Djafari-Rouhani, “Reso-nant and nonlocal properties of phononic metasolids,”Phys. Rev. B , 174110 (2015). A. A. Kutsenko, A. L. Shuvalov, and A. N. Norris,“Evaluation of the effective speed of sound in phononiccrystals by the monodromy matrix method (L),” J.Acoust. Soc. Am. , 3553–3557 (2011). A. A. Kutsenko, A. L. Shuvalov, A. N. Norris, andO. Poncelet, “Effective shear speed in two-dimensionalphononic crystals,” Phys. Rev. B , 064305 (2011). A. A. Kutsenko, A. L. Shuvalov, and A. N. Norris, “Onthe quasistatic effective elastic moduli for elastic wavesin three-dimensional phononic crystals,” J. Mech. Phys.Solids , 2260–2272 (2013). V. Laude,
Phononic Crystals: Artificial Crystals forSonic, Acoustic, and Elastic Waves , 2nd ed. (DeGruyter, Berlin, 2020). A. Authier, “International tables for crystallography,”(Kluwer Academic Publishers, 2006) Chap. 1.1. Intro-duction to the properties of tensors, p. 3–33. G. W. Postma, “Wave propagation in a stratifiedmedium,” Geophysics20