Higher-Order Dispersion, Stability, and Waveform Invariance in Nonlinear Monoatomic and Diatomic Systems

Fronk, Matthew D.; Leamy, Michael J.

doi:10.1115/1.4036501

Recent studies have presented first-order multiple time scale approaches for exploring amplitude-dependent plane-wave dispersion in weakly nonlinear chains and lattices characterized by cubic stiffness. These analyses have yet to assess solution stability, which requires an analysis incorporating damping. Furthermore, due to their first-order dependence, they make an implicit assumption that the cubic stiffness influences dispersion shifts to a greater degree than the quadratic stiffness, and they thus ignore quadratic shifts. This paper addresses these limitations by carrying-out higher-order, multiple scales perturbation analyses of linearly damped nonlinear monoatomic and diatomic chains. The study derives higher-order dispersion corrections informed by both quadratic and cubic stiffness and quantifies plane wave stability using evolution equations resulting from the multiple scales analysis and numerical experiments. Additionally, by reconstructing plane waves using both homogeneous and particular solutions at multiple orders, the study introduces a new interpretation of multiple scales results in which predicted waveforms are seen to exist over all space and time, constituting an invariant, multiharmonic wave of infinite extent analogous to cnoidal waves in continuous systems. Using example chains characterized by dimensionless parameters, numerical studies confirm that the spectral content of the predicted waveforms exhibits less growth/decay over time as higher-order approximations are used in defining the simulations' initial conditions. Thus, the study results suggest that the higher-order multiple scales perturbation analysis captures long-term, nonlocalized invariant plane waves, which have the potential for propagating coherent information over long distances.

Introduction

Wave propagation in linear monoatomic, diatomic, and other chains is of contemporary interest due to their nontrivial dispersion, filtering, and frequency bandgap behavior. These systems may arise in the modeling of one-dimensional (1D) waveguides, or wave propagation in three-dimensional crystals along preferred directions, such as the directions along the unit cell edges of zincblende crystals (e.g., GaAs and NaCl). Nonlinearities in these systems introduce behavior of interest to the engineering community such as amplitude-dependent dispersion and group velocity associated with weak nonlinearities and solitary wave propagation (e.g., solitons) associated with strong nonlinearities. Specific engineering technology that has been proposed from studying wave propagation in periodic structures are switches [1], filters [2], and diodes [3].

Several recent studies have analyzed amplitude-dependent dispersion (nontrivial temporal/spatial frequency relationships) and other nonlinear behavior in discrete chains. Many of such studies have applied perturbation theory to uncover amplitude-dependent phenomena. Vakakis and King [4] employed perturbation techniques to study the monocoupled nonlinear media. They identified amplitude-dependent attenuation and propagation zones, which were associated with synchronous motion of standing waves and nonsynchronous motion of traveling waves, respectively. Chakraborty and Mallik [5] obtained first-order perturbation results for a cubically nonlinear monoatomic chain, with a discussion on the bounding frequencies of the propagation zone. The interaction of oppositely traveling waves and the effect of boundary conditions for semi-infinite and finite chains were also examined. In Ref. [6], Dubus et al. applied higher-order multiple scales to a linear monoatomic chain coupled (linearly) to another monoatomic chain with quadratic nonlinearities. Tuning the linear interchain coupling stiffness shifted the magnitude of a critical wavenumber associated with resonant energy transfer between the chains. Narisetti et al. [7] and Manktelow et al. [8,9] investigated wave propagation in cubically nonlinear monoatomic and diatomic chains using Lindstedt–Poincaré and multiple scales analyses, respectively, with an emphasis on amplitude-dependent dispersion shifts. Further, they identified wave-based devices which exploit bandgap shifting to enable tunable filtering and wave-guiding. In Ref. [10], dispersion shifts of a monoatomic chain governed by quadratic and cubic interactions are obtained using a Lindstedt–Poincaré approach. Furthermore, particular solutions are obtained at the first-order only and waveform invariance is not considered. For a more comprehensive review of tunable dispersion behavior in discrete periodic systems, see Ref. [11].

The study of higher harmonics arising in these discrete systems has centered on their generation, energy exchange, and/or decay. Cabaret et al. [12] analytically and experimentally studied second harmonic generation in a precompressed weakly nonlinear diatomic chain, reporting amplitude-dependent self-induced attenuation and self-induced transparency at the band edges. Sánchez-Morcillo et al. [13] conducted a study on the propagation of spectral content of harmonically driven compressional waves in quadratically nonlinear 1D granular chains subject to a static preload. An interplay between the fundamental and second harmonic, governed by dispersion, was observed for nonevanescent modes. Swinteck et al. [14] employed second-order perturbation studies to examine the self-interaction and wave–wave interactions in a monoatomic chain with quadratically interacting masses. Molecular dynamics simulations affirmed that vibrational modes exchange energy and have lifetimes that are inherently long. These studies discuss the generation, energy exchange, and decay/lifetime of the spectral content in the chains, respectively, different from the invariant spectral content studied herein.

In the long wavelength limit, discrete media behave as if continuous, and within such an approximation, studies have identified the presence of invariant waveforms. Invariant waveforms, or specific harmonic content that remains unaltered over the course of propagation, are of interest to the engineering community because of their ability to maintain information over long distances and times. They arise due to an interplay between dispersion and nonlinearity, the former of which assigns different wave speeds for different harmonics, and the latter of which can shift such wave speeds to maintain the wave's spatial waveform [15]. Higher-order asymptotic analysis is often employed when analyzing these waveforms. Laitone [16] extended an expansion by Friedrichs [17] to higher-orders to capture the effect of vertical motions and pressure variation of finite-amplitude, long water waves. At the first-order, the solution by Friedrichs is shown to be a cnoidal wave. Cnoidal waves appear to be less-studied by the engineering community as compared to other nonlinear invariant solutions. Figure 1 provides a classical example of cnoidal water waves captured in 1933 [18]. These waves are characterized by their sharp crests formed by nontrivial spectral content that propagates without changing. It is known that spatially extended, long-crested cnoidal waves satisfy the nonlinear Kortegweg de-Vries (KdV) equation [19] and travel unaltered, much like solitons but with infinite, periodic extent (vis-à-vis compact support). Kraenkel and Manna [20] applied higher-order multiple scales perturbation theory to the KdV equations to yield solutions for solitary waves. Some studies have presented continuum approximations to nonlinear lattices and recovered the KdV equation as a result. Friesecke and Mikikits-Leitner [21] demonstrated how the nonlinear interaction of particles in a discrete chain uncovered, in the continuum limit, continuous partial differential equations, such as the KdV equation and thus cnoidal waves. The dynamics of discrete particles are not considered for the duration of their proofs. Gaison et al. [22] applied homogenization techniques to the long wavelength limit of a nonlinear monoatomic lattice with periodic material properties. They uncovered solutions that approximately satisfy a KdV equation. Friesecke and Pego [23] presented that sonic speed wave pulses in nonlinear lattices are governed by the continuum limit regardless of the length scale. The studies in Refs. [21–23] provide an insight into how invariant waveforms may arise in these lattices in the long wavelength limit. However, no application of multiple scales to demonstrate invariance of waves in discrete lattices, for all wavelengths, has been found in the prior literature.

Fig. 1

View large Download slide

Cnoidal waves captured by the U.S. Army bombers over the Panama coast in 1933

Waveform invariance has received attention in recent years, primarily focused on the propagation of solitary waves and solitons. This is in contrast to single-frequency plane waves in nonlinear media, which are known to induce higher harmonics as they propagate. Sen et al. [24] provide an overview on the phenomenon of solitary waves in discrete media, highlighting differences from their characteristics in continuous media. They also discuss ways to manipulate solitary waves in granular chains for shock mitigation. Another review of solitons in nonlinear lattices is provided by Kartashov et al. [25] with an emphasis on theoretical implications. Their existence and stability in one, two, and three dimensions are considered. Experimental studies have confirmed the existence of solitary waves and solitons in discrete chains. Molerón et al. [26] investigated a 1D chain of repelling magnets analytically and experimentally. The analysis discussed a solitary wave whose profile and speed depends on amplitude. It was shown that in the low energy limit these wave profiles take the form of a KdV soliton with a sech² shape, whereas in the high energy limit, energy becomes further localized, and the waveform exhibits a hat function with width on the order of the lattice spacing. Daraio et al. [27] experimentally observed solitary waves propagating down a strongly nonlinear chain of Teflon spheres, obtaining a good agreement between theoretically estimated values (based on Hertz interaction: 3/2 force–displacement characteristic) of the wave speed and experimental measurements.

The stability of wave propagation in nonlinear continuous and discrete systems has also been analyzed in a number of studies. Such studies are of practical interest because they define how systems may deviate from their expected behavior or define critical thresholds for system parameters. Newton and Keller [28] analytically studied the stability of continuous periodic plane waves by examining wave amplitude. They arrive at an eigenvalue problem in which positive eigenvalues correspond to exponential growth in amplitude (unstable), and negative eigenvalues correspond to decay in amplitude (stable). Bickham et al. [29] studied the stability of localized modes in a monoatomic chain with quadratic and cubic force interactions. Numerical simulations introduced local static distortions into chains. Instabilities of the localized vibration modes began to arise in simulations in which the cubic anharmonicity increased to the point where the lowest possible plane wave frequency was reached. Further increases in the cubic anharmonicity led to behavior reminiscent of a triatomic molecule until further increases created instabilities yet again. Friesecke and Pego [30] identified exponential growth in amplitude of waves in discrete lattices by applying Floquet theory and analyzing evolution equations. In Ref. [31], Flach and Gorbach uncovered a threshold for tangent bifurcations in discrete breathers. Huang and Hu [32] investigated the stability of diatomic lattices using a quasi-discreteness approach. They derived evolution equations and fixed points for acoustic and optical modes and identified the presence of asymmetrical gap solitons. Gorbach and Johansson [33] presented findings on the linear stability of discrete gap breathers in a 1D diatomic chain with harmonic intersite potential and nonlinear external potential. Six types of instabilities were identified, which took on either oscillatory or nonoscillatory forms. An “inversion of stability” regime was found characterized by practically radiationless mobility. The authors would like to emphasize that discrete breathers studied in Refs. [31–33] are spatially localized and thus in contrast to the plane waves of infinite extent studied herein.

This paper builds upon recent studies of plane waves in nonlinear monoatomic and diatomic chains and other discrete systems by developing a higher-order multiple scales procedure to inform dispersion, stability, and waveform invariance. Previous efforts have focused primarily on first-order multiple scales analysis and have not yet considered both quadratic and cubic nonlinearities, which are of interest since they both necessarily arise in Taylor Series expansions of nonlinear interactions. With the higher-order analysis applied to these discrete systems, former assumptions neglecting quadratic nonlinearities on their amplitude-dependent dispersion shifts can now be rigorously justified. The reviewed literature also has not presented an analytical approach based on multiple scales to describe the loss of plane wave stability in nonlinear monoatomic and diatomic chains. This paper formulates a local fixed point analysis based on the multiple scales-derived evolution equations, and thus provides stability information in analytical form. It also numerically simulates nondimensionalized chains to identify and characterize their amplitude-dependent stability. Waveform invariance in discrete media has generally focused on localized waveforms such as solitary waves and solitons. To the best of our knowledge, no analysis has been put forth showing that multiple scales theory, applied to nonlinear monoatomic and diatomic chains, can capture solutions which asymptotically approach periodic, invariant waveforms akin to cnoidal waves. This paper provides both analytical and numerical evidences that a higher-order multiple scales approach yields such invariant waveforms. This contrasts, for example, with other equally valid interpretations [14] in which similar higher-order perturbation analyses are used instead to provide information about multiple wave scattering and phonon lifetimes.

System Descriptions

This section introduces the two systems studied herein, namely nonlinear monoatomic and diatomic chains, defining their equations of motion to which the asymptotic analyses are applied and presenting schematics for the ease of system visualization. For the monoatomic chain, the smallest repeatable subsystem, or unit cell, consists of a single mass and its connecting springs and dampers, as presented in Fig. 2.

Fig. 2

View large Download slide

Schematic of the unit cell for the monoatomic chain

A force balance on the jth mass yields its equation of motion

\begin{matrix} m {\ddot{x}}_{j} + k_{1} (2 x_{j} - x_{j + 1} - x_{j - 1}) - ε k_{2} {(x_{j + 1} - x_{j})}^{2} + ε k_{2} {(x_{j - 1} - x_{j})}^{2} \\ - ε k_{3} {(x_{j + 1} - x_{j})}^{3} - ε k_{3} {(x_{j - 1} - x_{j})}^{3} + ε c (2 {\dot{x}}_{j} - {\dot{x}}_{j + 1} - {\dot{x}}_{j - 1}) = 0 \end{matrix}

(1)

where $x_{j} = x (j, t)$ denotes the displacement from equilibrium of the jth mass, $m$ its mass, $k_{1}$ ⁠, $k_{2}$ ⁠, and $k_{3}$ are the linear, quadratic, and cubic stiffnesses, respectively, and $c$ denotes the linear damping coefficient. Time derivatives of $x_{j}$ — ${\dot{x}}_{j}$ and ${\ddot{x}}_{j}$ —denote the velocity and acceleration of the jth mass, respectively. The parameter $ε$ is assumed to be small and introduced herein as a bookkeeping device [34].

By contrast, the unit cell of the diatomic chain, depicted in Fig. 3, contains two masses coupled with identical springs and dampers.

Fig. 3

View large Download slide

Schematic of the unit cell for the diatomic chain

Summing the forces on each degree-of-freedom in the jth unit cell leads to its equation of motion

M {\ddot{\underline{x}}}_{j} + \sum_{l = j - 1}^{j + 1} K_{l} {\underline{x}}_{l} + \sum_{l = j - 1}^{j + 1} ε C_{l} \underline{\dot{x}} - ε {\underline{f}}_{quadratic} - ε {\underline{f}}_{cubic} = \underline{0}

(2)

where

{\underline{x}}_{j} = [\begin{matrix} x_{a} (j, t) \\ x_{b} (j, t) \end{matrix}]

(3)

M = (\begin{matrix} m_{a} & 0 \\ 0 & m_{b} \end{matrix})

(4)

K_{j - 1} = (\begin{matrix} 0 & - k_{1} \\ 0 & 0 \end{matrix}), K_{j} = (\begin{matrix} 2 k_{1} & - k_{1} \\ - k_{1} & 2 k_{1} \end{matrix}), K_{j + 1} = (\begin{matrix} 0 & 0 \\ - k_{1} & 0 \end{matrix})

(5)

C_{j - 1} = (\begin{matrix} 0 & - c \\ 0 & 0 \end{matrix}), C_{j} = (\begin{matrix} 2 c & - c \\ - c & 2 c \end{matrix}), C_{j + 1} = (\begin{matrix} 0 & 0 \\ - c & 0 \end{matrix})

(6)

{\underline{f}}_{quadratic} = [\begin{matrix} - k_{2} {(x_{b} (j - 1, t) - x_{a} (j, t))}^{2} + k_{2} {(x_{b} (j, t) - x_{a} (j, t))}^{2} \\ - k_{2} {(x_{a} (j, t) - x_{b} (j, t))}^{2} + k_{2} {(x_{a} (j + 1, t) - x_{b} (j, t))}^{2} \end{matrix}]

(7)

{\underline{f}}_{cubic} = [\begin{matrix} k_{3} {(x_{b} (j - 1, t) - x_{a} (j, t))}^{3} + k_{3} {(x_{b} (j, t) - x_{a} (j, t))}^{3} \\ k_{3} {(x_{a} (j, t) - x_{b} (j, t))}^{3} + k_{3} {(x_{a} (j + 1, t) - x_{b} (j, t))}^{3} \end{matrix}]

(8)

It is noted here that the presence of linear and small cubic and quadratic force–displacement relationships describes many systems whose nonlinear interactions are expanded in a Taylor series.

Analysis Approach

Method of Multiple Scales: Monoatomic Chain.

In this section, a higher-order method of multiple scales (MMS) is developed for the monoatomic chain presented in Sec. 2, extending the analysis through the second-order to capture harmonic content and time evolution of amplitude and phase. It is shown in Sec. 5.3, that the inclusion of these higher-order features increases the invariance of the cnoidal-like wave. Higher-order results also facilitate the stability and dispersion analyses, providing insight to behavior not evident at solely the first-order.

Multiple time scales are assumed to exist such that

T_{0} = t, T_{1} = ε t, T_{2} = ε^{2} t, \dots, T_{n} = ε^{n} t

(9)

where $t$ denotes the original time scale, and $T_{n}$ represents the nth time scale. Since $ε ≪ 1$ ⁠, each time scale advances more slowly than its predecessor. Introducing these times scales plays an important role in detecting amplitude-dependent dispersion shifts since phase change at each slow time scale can be interpreted as dispersion corrections.

In accordance with these time scales, differentiation with respect to time contains multiple orders

\dot{()} = D_{0} () + ε D_{1} () + ε^{2} D_{2} () + \dots + ε^{n} D_{n} ()

(10)

where $D_{n} ()$ denotes an operator defined to represent differentiation with respect to $T_{n}$ ⁠.

It follows that

\begin{matrix} \ddot{()} = {(D_{0} () + ε D_{1} () + ε^{2} D_{2} () + \dots + ε^{n} D_{n} ())}^{2} \\ = D_{0}^{2} () + 2 ε D_{0} D_{1} () + ε^{2} D_{1}^{2} () + 2 ε^{2} D_{0} D_{2} () + O (ε^{3}) \end{matrix}

(11)

The MMS seeks a series expansion of the solution of the form

x_{j} = x_{j}^{(0)} + ε x_{j}^{(1)} + ε^{2} x_{j}^{(2)} + \dots + ε^{n} x_{j}^{(n)}

(12)

where

x_{j}^{(n)}

denotes the nth-order solution. Thus, the asymptotic solution is expanded about a known zeroth-order solution

x_{j}^{(0)}

⁠. Substituting Eqs. – into Eq. and collecting matching orders of

ε

yields

n + 1

equations. The first three are presented here

ε^{0} : m D_{0}^{2} x_{j}^{(0)} + k_{1} (2 x_{j}^{(0)} - x_{j + 1}^{(0)} - x_{j - 1}^{(0)}) = 0

(13)

\begin{array}{l} ε^{1} : m D_{0}^{2} x_{j}^{(1)} + k_{1} (2 x_{j}^{(1)} - x_{j + 1}^{(1)} - x_{j - 1}^{(1)}) = - 2 m D_{0} D_{1} x_{j}^{(0)} - c D_{0} (2 x_{j}^{(0)} - x_{j + 1}^{(0)} - x_{j - 1}^{(0)}) \\ + k_{2} {(x_{j + 1}^{(0)} - x_{j}^{(0)})}^{2} - k_{2} {(x_{j - 1}^{(0)} - x_{j}^{(0)})}^{2} + k_{3} {(x_{j + 1}^{(0)} - x_{j}^{(0)})}^{3} + k_{3} {(x_{j - 1}^{(0)} - x_{j}^{(0)})}^{3} \end{array}

(14)

\begin{array}{l} ε^{2} : m D_{0}^{2} x_{j}^{(2)} + k_{1} (2 x_{j}^{(2)} - x_{j + 1}^{(2)} - x_{j - 1}^{(2)}) = - 2 m D_{0} D_{2} x_{j}^{(0)} - m D_{1}^{2} x_{j}^{(0)} - 2 m D_{0} D_{1} x_{j}^{(1)} - c D_{1} (2 x_{j}^{(0)} - x_{j + 1}^{(0)} - x_{j - 1}^{(0)}) - c D_{0} (2 x_{j}^{(1)} - x_{j + 1}^{(1)} - x_{j - 1}^{(1)}) + f_{quadratic} (x_{j - 1}, x_{j}, x_{j + 1}) + f_{cubic} (x_{j - 1}, x_{j}, x_{j - 1}) \end{array}

(15)

where

\begin{matrix} f_{quadratic} (x_{j - 1}, x_{j}, x_{j + 1}) = - 2 k_{2} x_{j - 1}^{(0)} x_{j - 1}^{(1)} + 2 k_{2} x_{j - 1}^{(0)} x_{j}^{(1)} + 2 k_{2} x_{j - 1}^{(1)} x_{j}^{(0)} - 2 k_{2} x_{j + 1}^{(1)} x_{j}^{(0)} \\ - 2 k_{2} x_{j + 1}^{(0)} x_{j}^{(1)} + 2 k_{2} x_{j + 1}^{(0)} x_{j + 1}^{(1)} \end{matrix}

(16)

\begin{matrix} f_{cubic} (x_{j - 1}, x_{j}, x_{j + 1}) = 3 k_{3} x_{j - 1}^{(0)}^{2} x_{j - 1}^{(1)} - 3 k_{3} x_{j - 1}^{(0)}^{2} x_{j}^{(1)} - 6 k_{3} x_{j - 1}^{(0)} x_{j - 1}^{(1)} x_{j}^{(0)} + 6 k_{3} x_{j - 1}^{(0)} x_{j}^{(0)} x_{j}^{(1)} + 3 k_{3} x_{j - 1}^{(1)} x_{j}^{(0)}^{2} - 6 k_{3} x_{j}^{(0)}^{2} x_{j}^{(1)} + 3 k_{3} x_{j + 1}^{(1)} x_{j}^{(0)}^{2} + 6 k_{3} x_{j + 1}^{(0)} x_{j}^{(0)} x_{j}^{(1)} - 6 k_{3} x_{j + 1}^{(0)} x_{j + 1}^{(1)} x_{j}^{(0)} - 3 k_{3} x_{j + 1}^{(0)}^{2} x_{j}^{(1)} + 3 k_{3} x_{j + 1}^{(0)}^{2} x_{j + 1}^{(1)} \end{matrix}

(17)

Note that the left-hand sides of Eq. (13) through (15) are similar, i.e., the linear kernels of the collected equations take the same form. Also note that each order is effectively forced by the solution obtained from the previous order.

By examination of Eq. , it is apparent that

x_{j}^{(0)}

takes the form of a Bloch wave since its governing equation is linear with periodic coefficients

x_{j}^{(0)} = \frac{1}{2} A e^{i ω_{0} T_{0}} e^{- i μ j} + c . c .

(18)

where $A$ denotes the wave's amplitude, $ω_{0}$ and $μ$ are the wave's fundamental temporal and spatial frequencies, respectively, and $c . c .$ denotes the complex conjugate of all preceding terms.

Since

A

is generally a complex quantity, it is useful to express it in polar form

A = α e^{i β}

(19)

where $α = α (T_{1}, T_{2}, \dots, T_{n})$ and $β = β (T_{1}, T_{2}, \dots, T_{n})$ are real quantities. Wave amplitude is a critical parameter as it determines the strength of the nonlinearities and consequently the amount of dispersion shifting and the relevancy of the MMS results.

The zeroth-order dispersion relationship follows from the Bloch solution form Eq. and the governing equation

ω_{0} = \sqrt{\frac{2 k_{1}}{m} (1 - cos μ)}

(20)

This relationship, which stems from the linear, undamped chain (⁠ $ε = 0)$ ⁠, shifts for waves at different amplitudes in the presence of nonlinearities.

With

x_{j}^{(0)}

known,

x_{j}^{(1)}

is found by updating Eq. with Eq. . Secular terms, i.e., those that contain

e^{i ω_{0} T_{0}} e^{- i μ j}

⁠, or its complex conjugate, induce resonant responses and thus cause the series solution to become nonuniform over the time scales of interest. Consequently, they are set to zero

0 = k_{3} α^{3} (\frac{3}{4} cos 2 μ - 3 cos μ + \frac{9}{4}) + ci ω_{0} α (1 - cos μ) + i m ω_{0} D_{1} α - m ω_{0} α D_{1} β

(21)

which, after separating real and imaginary parts and setting to zero, restricts the evolution of

α

and

β

with respect to

T_{1}

D_{1} α = - γ α

(22)

D_{1} β = δ α^{2}

(23)

where

γ = \frac{c}{m} (1 - cos μ)

(24)

δ = \frac{3 k_{3}}{4 m ω_{0}} (cos 2 μ - 4 cos μ + 3)

(25)

Equations (22) and (23) represent evolution equations describing the dynamics of the wave on a slower time scale. Phase evolution governs dispersion shifts, whereas amplitude evolution affects wave stability as will be presented in Secs. 5.1 and 5.2.1, respectively.

With secular terms removed, Eq. reduces to a linear, first-order ordinary differential equation, with constant coefficients and harmonic forcing at

2 ω_{0}

and

3 ω_{0}

⁠. The homogeneous solution is ignored since it contains no information than that already given by

x_{j}^{(0)}

⁠, so a particular solution is sought of the form

x_{j}^{(1)} = \frac{1}{2} B_{1} e^{2 i ω_{0} T_{0}} e^{- 2 i μ j} + \frac{1}{2} C_{1} e^{3 i ω_{0} T_{0}} e^{- 3 i μ j} + c . c .

(26)

where

B_{1}

and

C_{1}

can be found using the method of undetermined coefficients

B_{1} = \frac{i k_{2} A^{2} (sin 2 μ - 2 sin μ)}{4 m ω_{0}^{2} + 2 k_{1} cos 2 μ - 2 k_{1}} = b_{1} A^{2}

(27)

C_{1} = \frac{k_{3} A^{3} (3 cos 2 μ - cos 3 μ - 3 cos μ + 1)}{18 m ω_{0}^{2} + 4 k_{1} cos 3 μ - 4 k_{1}} = c_{1} A^{3}

(28)

The amplitudes in Eqs. (27) and (28) describe the amount of second- and third-harmonic content that make-up the first-order plane wave solution.

Turning attention to the second-order problem, Eq. , with known solutions for

x_{j}^{(0)}

and

x_{j}^{(1)}

⁠, the same procedure is carried out. Removal of secular terms now determines the evolution of

α

and

β

with respect to

T_{2}

D_{2} α = \frac{γ δ}{ω_{0}} α^{3}

(29)

D_{2} β = X α^{4} + Y α^{2} + Z

(30)

where

X = - \frac{1}{4 m ω_{0}} (k_{3} c_{1} (9 cos μ - 9 cos 2 μ + 3 cos 3 μ - 3) + 2 δ^{2} m)

(31)

Y = - \frac{k_{2} b_{1}}{i m ω_{0}} (2 sin μ - sin 2 μ)

(32)

Z = - \frac{γ^{2}}{2 ω_{0}}

(33)

The particular solution for

x_{j}^{(2)}

takes the form

x_{j}^{(2)} = \frac{1}{2} B_{2} e^{2 i ω_{0} T_{0}} e^{- 2 i μ j} + \frac{1}{2} C_{2} e^{3 i ω_{0} T_{0}} e^{- 3 i μ j} + \frac{1}{2} E_{2} e^{4 i ω_{0} T_{0}} e^{- 4 i μ j}

+ \frac{1}{2} F_{2} e^{5 i ω_{0} T_{0}} e^{- 5 i μ j} + c . c .

(34)

where

B_{2} = b_{21} α^{4} e^{2 i β} + b_{22} α^{2} e^{2 i β}

(35)

C_{2} = c_{21} α^{5} e^{3 i β} + c_{22} α^{3} e^{3 i β}

(36)

E_{2} = e_{2} A^{4}

(37)

F_{2} = f_{2} A^{5}

(38)

and

\begin{matrix} b_{21} = \frac{1}{2 e^{i μ} (4 m ω_{0}^{2} e^{2 i μ} + k_{1} e^{4 i μ} - 2 k_{1} e^{2 i μ} + k_{1})} (3 k_{3} b_{1} e^{6 i μ} - 2 k_{2} c_{1} e^{6 i μ} - 6 k_{3} b_{1} e^{5 i μ} + 2 k_{2} c_{1} e^{5 i μ} - \\ 16 m b_{1} ω_{0} δ e^{3 i μ} - 3 k_{3} b_{1} e^{4 i μ} + 2 k_{2} c_{1} e^{4 i μ} + 12 k_{3} b_{1} e^{3 i μ} - 3 k_{3} b_{1} e^{2 i μ} - 2 k_{2} c_{1} e^{2 i μ} - \\ 6 k_{3} b_{1} e^{i μ} - 2 k_{2} c_{1} e^{i μ} + 3 k_{3} b_{1} + 2 k_{2} c_{1}) \end{matrix}

(39)

b_{22} = \frac{4 i c b_{1} ω_{0} (1 - cos 2 μ) - 8 i b_{1} γ m ω_{0}}{4 (m ω_{0}^{2} + k_{1} {cos}^{2} μ - k_{1})}

(40)

\begin{matrix} c_{21} = - \frac{3}{2 (k_{1} e^{6 i μ} + 9 m ω_{0}^{2} e^{3 i μ} - 2 k_{1} e^{3 i μ} + k_{1})} c_{1} (\begin{matrix} - 16 k_{3} e^{3 i μ} {cos}^{4} μ + 2 k_{3} e^{6 i μ} - k_{3} e^{5 i μ} + 16 k_{3} e^{3 i μ} {cos}^{2} μ + \\ 12 e^{3 i μ} δ m ω_{0} + 2 k_{3} e^{4 i μ} - 6 k_{3} e^{3 i μ} + 2 k_{3} e^{2 i μ} - k_{3} e^{i μ} + 2 k_{3} \end{matrix}) \end{matrix}

(41)

c_{22} = - \frac{36 i γ c_{1} m ω_{0} + i c c_{1} ω_{0} (12 cos 3 μ - 12) + i b_{1} k_{2} (- 4 sin 3 μ + 4 sin 2 μ + 4 sin μ)}{2 (9 m ω_{0}^{2} + 2 k_{1} cos 3 μ - 2 k_{1})}

(42)

e_{2} = - \frac{b_{1} k_{3} (12 cos μ + 6 cos 4 μ - 12 cos 3 μ - 6) + i c_{1} k_{2} (8 sin μ - 8 sin 4 μ + 8 sin 3 μ)}{8 (8 m ω_{0}^{2} + k_{1} cos 4 μ - k_{1})}

(43)

f_{2} = \frac{c_{1} k_{3} (3 cos 2 μ - 6 cos μ + 6 cos 4 μ - 3 cos 5 μ - 3 cos 3 μ + 3)}{50 m ω_{0}^{2} + 4 k_{1} cos 5 μ - 4 k_{1}}

(44)

These amplitudes build upon the higher harmonic content introduced at the first-order to describe more completely the multiharmonic solutions that exist in this nonlinear media. As will be discussed in Sec. 5.3, this higher-order solution, when added to the first- and zeroth-order solutions, propagates more invariantly, i.e., with less growth and decay in its spectral content.

After substituting the results from Eqs. (18), (26), and (34) into Eq. (12), the solution for $x_{j}$ is now known up to $O (ε^{3})$ ⁠.

Method of Multiple Scales: Diatomic Chain.

Next, the MMS is developed for the diatomic chain. As with the monoatomic chain, higher-order results reveal the magnitudes of the spatially and temporally invariant harmonics induced by the quadratic and cubic nonlinearities as well as amplitude and phase evolution equations utilized in the stability and dispersion shift analyses. It is again assumed that there are multiple time scales at which the system evolves, and a series expansion of the solution is assumed

{\underline{x}}_{j} = {\underline{x}}_{j}^{(0)} + ε {\underline{x}}_{j}^{(1)} + ε^{2} {\underline{x}}_{j}^{(2)} + \dots + ε^{n} {\underline{x}}_{j}^{(n)}

(45)

where

{\underline{x}}_{j}^{(n)}

denotes the nth-order solution vector

{\underline{x}}_{j}^{(n)} = [\begin{matrix} x_{a}^{(n)} (j, t) \\ x_{b}^{(n)} (j, t) \end{matrix}]

(46)

Collecting matching orders of

ε

in Eq. yields

n + 1

vector equations. The first two are presented here

ε^{0} : M D_{0}^{2} {\underline{x}}_{j}^{(0)} + \sum_{l = j - 1}^{j + 1} K_{l} {\underline{x}}_{l}^{(0)} = \underline{0}

(47)

ε^{1} : M D_{0}^{2} {\underline{x}}_{j}^{(1)} + \sum_{l = j - 1}^{j + 1} K_{l} {\underline{x}}_{l}^{(1)} = - 2 M D_{0} D_{1} {\underline{x}}_{j}^{(0)} + {\underline{f}}_{damping, 1} + {\underline{f}}_{quadratic, 1} + {\underline{f}}_{cubic, 1}

(48)

{\underline{f}}_{damping, 1} = [\begin{matrix} c D_{0} (x_{b}^{(0)} (j, t) - x_{a}^{(0)} (j, t)) + c D_{0} (x_{b}^{(0)} (j - 1, t) - x_{a}^{(0)} (j, t)) \\ c D_{0} (x_{a}^{(0)} (j + 1, t) - x_{b}^{(0)} (j, t)) + c D_{0} (x_{a}^{(0)} (j, t) - x_{b}^{(0)} (j, t)) \end{matrix}]

(49)

{\underline{f}}_{quadratic, 1} = [\begin{matrix} k_{2} {(x_{b}^{(0)} (j, t) - x_{a}^{(0)} (j, t))}^{2} - k_{2} {(x_{b}^{(0)} (j - 1, t) - x_{a}^{(0)} (j, t))}^{2} \\ k_{2} {(x_{a}^{(0)} (j + 1, t) - x_{b}^{(0)} (j, t))}^{2} - k_{2} {(x_{a}^{(0)} (j, t) - x_{b}^{(0)} (j, t))}^{2} \end{matrix}]

(50)

{\underline{f}}_{cubic, 1} = [\begin{matrix} k_{3} {(x_{b}^{(0)} (j, t) - x_{a}^{(0)} (j, t))}^{3} + k_{3} {(x_{b}^{(0)} (j - 1, t) - x_{a}^{(0)} (j, t))}^{3} \\ k_{3} {(x_{a}^{(0)} (j + 1, t) - x_{b}^{(0)} (j, t))}^{3} + k_{3} {(x_{a}^{(0)} (j, t) - x_{b}^{(0)} (j, t))}^{3} \end{matrix}]

(51)

By inspection of Eq. , the zeroth-order solution takes the Bloch form

{\underline{x}}_{j}^{(0)} = \frac{1}{2} A \underline{ϕ} e^{i ω_{0} T_{0}} e^{- i μ j} + c . c .

(52)

where $A$ denotes the complex amplitude and $\underline{ϕ}$ a wave mode shape.

Substituting Eq. into Eq. leads to an eigenvalue problem for two wave modes and their dispersion relationship,

ω_{0} (μ)

⁠. There are two independent solutions, corresponding to the acoustic and optical branches

\underline{ϕ_{a c}} = [\begin{matrix} \frac{(1 + e^{i μ}) m_{b}}{m_{b} - m_{a} + \sqrt{m_{a}^{2} + m_{b}^{2} + 2 m_{a} m_{b} cos μ}} \\ 1 \end{matrix}]

(53)

\underline{ϕ_{opt}} = [\begin{matrix} 1 \\ \frac{m_{a} - m_{b} + \sqrt{m_{a}^{2} + m_{b}^{2} + 2 m_{a} m_{b} cos μ}}{- (1 + e^{i μ}) m_{b}} \end{matrix}]

(54)

ω_{0, a c} = \sqrt{\frac{k_{1} (m_{b} + m_{a})}{m_{a} m_{b}} - \frac{k_{1}}{2} \sqrt{\frac{4 m_{b}^{2} + 4 m_{a}^{2} + 8 m_{a} m_{b} cos μ}{m_{a}^{2} m_{b}^{2}}}}

(55)

ω_{0, opt} = \sqrt{\frac{k_{1} (m_{b} + m_{a})}{m_{a} m_{b}} + \frac{k_{1}}{2} \sqrt{\frac{4 m_{b}^{2} + 4 m_{a}^{2} + 8 m_{a} m_{b} cos μ}{m_{a}^{2} m_{b}^{2}}}}

(56)

As with the monoatomic chain, these zeroth-order dispersion relationships completely characterize the linear, undamped chain (⁠ $ε = 0$ ⁠). Higher-order results add amplitude-dependent corrections to these relationships due to the nonlinearities.

It is worth noting that at $μ = π$ ⁠, when $m_{b} > m_{a}$ ⁠, the optical mode shape simplifies from an indeterminate form to $[\begin{matrix} 1 \\ 0 \end{matrix}]$ through application of L'Hopital's rule. Also, the mode shapes presented here are the complex conjugates of those in Ref. [8] due to a difference in sign convention for the zeroth-order solution in Eq. (52).

The known form of

{\underline{x}}_{j}^{(0)}

is then substituted into Eq. . Projecting the secular terms against the Hermitian conjugate of the corresponding mode shape leads to evolution expressions for

D_{1} α

and

D_{1} β

D_{1} α = - γ_{d} α

(57)

D_{1} β = δ_{d} α^{2}

(58)

where

γ_{d} = γ_{d} (μ, c, \underline{ϕ_{a}}, m_{a}, m_{b})

(59)

δ_{d} = δ_{d} (μ, k_{3}, \underline{ϕ_{a}}, \underline{ϕ_{b}}, m_{a}, m_{b}, ω_{0})

(60)

The functional dependence only is provided although complete expressions can be found in a straight-forward manner using a symbolic manipulator. Such dependence is significant because of arguments made about the influence of different chain parameters on dispersion shifts and invariance in Secs. 5.1 and 5.3, respectively.

A particular solution for

{\underline{x}}_{j}^{(1)}

can be found using the form

{\underline{x}}_{j}^{(1)} = \frac{1}{2} \underline{B_{1}} e^{2 i ω_{0} T_{0}} e^{- 2 i μ j} + \frac{1}{2} \underline{C_{1}} e^{3 i ω_{0} T_{0}} e^{- 3 i μ j} + c . c .

(61)

where the vectors

\underline{B_{1}}

and

\underline{C_{1}}

are found using the method of undetermined coefficients

\underline{B_{1}} = [\begin{matrix} b_{1, a} A^{2} \\ b_{1, b} A^{2} \end{matrix}]

(62)

\underline{C_{1}} = [\begin{matrix} c_{1, a} A^{3} \\ c_{1, b} A^{3} \end{matrix}]

(63)

where

b_{1, a} = b_{1, a} (μ, m_{a}, m_{b}, k_{1}, k_{2}, \underline{ϕ_{a}}, \underline{ϕ_{b}}, ω_{0})

(64)

b_{1, b} = b_{1, b} (μ, m_{a}, m_{b}, k_{1}, k_{2}, \underline{ϕ_{a}}, \underline{ϕ_{b}}, ω_{0})

(65)

c_{1, a} = c_{1, a} (μ, m_{a}, m_{b}, k_{1}, k_{3}, \underline{ϕ_{a}}, \underline{ϕ_{b}}, ω_{0})

(66)

c_{1, b} = c_{1, b} (μ, m_{a}, m_{b}, k_{1}, k_{3}, \underline{ϕ_{a}}, \underline{ϕ_{b}}, ω_{0})

(67)

These amplitudes detail how much second and third harmonics arise due to the quadratic and cubic nonlinearities, respectively, and will be shown to make up an invariant waveform.

Note that this solution may not exist very close to $μ$ values associated with internal resonance—i.e., where $(2 ω_{0}, 2 μ)$ or $(3 ω_{0}, 3 μ)$ satisfy the zeroth-order dispersion relationship in Eq. (55) or (56), inducing a resonant response. For these cases, the particular solution becomes unbounded and must be removed using an additional (e.g., detuning) parameter. Chains exhibiting this behavior are omitted in this work, but this phenomenon can be explored in greater detail in future studies.

Reapplying the procedure on the second-order yields evolution equations for

α

and

β

of the following form:

D_{2} α = D_{2} α (α, μ, m_{a}, m_{b}, k_{2}, k_{3}, c, γ_{d}, δ_{d}, \underline{ϕ_{a}}, b_{1, a}, b_{1, b}, c_{1, a}, c_{1, b}, ω_{0})

(68)

D_{2} β = D_{2} β (α, μ, m_{a}, m_{b}, k_{2}, k_{3}, c, γ_{d}, δ_{d}, \underline{ϕ_{a}}, b_{1, a}, b_{1, b}, c_{1, a}, c_{1, b}, ω_{0})

(69)

The second-order particular solution then takes the form

\begin{matrix} {\underline{x}}_{j}^{(2)} = \frac{1}{2} \underline{B_{2}} e^{2 i ω_{0} T_{0}} e^{- 2 i μ j} + \frac{1}{2} \underline{C_{2}} e^{3 i ω_{0} T_{0}} e^{- 3 i μ j} + \frac{1}{2} \underline{E_{2}} e^{4 i ω_{0} T_{0}} e^{- 4 i μ j} \\ + \frac{1}{2} \underline{F_{2}} e^{5 i ω_{0} T_{0}} e^{- 5 i μ j} + c . c \end{matrix}

(70)

where the vectors $\underline{B_{2}}$ ⁠, $\underline{C_{2}}$ ⁠, $\underline{E_{2}}$ ⁠, and $\underline{F_{2}}$ can be found using the method of undetermined coefficients.

As with the monoatomic chain, second-order analysis adds to the spectral content at the second and third harmonics and introduces new spectral content at the fourth and fifth harmonics further defining the multiharmonic invariant wave. It is also noted here that the diatomic results can recover the monoatomic results by setting $m_{a} = m_{b} = m$ and replacing $μ$ with $μ / 2$ ⁠.

Numerical Simulation

This section briefly presents the procedure for numerically injecting perturbation-based waveforms into numerical simulations of Eqs. (1) and (2), forming the foundation to the numerical results presented in the remaining sections. They rely on initial conditions associated with the multiple scales solutions, assigning the order of terms at which to truncate the multiharmonic wave solution. In doing so, this procedure exhibits more control over the ensuing wave propagation than the approach used in Refs. [7] and [9], whereby the chain is excited with an arbitrary forcing until a measurable response is generated. It also allows for measurement of the validity of the MMS solutions since waveforms existing in the medium should persist after they are initially set, as detailed in Sec. 5.3. This section also introduces dimensionless parameters that combine the effects of nonlinear stiffness and wave amplitude into convenient expressions useful for simulating chains and interpreting results.

Displacement and velocity initial conditions required to numerically integrate Eqs. (1) and (2) are found using the series solutions in Eqs. (12) and (45), respectively. Thus, $x (j, 0) = {x_{j} |}_{t = 0}$ and $\dot{x} (j, 0) = {(d / d t) (x_{j}) |}_{t = 0}$ ⁠, where ${x_{j} |}_{t = 0}$ is truncated to a desired order $n$ ⁠, and $(d / d t) ()$ must correspond to the total time derivative defined in Eq. (10), truncated at the same order $n$ ⁠.

By carrying the expansions out and truncating, first- and second-order initial conditions can be determined. To emulate an infinite domain, the number of unit cells $N$ should be sufficiently large—generally on the order of 1800–2000—such that $μ N ≫ 1$ ⁠. Within this domain, analysis of a chain's behavior remains at a central region—usually the middle 800–1000 unit cells—such that reflection from the boundaries does not propagate to the region within the time scale of the simulation. In these studies, either fixed–fixed or free–free boundary conditions are used.

A given simulation is characterized by the following dimensionless parameters introduced here:

τ = ω_{0} t, Π_{1} = \frac{ε c ω_{0}}{k_{1}}, Π_{2} = \frac{ε k_{2} α_{0}}{k_{1}}, Π_{3} = \frac{ε k_{3} α_{0}^{2}}{k_{1}}, Π_{4} = \frac{m_{b}}{m_{a}}

(71)

where $Π_{1}, Π_{2},$ and $Π_{3}$ can be interpreted as relative strengths of the damping, quadratic nonlinearity, and cubic nonlinearity, respectively, and $Π_{4}$ ⁠, a mass ratio, applies only to the diatomic chain. The parameters $Π_{2}$ and $Π_{3}$ are particularly useful since they combine the effects of the nonlinear stiffness and wave amplitude into single parameters, both of which affect the magnitude of dispersion shifts. Note, however, that $Π_{2}$ and $Π_{3}$ are proportional to—not exact measurements of—the strength of the quadratic and cubic nonlinearities, respectively. That is to say, these parameters are not an exact ratio of spring forces because $α_{0}$ is the zeroth-order wave amplitude and not necessarily the spring stretch, $x_{j + 1} - x_{j}$ ⁠. Thus, $Π_{2}, Π_{3} < 0.1$ are conservative bounds on the MMS applicability (weak nonlinearity assumption). In general, the parameters in Eq. (71) facilitate parameter definition as outlined in the following procedure outline.

These steps were carried out for each simulation:

Specify $Π_{1}, Π_{2},$ and $Π_{3}$ values. Set $m$ and $k_{1}$ equal to 1, without loss of generality. For diatomic simulations, prescribe $m_{a}, m_{b}$ values (based on a desired $Π_{4}$ value) as well as a branch of excitation (acoustic or optical).
Specify a propagation constant $μ$ restricted to the first irreducible Brillouin zone. A corresponding $ω_{0}$ value can then be assigned to it according to Eqs. (20), (55), or (56). Note that this is not guaranteed to be the ensuing temporal frequency of the chain but merely a parameter used to evaluate initial conditions. For diatomic simulations, $\underline{ϕ_{a c}}$ or $\underline{ϕ_{opt}}$ can now be evaluated, depending on the branch of excitation that was selected.
Set $k_{3}$ equal to 1 (hardening) or −1 (softening) and then solve for $c, k_{2},$ and $α_{0}$ values which satisfy the desired $Π_{1}, Π_{2},$ and $Π_{3}$ values, respectively.
Define the number of unit cells to emulate an infinite system (Usually $N = 1800$ or $2000$ such that $μ N ≫ 1$ ⁠)
Assign an initial displacement and velocity distribution to the chain according to Eqs. (12) or (45) evaluated at time $t = 0$ and truncated at a specified order $n$ (⁠ $n = 0, 1,$ or $2$ ⁠). The $ε$ parameter is chosen to be small (e.g., $ε = 0.1$ ⁠).
Set two boundary conditions. For each end, use either fixed—zero displacement—or free—zero force.
Simulate the equations of motion by numerically integrating Eq. (1) or (2) using matlab's ODE45 routine.
For analyzing the results, study only a fixed middle portion of the chain (e.g., the middle 800–1000 unit cells) such that boundary effects can be ignored. Visually inspect the results to make sure that the boundary effects have not propagated into this central region during the time range of interest.

Results

Dispersion Shifts.

This section details amplitude-dependent dispersion shifts arising from both cubic and quadratic nonlinearities, a finding useful for device design [7]. It will be shown that the MMS analysis uncovers dispersion shifts to the zeroth-order (⁠ $ε = 0$ ⁠, linear undamped) relationship. The measurement of these predicted shifts in numerical simulations affirms the validity of the MMS results at the first and second-orders and demonstrates that quadratic nonlinearities govern shifts at higher-orders only, i.e., an effect not captured by a first-order analysis.

At the zeroth-order, the temporal frequency $ω_{0}$ relates to the spatial frequency $μ$ through dispersion relationships represented by Eqs. (20), (55), and (56). Figure 4 displays these zeroth-order expressions for the monoatomic and diatomic chains. Note that the monoatomic chain has a cutoff for $ω_{0} > 2 \sqrt{(k_{1} / m)}$ ⁠, and the diatomic chain exhibits a bandgap.

Fig. 4

View large Download slide

Zeroth-order dispersion relationships for the monoatomic and diatomic chains

Closed-form, amplitude-dependent dispersion shifts will now be derived from the MMS results. For undamped chains,

D_{1} α = D_{2} α = 0

⁠. Thus, the wave amplitude is constant (⁠

α = α_{0}

⁠). For the monoatomic chain, the evolution equations for phase, Eqs. and , reduce to

D_{1} β = δ α_{0}^{2}, c = 0

(72)

D_{2} β = X α_{0}^{4} + Y α_{0}^{2} + Z, c = 0

(73)

Similar expressions result for the diatomic chain by replacing $δ$ in Eq. (72) with $δ_{d}$ from Eq. (60), and replacing the left side of Eq. (73) with the expression in Eq. (69) evaluated at $α = α_{0}$ ⁠.

Reconstituting the phase expression using

D_{1} β

and

D_{2} β

in Eq. and returning to the original time scale, leads to the wave's temporal phase behavior

\dot{β} = ε D_{1} β + ε^{2} D_{2} β = ε δ α_{0}^{2} + ε^{2} (X α_{0}^{4} + Y α_{0}^{2} + Z)

(74)

for the monoatomic chain; a similar expression occurs for the diatomic chain.

Since

D_{1} β

and

D_{2} β

do not explicitly depend on time, they can be interpreted as corrections to the frequency at their respective order, yielding first- and second-order-accurate corrections

ω_{O (ε^{1})} = ω_{0} + ε D_{1} β

(75)

ω_{O (ε^{2})} = ω_{0} + ε D_{1} β + ε^{2} D_{2} β

(76)

Inspection of Eqs. (72), (73), (75), and (76) shows that the quadratic nonlinearity does not shift the frequency until the second-order correction is included, implying a small influence as compared to the cubic nonlinearity, which appears at both the first- and second-orders. For positive coefficients $k_{2}$ and $k_{3}$ ⁠, the shifts further increase the frequency (for a given $μ$ ⁠) with increasing wave amplitude $α_{0}$ ⁠. Such effects can inform the design of tunable devices (filters, waveguides, acoustic logic devices, etc.) which selectively include or exclude waves based on wave gain and subsequent bandgap shifting [7,8].

For both monoatomic and diatomic chains, numerical studies of undamped systems confirm the theoretical shifts as documented in Fig. 5. For different $μ$ values prescribed in simulation initial conditions, the first- and second-order MMS-predicted temporal frequencies are compared to the measured dominant temporal frequency of a central unit cell in the chain. Recall that assigning a $μ$ value does not guarantee the ensuing temporal frequency but rather is used to define all initial parameters. Thus, the simulations given first- and second-order initial conditions (⁠ $n = 1, 2$ ⁠) were measured to have dominant temporal frequencies close to the first- and second-order perturbation corrections (⁠ $ω_{O (ε^{1})}$ ⁠, $ω_{O (ε^{2})}$ ⁠), respectively. Theoretical expressions for the zeroth-order, unshifted temporal frequency, $ω_{0}$ (Eq. (20)—monoatomic and Eqs. (55) and (56)—diatomic) are also plotted as references to visualize the shifts that occur. For each simulation, the dominant temporal frequency is computed by taking an fast Fourier transform (FFT) of the time history of a central mass in the chain. To collect accurate data, temporal sampling frequencies are many integer multiple times higher than the expected shifted frequency, e.g., $f_{sampling} = 600 ((ω_{O (ε^{1})}) / 2 π)$ or $f_{sampling} = 600 ((ω_{O (ε^{2})}) / 2 π)$ ⁠.

Fig. 5

View large Download slide

First- and second-order dispersion shifts in the monoatomic (left) and diatomic (right) chains. Monoatomic: $Π_{1} = 0, Π_{2} = 0.04$ (top), $Π_{2} = 0$ (bottom), $Π_{3} = 0.04$ ⁠, Diatomic: $Π_{1} = 0, Π_{2} = 0.04$ (top), $Π_{2} = 0$ (bottom), $Π_{3} = 0.04, Π_{4} = 1.5$ ⁠.

Inspection of the bottom two subplots in Fig. 5 reveals that $k_{2}$ bears little to no effect on the first-order dispersion correction. As $Π_{2}$ changes from 0.04 to 0, which corresponds to the no quadratic nonlinearity present, the first-order perturbation and simulation results change by an average amount of 0.043%. For lightly damped chains $D_{1} α, D_{2} α \neq 0$ ⁠, corrections can be interpreted as instantaneous dispersion shifts since the amplitude decays slowly. The nonconstant reconstituted expressions for $\dot{β}$ are then interpreted similar to Eqs. (75) and (76) at each instance of time.

Stability Considerations.

The stability of wave propagation is of primary concern for designing devices informed by these discrete systems. Herein, instability is defined as any qualitative shift from the multiharmonic plane wave solution predicted by the MMS. Thus, it is important to characterize the types of instabilities that take place in these systems and the conditions at which they occur. An analytical stability analysis is first presented to see if the evolution equations for the amplitude reveal any information on stability loss. Informed by this analysis, numerical studies are carried out to determine the onset of instabilities.

Local Stability Analysis.

Because it is a critical parameter that determines the strength of the nonlinearities in the system, amplitude is chosen to be analyzed to assess the stability of the wave. It will be shown that the first-order results reveal only a damping-induced decay of amplitude (stable solution), whereas the inclusion of second-order results reveals a possible growth in amplitude (unstable solution). Though it is outside the conservative bounds of the MMS solution, the unstable solution provides a foundation to the numerical stability studies presented in Sec. 5.2.2.

Applying Eq. to wave amplitude

α

and recalling that

D_{0} α = 0

leads to its reconstituted evolution equation

\dot{α} = ε D_{1} α + ε^{2} D_{2} α

(77)

which suggests the existence of fixed points, or

α^{*}

such that

\dot{α} = 0

⁠. These fixed points are found from

{\dot{α} |}_{α^{*}} = 0

(78)

A local stability analysis can then be performed near each fixed point by introducing small perturbations from

α^{*}

⁠. The result of doing so is that the stability of each fixed point depends on the sign of

λ

arising from

\frac{d}{d α} {\dot{α} |}_{α^{*}} \equiv λ

(79)

where

λ > 0

yields an unstable fixed point,

λ < 0

a stable fixed point, and

λ = 0

is neutrally stable. Thus, the fixed points for the monoatomic chain can be readily found by substituting in the values from Eqs. and into Eqs. and

α^{*} = 0, \pm \sqrt{\frac{ω_{0}}{ε δ}}

(80)

The fixed point

α^{*} = \pm \sqrt{(ω_{0} / ε δ)}

exists only when

δ > 0

(cubic hardening) since

α

is assumed to be strictly real. Application of Eq. yields

{λ |}_{α^{*} = 0} = - ε γ

(81)

{λ |}_{α^{*} = \pm \sqrt{\frac{ω_{0}}{ε δ}}} = 2 ε γ

(82)

The fixed point $α^{*} = 0$ is of the highest relevancy since (1) amplitudes near it are small-enough to satisfy the conservative weakly nonlinear ansatz, $Π_{2}, Π_{3} < 0.1$ ⁠; (2) this is the fixed point to which asymptotically predicted plane waves will be attracted to in the presence of damping, therefore, its stability establishes plane wave stability of interest herein. Since, ${λ |}_{α^{*} = 0} < 0$ for damped chains (⁠ $c > 0$ ⁠), the fixed point $α^{*} = 0$ is stable. Note that $α^{*} = 0$ is the only fixed point informed by a first-order analysis as $α^{*} = \pm \sqrt{(ω_{0} / ε δ)}$ originates from second-order terms. Also note that the presence of damping enables the stability of the fixed points to be classified: undamped chains yield ${λ |}_{α^{*} = 0} = {λ |}_{α^{*} = \pm \sqrt{(ω_{0} / ε δ)}} = 0$ ⁠, which is neutrally stable.

By similar arguments, since Eq. is positive for damped chains (⁠

c > 0

⁠), the fixed point

α^{*} = \pm \sqrt{(ω_{0} / ε δ)}

is unstable. However, closer inspection of this fixed point reveals that it is at a large-enough amplitude that its accuracy cannot be expected. After substituting Eqs. and into the expression for

α^{*}

⁠, a critical

Π_{3}

value results dependent on only

μ

Π_{3, crit} = \frac{ε k_{3} α^{*^{2}}}{k_{1}} = \frac{8 (1 - cos μ)}{3 (cos 2 μ - 4 cos μ + 3)}

(83)

Examining Eq. (83), the lowest value of $Π_{3, crit}$ in the first Brillouin zone, $0 \leq μ \leq π$ ⁠, is approximately 0.66, which is beyond the conservative weak nonlinearity ansatz used for MMS analysis. However, as evident in the results presented in Sec. 5.2.2, the existence of such a fixed point is likely since large-amplitude plane waves are observed to grow unboundedly in numerical simulations.

The stability of the monoatomic chain is summarized in Fig. 6. As indicated, only initial wave amplitude, and not initial phase, determines the basins of attraction.

Fig. 6

View large Download slide

Fixed points and basins of attraction for the monoatomic chain

Stability of the diatomic chain can also be studied analytically. Its reconstituted evolution equation for

α

can be written as

\dot{α} = P_{1} α^{5} + P_{2} α^{3} + P_{3} α

(84)

where

P_{1}, P_{2}, and P_{3}

denote constants known from Eqs. and . Fixed-points can be solved for by setting

\dot{α} = 0

⁠. Thus, a given system admits at most five fixed points; however, because

α

is assumed to be real, imaginary fixed points are ignored. Inspection of Eq. indicates that

α^{*} = 0, \pm {[- \frac{P_{2}}{2 P_{1}} \pm \frac{1}{2} \sqrt{\frac{P_{2}^{2}}{P_{1}^{2}} - \frac{4 P_{3}}{P_{1}}}]}^{\frac{1}{2}}

(85)

where it is again emphasized that

α^{*}

must be strictly real. The stability of a given fixed point can be studied by examining the sign of each value of

λ

λ = 5 P_{1} α^{* 4} + 3 P_{2} α^{* 2} + P_{3}

(86)

Results are similar to the monoatomic chain. For $α^{*} \neq 0$ ⁠, many values of $λ$ can be computed using the $P_{1}, P_{2},$ and $P_{3}$ values associated with a prescribed set of $Π_{1}, Π_{2}, Π_{3},$ and $Π_{4}$ values and branch. The fixed point $α^{*} = 0$ is always stable for $c > 0$ ⁠. Other fixed points are unstable but have corresponding $Π_{2}$ and $Π_{3}$ values never simultaneously less than 0.1. Therefore, they are of large-enough amplitude to not be of interest to plane waves predicted herein.

Numerical Stability Studies.

This section details numerical studies that investigate stability beyond the applicable range of the MMS. It also confirms the general trends predicted by the local stability analysis, which is that at low amplitudes, multiharmonic plane waves are stable, and at high amplitudes (or, equivalently, large nonlinearities) the plane wave solutions predicted may lose stability. Note that such a conclusion cannot be drawn from only first-order information. These studies quantify the stability of the predicted waveforms over a wide range of $Π_{2}$ and $Π_{3}$ values for a given value of $Π_{1}$ and, if diatomic, $Π_{4}$ ⁠. Simulations begin with the initial conditions corresponding to the series solutions in Eqs. (12) and (45) for the monoatomic and diatomic chain, respectively. An unstable response is defined to be any deviation from the nature of the plane wave propagation predicted by the perturbation analysis. Thus, such studies aim to identify the stability limits of the MMS-predicted multiharmonic solutions. Although different types of instability are observed, it is not the primary concern of this work to explore further the ensuing solutions.

One such type of instability that occurred in numerical simulations is characterized by a spreading of the wave's spectral content. Chains with large $Π_{3}$ values exhibit this behavior. The spectral content of these waves initially contains integer multiples of $ω_{0}$ and $μ$ ⁠, but following longer simulation times, the plane wave nature of the response breaks down into an assortment of incommensurate frequencies. Figure 7 depicts this phenomenon.

Fig. 7

View large Download slide

An example chain exhibiting a breakdown in spectral content: $Π_{1} = 0.01, Π_{2} = 0, Π_{3} = 3.7$ ⁠, second-order initial conditions (ICs)

The other type of instability observed in simulations is characterized by an unbounded growth in the wave's amplitude. Chains with large positive or large negative $Π_{2}$ values, as well as chains with large negative $Π_{3}$ values, exhibit this behavior. The initial amplitude is unstable such that, as time progresses, the wave propagates with increasingly larger amplitudes. Note that this is similar to the exponential growth in amplitude instability presented in Ref. [28]. Figure 8 depicts this phenomenon. It is noted that all oscillations cease when the amplitude starts to grow.

Fig. 8

View large Download slide

An example chain exhibiting a growth in amplitude: $Π_{1} = 0.01, Π_{2} = 1, Π_{3} = 0$ ⁠, second-order ICs

For automating instability detection in the simulated results, criteria are developed which determine when one of the two types occur, thus signifying a loss of stability, or deviation from the multiharmonic plane wave solution. When either one of the two instability types occur, the parameter set is considered “unstable,” without distinguishing exactly which type of instability developed. Thus, an instability occurs when a set number of peaks in a spatial FFT exceeds a certain fraction of the maximum peak at that instant or when the instantaneous amplitude grows to a certain multiple of the simulation's original amplitude. A time limit is prescribed for simulations: if an instability does not arise within the prescribed time, or $τ = 12 (2 π)$ for these simulations, the parameter set is considered stable. For the diatomic chain, instability criteria are applied to $x_{a}$ and $x_{b}$ separately.

Figure 9 presents the stability results for the monoatomic chain as a function of $Π_{2}$ and $Π_{3}$ ⁠. The parameters $Π_{2, crit}$ and $Π_{3, crit}$ are the values of $Π_{2}$ and $Π_{3}$ ⁠, accurate to the nearest tenth, at which an unstable response is first observed when the other nonlinear $Π$ -value is set to zero (i.e., $Π_{2, crit}$ occurs when $Π_{3} = 0$ and $Π_{3, crit}$ occurs when $Π_{2} = 0$ ⁠). Each plot verifies the local stability analysis; namely, that plane waves characterized by small values of $Π_{2} / Π_{2, crit}$ and $Π_{3} / Π_{3, crit}$ are stable. Figure 10 presents similar stability results for the diatomic chain applied to $x_{a}$ ⁠. Plots for the stability study of $x_{b}$ depict similar trends and are shown in Figs. 13–15 in the Appendix.

Fig. 9

View large Download slide

Monoatomic chain stability study: second-order ICs, $Π_{1} = 0.01, μ = (π / 4)$

Fig. 10

View large Download slide

Diatomic chain stability study: $x_{a}$ ⁠, second-order ICs, $Π_{1} = 0.01, Π_{4} = 1.5, μ = (π / 4)$ ⁠, acoustic (a), and optical (b)

For some chains simulated using free–free boundary conditions, a rigid body (DC) harmonic appears that causes the two masses in each unit cell to separate by a constant amount, in addition to their expected oscillation. The location of the center of mass remains unchanged, as expected to satisfy the conservation of momentum. However, since the magnitude of this separation is generally small for the parameters being studied, such responses are considered stable.

Although the amplitudes at which point stability is lost in the numerical simulations do not align with the nonzero fixed points predicted by multiple scales approach, the results of these numerical studies qualitatively confirm the findings. Wave amplitudes below a critical value decay over time toward zero due to the presence of damping. Above this critical amplitude value, the wave will become unstable by either growing unboundedly or decomposing its original multiharmonic content into incommensurate frequencies.

Invariance.

This section describes how higher-order MMS results converge to an invariant waveform, a finding that offers new insight for predicting and modeling the wave propagation in these weakly nonlinear discrete systems and further validation of the MMS approach. In the absence of damping, the presented perturbation approach does not distinguish between any two points in space, nor any two points in time. Thus, in the limit as the asymptotic expansions are taken to higher orders, and assuming convergent vis-à-vis divergent behavior, it can be hypothesized that the waveforms should propagate for all space and time without change—i.e., they should be an invariant form akin to cnoidal waves studied in other contexts. It is worth noting that such an argument does not rely on a continuum approximation as in other studies [21–23]. Note that this is a very different perspective taken from other works carrying-out perturbation approaches in wave propagation problems [14,35], which more conventionally analyze higher harmonic generation and/or the closely related concept of phonon generation and lifetimes.

To test the invariant hypothesis, numerical integration of Eqs. (1) and (2) is carried out using the procedure presented in Sec. 4. Thus, plane waves corresponding to the asymptotic expansions in Eqs. (12) and (45) truncated at different orders are injected into the simulated chains. During simulation, the spatial harmonic content of the waves, is quantified using FFTs at different instances in time to track the generation (or loss) of the initially prescribed harmonic content. As mentioned in Sec. 4, to emulate an infinite chain, the analysis was restricted to a central region of the chain so that boundary effects could be ignored. According to the hypothesis, waves incorporating higher-order terms should exhibit measurably less variation, and hence propagate more invariantly. Only $Π_{3} > 0$ chains were sampled since it was found while conducting the numerical stability studies in Sec. 5.2.2, that the negative $Π_{3}$ values become unstable at much lower magnitudes than positive ones.

For monoatomic chains, a parameter that measures the spatial variance of a simulated waveform is proposed as

σ_{p}^{(n)} \equiv \frac{1}{L} \sum_{l = 1}^{L} | 1 - a_{l, p} |

(87)

where $a_{l, p}$ denotes the amplitude of the pth harmonic at the lth instant in time normalized by its initially prescribed value, $L$ is the number of time instances that make up the numerical solution, and $n$ denotes the order of perturbation terms retained in Eq. (12). Since $σ_{p}^{(n)}$ measures variance, lower values of $σ_{p}^{(n)}$ indicate more invariant waveforms, and thus, it is expected (if the hypothesis holds) that $σ_{p}^{(n)}$ decreases as $n$ increases. Figure 11 verifies the expected reduction in variance for monoatomic chains given second versus first-order initial conditions—all values of (⁠ $Π_{2}$ ⁠, $Π_{3}$ ⁠) simulated experience decreases in variance. Because Fig. 11 plots the normalized reduction in variance— $σ_{p}^{(1)} - σ_{p}^{(2)} / σ_{p, \max}^{(1)} (p = 2, 3)$ —a value of 1.0 indicates elimination of all variance using second-order initial conditions, while a value of 0.0 indicates no variance reduction.

Fig. 11

View large Download slide

Reduction of variance of the (a) second and (b) third harmonics, monoatomic chain, $μ = (π / 4), Π_{1} = 0$

For diatomic chains, the variance parameter specializes to each degree-of-freedom: $σ_{p, a}^{(n)}$ and $σ_{p, b}^{(n)}$ for $x_{a}$ and $x_{b}$ ⁠, respectively. Similar numerical simulations are carried out and the second and third spatial harmonics are tracked over time for waves on the acoustic and optical branches. Figure 12 presents the variance results for $x_{a}$ ⁠; the Appendix provides the results for $x_{b}$ ⁠. The trend in variance reduction is again present.

Fig. 12

View large Download slide

Reduction of variance in $x_{a}$ ⁠, second (a, c) and third (b,d) harmonics, diatomic chain. $Π_{4} = 2, μ = (π / 4), Π_{1} = 0$ acoustic (a and b), and optical (c and d).

Fig. 13

View large Download slide

Diatomic chain stability study: $x_{b}$ ⁠, second-order ICs, $Π_{1} = 0.01, Π_{4} = 1.5, μ = (π / 4),$ acoustic (a), and optical (b)

Fig. 14

View large Download slide

Reduction of variance in $x_{b}$ ⁠, second (a) and third (b) harmonics, acoustic branch, $Π_{4} = 2$ ⁠, $μ = (π / 4), Π_{1} = 0$

Fig. 15

View large Download slide

Reduction of variance in $x_{b}$ ⁠, second (a) and third (b) harmonics, optical branch, $Π_{4} = 2$ ⁠, $μ = (π / 4)$ ⁠, $Π_{1} = 0$

Figures 11 and 12 also show that an increase in the magnitude of $Π_{2}$ leads to a large variance reduction in the third harmonic, and conversely, an increase in $Π_{3}$ leads to a large variance reduction of the second harmonic. This counter-intuitive crossing effect can be understood by close examination of the second-order results in Eqs. (39)–(42) for the monoatomic chain (the same argument can be made for the diatomic chain). There, $k_{2}$ influences the particular solution governing the third harmonic, and $k_{3}$ influences the particular solution governing the second harmonic. Note that, the first-order solution does not contain such a crossing effect and thus injected waveforms based solely on it will not be sufficient for preventing an exchange of energy between different harmonics. In contrast, because the presence of the crossing terms in second-order injected waves keeps spectral energy, where it is initially assigned, the crossing effect becomes necessary for increasing invariance. This is yet another example of the richness in system behavior uncovered by higher-order perturbation analysis.

Conclusions

Higher-order multiple scales solutions have been detailed for predicting multiharmonic plane wave propagation in cubically and quadratically nonlinear monoatomic and diatomic chains. The method yields higher-order amplitude-dependent dispersion relationships, the result of which clearly show the importance of cubic terms vis-à-vis quadratic terms in frequency/wavenumber shifts. Wave stability has also been assessed through introduction of linear damping and a local stability analysis, and the notion of wave invariance has been explored and demonstrated for waveforms composed using higher-order terms in the asymptotic expansions. Finally, trends in dispersion shifts, wave stability, and wave invariance have been verified using dimensionless parameters and direct numerical integration.

Acknowledgment

The authors would like to thank the National Science Foundation for support of this research under Grant No. CMMI 1332862.

Appendix

References

1.

Li

,

F.

,

Anzel

,

P.

,

Yang

,

J.

,

Kevrekidis

,

P. G.

, and

Daraio

,

C.

,

2014

, “

Granular Acoustic Switches and Logic Elements

,”

Nat. Commun.

,

5

, pp.

1

–

6

.

2.

Jensen

,

J. S.

,

2003

, “

Phononic Band Gaps and Vibrations in One- and Two-Dimensional Mass–Spring Structures

,”

J. Sound Vib.

,

266

(

5

), pp.

1053

–

1078

.

Google Scholar

Crossref

3.

Sun

,

H.-X.

,

Zhang

,

S.-Y.

, and

Shui

,

X.-J.

,

2012

, “

A Tunable Acoustic Diode Made by a Metal Plate With Periodical Structure

,”

Appl. Phys. Lett.

,

100

(

10

), p.

103507

.

Google Scholar

Crossref

4.

Vakakis

,

A. F.

, and

King

,

M. E.

,

1995

, “

Nonlinear Wave Transmission in a Monocoupled Elastic Periodic System

,”

J. Acoust. Soc. Am.

,

98

(

3

), pp.

1534

–

1546

.

Google Scholar

Crossref

5.

Chakraborty

,

G.

, and

Mallik

,

A. K.

,

2001

, “

Dynamics of a Weakly Non-Linear Periodic Chain

,”

Int. J. Nonlinear Mech.

,

36

(

2

), pp.

375

–

389

.

Google Scholar

Crossref

6.

Dubus

,

B.

,

Swinteck

,

N.

,

Muralidharan

,

K.

,

Vasseur

,

J. O.

, and

Deymier

,

P. A.

,

2016

, “

Nonlinear Phonon Modes in Second-Order Anharmonic Coupled Monoatomic Chains

,”

ASME J. Vib. Acoust.

,

138

(

4

), p.

041016

.

Google Scholar

Crossref

7.

Narisetti

,

R. K.

,

Leamy

,

M. J.

, and

Ruzzene

,

M.

,

2010

, “

A Perturbation Approach for Predicting Wave Propagation in One-Dimensional Nonlinear Periodic Structures

,”

ASME J. Vib. Acoust.

,

132

(

3

), p.

031001

.

Google Scholar

Crossref

8.

Manktelow

,

K. L.

,

Leamy

,

M. J.

, and

Ruzzene

,

M.

,

2014

, “

Weakly Nonlinear Wave Interactions in Multi-Degree of Freedom Periodic Structures

,”

Wave Motion

,

51

(

6

), pp.

886

–

904

.

Google Scholar

Crossref

9.

Manktelow

,

K.

,

Leamy

,

M. J.

, and

Ruzzene

,

M.

,

2010

, “

Multiple Scales Analysis of Wave-Wave Interactions in a Cubically Nonlinear Monoatomic Chain

,”

Nonlinear Dyn.

,

63

(

1

), pp.

193

–

203

.

10.

Narisetti

,

R. K.

,

2010

, “

Wave Propagation in Nonlinear Periodic Structures

,”

Ph.D. thesis

, Georgia Institute of Technology, Atlanta, GA.

11.

Hussein

,

M. I.

,

Leamy

,

M. J.

, and

Ruzzene

,

M.

,

2014

, “

Dynamics of Phononic Materials and Structures: Historical Origins, Recent Progress, and Future Outlook

,”

ASME Appl. Mech. Rev.

,

66

(

4

), p.

040802

.

Google Scholar

Crossref

12.

Cabaret

,

J.

,

Tournat

,

V.

, and

Béquin

,

P.

,

2012

, “

Amplitude-Dependent Phononic Processes in a Diatomic Granular Chain in the Weakly Nonlinear Regime

,”

Phys. Rev. E

,

86

(

4

), p.

041305

.

Google Scholar

Crossref

13.

Sánchez-Morcillo

,

V. J.

,

Pérez-Arjona

,

I.

,

Romero-García

,

V.

,

Tournat

,

V.

, and

Gusev

,

V. E.

,

2013

, “

Second-Harmonic Generation for Dispersive Elastic Waves in a Discrete Granular Chain

,”

Phys. Rev. E

,

88

(

4

), p.

043203

.

Google Scholar

Crossref

14.

Swinteck

,

N. Z.

,

Muralidharan

,

K.

, and

Deymier

,

P. A.

,

2013

, “

Phonon Scattering in One-Dimensional Anharmonic Crystals and Superlattices: Analytical and Numerical Study

,”

ASME J. Vib. Acoust.

,

135

(

4

), p.

041016

.

Google Scholar

Crossref

15.

Hasegawa

,

A.

, and

Tappert

,

F.

,

1973

, “

Transmission of Stationary Nonlinear Optical Pulses in Dispersive Dielectric Fibers—I: Anomalous Dispersion

,”

Appl. Phys. Lett.

,

23

(

3

), pp.

142

–

144

.

Google Scholar

Crossref

16.

Laitone

,

E.

,

1960

, “

The Second Approximation to Cnoidal and Solitary Waves

,”

J. Fluid Mech.

,

9

(

03

), pp.

430

–

444

.

Google Scholar

Crossref

17.

Friedrichs

,

K. O.

,

1948

, “

On the Derivation of the Shallow Water Theory

,”

Commun. Pure Appl. Math.

,

1

, pp.

81

–

85

.

18.

Grosvenor

,

G.

,

1933

, “

Flying

,”

Natl. Geogr. Mag.

,

63

(

5

), pp.

585

–

630

.

19.

Fenton

,

J. D.

,

1999

, “

The Cnoidal Theory of Water Waves

,”

Developments in Offshore Engineering

,

J. B. Herbich and K. A. Ansari K.

,

eds.

, Gulf Publisher Company, Houston, TX, pp. 2–32.

20.

Kraenkel

,

R. A.

,

Manna

,

M. A.

,

Merle

,

V.

,

Montero

,

J. C.

, and

Pereira

,

J. G.

,

1996

, “

Multiple-Time Higher-Order Perturbation Analysis of the Regularized Long-Wavelength Equation

,”

Phys. Rev. E.

,

54

(3), pp.

2976

–

2981

.

Google Scholar

Crossref

21.

Friesecke

,

G.

, and

Mikikits-Leitner

,

A.

,

2015

, “

Cnoidal Waves on Fermi-Pasta-Ulam Lattices

,”

J. Dyn. Differ. Equations

,

27

(

3

), pp.

627

–

652

.

Google Scholar

Crossref

22.

Gaison

,

J.

,

Moskow

,

S.

,

Wright

,

J. D.

, and

Zhang

,

Q.

,

2014

, “

Approximation of Polyatomic FPU Lattices by KdV Equations

,”

Multiscale Model. Simul.

,

12

(

3

), pp.

953

–

995

.

Google Scholar

Crossref

23.

Friesecke

,

G.

, and

Pego

,

R. L.

,

1999

, “

Solitary Waves on FPU Lattices: I. Qualitative Properties, Renormalization and Continuum Limit

,”

Nonlinearity

,

12

(

6

), pp.

1601

–

1627

.

Google Scholar

Crossref

24.

Sen

,

S.

,

Hong

,

J.

,

Bang

,

J.

,

Avalos

,

E.

, and

Doney

,

R.

,

2008

, “

Solitary Waves in the Granular Chain

,”

Phys. Rep.

,

462

(

2

), pp.

21

–

66

.

Google Scholar

Crossref

25.

Kartashov

,

Y. V.

,

Malomed

,

B. A.

, and

Torner

,

L.

,

2011

, “

Solitons in Nonlinear Lattices

,”

Rev. Mod. Phys.

,

83

(

1

), pp.

247

–

305

.

Google Scholar

Crossref

26.

Molerón

,

M.

,

Leonard

,

A.

, and

Daraio

,

C.

,

2014

, “

Solitary Waves in a Chain of Repelling Magnets

,”

J. Appl. Phys.

,

115

(

18

), p.

184901

.

Google Scholar

Crossref

27.

Daraio

,

C.

,

Nesterenko

,

V. F.

,

Herbold

,

E. B.

, and

Jin

,

S.

,

2005

, “

Strongly Nonlinear Waves in a Chain of Teflon Beads

,”

Phys. Rev. E

,

72

(

1

), p.

016603

.

Google Scholar

Crossref

28.

Newton

,

P. K.

, and

Keller

,

J. B.

,

1987

, “

Stability of Periodic Plane Waves

,”

SIAM J. Appl. Math.

,

47

(

5

), pp.

959

–

964

.

Google Scholar

Crossref

29.

Bickham

,

S. R.

,

Kiselev

,

S. A.

, and

Sievers

,

A. J.

,

1993

, “

Stationary and Moving Intrinsic Localized Modes in One-Dimensional Monatomic Lattices With Cubic and Quartic Anharmonicity

,”

Phys. Rev. B

,

47

(

21

), pp.

14206

–

14211

.

Google Scholar

Crossref

30.

Friesecke

,

G.

, and

Pego

,

R. L.

,

2004

, “

Solitary Waves on Fermi–Pasta–Ulam Lattices: III. Howland-Type Floquet Theory

,”

Nonlinearity

,

17

(

1

), pp.

207

–

227

.

Google Scholar

Crossref

31.

Flach

,

S.

, and

Gorbach

,

A. V.

,

2008

, “

Discrete Breathers—Advances in Theory and Applications

,”

Phys. Rep.

,

467

(1–3), pp.

1

–

116

.

Google Scholar

Crossref

32.

Huang

,

G.

, and

Hu

,

B.

,

1998

, “

Asymmetric Gap Soliton Modes in Diatomic Lattices With Cubic and Quartic Nonlinearity

,”

Phys. Rev. B

,

57

(

10

), pp.

5746

–

5757

.

Google Scholar

Crossref

33.

Gorbach

,

A. V.

, and

Johansson

,

M.

,

2003

, “

Discrete Gap Breathers in a Diatomic Klein-Gordon Chain: Stability and Mobility

,”

Phys. Rev. E

,

67

(

6

), p.

066608

.

Google Scholar

Crossref

34.

Nayfeh

,

A. H.

, and

Mook

,

D. T.

,

1979

,

Nonlinear Oscillations

,

Wiley

, New York.

35.

Hamilton

,

M. F.

, and

Blackstock

,

D. T.

,

1998

,

Nonlinear Acoustics

,

Academic Press

,

San Diego, CA

.

Higher-Order Dispersion, Stability, and Waveform Invariance in Nonlinear Monoatomic and Diatomic Systems

Introduction

System Descriptions

Analysis Approach

Method of Multiple Scales: Monoatomic Chain.

Method of Multiple Scales: Diatomic Chain.

Numerical Simulation

Results

Dispersion Shifts.

Stability Considerations.

Local Stability Analysis.

Numerical Stability Studies.

Invariance.

Conclusions

Acknowledgment

Appendix

References

Get Email Alerts

Cited By

ASME Journals

ASME Conference Proceedings

ASME eBooks

Resources

Opportunities

Higher-Order Dispersion, Stability, and Waveform Invariance in Nonlinear Monoatomic and Diatomic Systems

Introduction

System Descriptions

Analysis Approach

Method of Multiple Scales: Monoatomic Chain.

Method of Multiple Scales: Diatomic Chain.

Numerical Simulation

Results

Dispersion Shifts.

Stability Considerations.

Local Stability Analysis.

Numerical Stability Studies.

Invariance.

Conclusions

Acknowledgment

Appendix

References

Get Email Alerts

Cited By

Related Articles

Related Proceedings Papers

Related Chapters

ASME Journals

ASME Conference Proceedings

ASME eBooks

Resources

Opportunities

This Feature Is Available To Subscribers Only