Doubly-nonlinear waveguides with self-switching functionality selection capabilities

TL;DR

This paper explores how combined quadratic and cubic nonlinearities in doubly-nonlinear waveguides enable tunable elastic wave control, acting as internal switches for high-frequency response features.

Contribution

It introduces a novel framework analyzing the interplay of two nonlinearities in elastic waveguides, revealing new wave control strategies and design principles for tunable metamaterials.

Findings

Cubic nonlinearity can act as an internal switch for wave modes.

Interplay of nonlinearities modifies dynamic response signatures.

Demonstrated control of high-frequency features via nonlinearity tuning.

Abstract

In this article, we investigate the effects of the interplay between quadratic and cubic nonlinearities on the propagation of elastic waves in periodic waveguides. Through this framework, we unveil an array of wave control strategies that are intrinsically available in the response of doubly-nonlinear systems and we infer some basic design principles for tunable elastic metamaterials. The objective is to simultaneously account for two sources of nonlinearity that are responsible for distinct and complementary phenomena and whose effects are therefore typically discussed separately in the literature. Our study explicitly targets the intertwined effects that the two types of nonlinearity exert on each other, which modify the way in which their respective signatures are observed in the dynamic response. Through two illustrative examples we show how the dispersion correction caused by cubic nonlinearity can be used as an internal switch, or mode selector, capable of tuning on/off certain high-frequency response features that are generated through quadratic mechanisms.

Equations

$$\begin{array}{rl}{f}_1& =J_2\delta \\ f_2& =G_2\delta +{ϵ}_1{G}_4{\delta }^3\\ f_3& =K_2\delta +{ϵ}_2{K}_3{\delta }^2\end{array}$$

$$\begin{array}{rl}& \left[\begin{array}{cc}{m}_1& 0\\ 0& m_2\end{array}\right]\left\{\begin{array}{c}{\ddot{u}}_n\\ {\ddot{v}}_n\end{array}\right\}+\left[\begin{array}{cc}2K_2+J_2& -J_2\\ -J_2& J_2+G_2\end{array}\right]\left\{\begin{array}{c}{u}_n\\ v_n\end{array}\right\}+\left[\begin{array}{cc}-K_2& 0\\ 0& 0\end{array}\right]\left\{\begin{array}{c}{u}_{n-1}\\ v_{n-1}\end{array}\right\}+\left[\begin{array}{cc}-K_2& 0\\ 0& 0\end{array}\right]\left\{\begin{array}{c}{u}_{n+1}\\ v_{n+1}\end{array}\right\}\\ & +{ϵ}_1\left\{\begin{array}{c}0\\ G_4{v}_n^3\end{array}\right\}+{ϵ}_2\left\{\begin{array}{c}{K}_3\left[{\left(u_n-u_{n-1}\right)}^2-{\left(u_{n+1}-u_n\right)}^2\right]\\ 0\end{array}\right\}=\left\{\begin{array}{c}0\\ 0\end{array}\right\}\end{array}$$

$$\begin{array}{rl}{\mathbf{u}}_n& =\left\{\begin{array}{c}{u}_n(t)\\ v_n(t)\end{array}\right\}\\ & =\left\{\begin{array}{c}{u}_n^0({\theta }_n,s,\tau )+{ϵ}_1{u}_n^1({\theta }_n,s,\tau )+{ϵ}_2{u}_n^2({\theta }_n,s,\tau )\\ v_n^0({\theta }_n,s,\tau )+{ϵ}_1{v}_n^1({\theta }_n,s,\tau )+{ϵ}_2{v}_n^2({\theta }_n,s,\tau )\end{array}\right\}\end{array}$$

$$\begin{array}{rl}O(1)& \text{mathchar}58{\omega }^2\mathbf{M}\frac{{\partial }^2{\mathbf{u}}_n^0}{\partial {\theta }_n^2}+{\mathbf{K}}_1{\mathbf{u}}_n+{\mathbf{K}}_2{\mathbf{u}}_{n-1}^0+{\mathbf{K}}_2{\mathbf{u}}_{n+1}^0=0\end{array}$$

$$\begin{array}{rl}O({ϵ}_1)& \text{mathchar}58{\omega }^2\mathbf{M}\frac{{\partial }^2{\mathbf{u}}_n^1}{\partial {\theta }_n^2}+{\mathbf{K}}_1{\mathbf{u}}_n^1+{\mathbf{K}}_2{\mathbf{u}}_{n-1}^1+{\mathbf{K}}_2{\mathbf{u}}_{n+1}^1={\mathbf{f}}^1\end{array}$$

$$\begin{array}{rl}O({ϵ}_2)& \text{mathchar}58{\omega }^2\mathbf{M}\frac{{\partial }^2{\mathbf{u}}_n^2}{\partial {\theta }_n^2}+{\mathbf{K}}_1{\mathbf{u}}_n^2+{\mathbf{K}}_2{\mathbf{u}}_{n-1}^2+{\mathbf{K}}_2{\mathbf{u}}_{n+1}^2={\mathbf{f}}^2\end{array}$$

$$\begin{array}{rl}{\mathbf{f}}^1& =2\omega \mathbf{M}\frac{{\partial }^2{\mathbf{u}}_n^0}{\partial {\theta }_n\partial \tau }\\ & -\left\{\begin{array}{c}0\\ G_4{\left(v_n^0\right)}^3\end{array}\right\}+\left\{\begin{array}{c}{K}_2\left(\frac{\partial u_{n+1}^0}{\partial s}-\frac{\partial u_{n-1}^0}{\partial s}\right)\\ 0\end{array}\right\}\end{array}$$

$$\begin{array}{rl}{\mathbf{f}}^2& =-\left\{\begin{array}{c}{K}_3\left[{\left(u_n^0-u_{n-1}^0\right)}^2-{\left(u_{n+1}^0-u_n^0\right)}^2\right]\\ 0\end{array}\right\}\end{array}$$

$${\mathbf{u}}_n^0=A(s,\tau )\mathbf{\varphi }{e}^{i{\theta }_n}+A^{\ast }(s,\tau ){\mathbf{\varphi }}^{\ast }{e}^{-i{\theta }_n}$$

$${\mathbf{u}}_{n\pm 1}=A(s,\tau )\mathbf{\varphi }{e}^{i{\theta }_n}{e}^{\pm i\xi }+A^{\ast }(s,\tau ){\mathbf{\varphi }}^{\ast }{e}^{-i{\theta }_n}{e}^{\mp i\xi }$$

$$\left(-{\omega }^2\mathbf{M}+\mathbf{K}(\xi )\right)\mathbf{\varphi }=0$$

$${\omega }_{L/U}=\sqrt{\frac{b\mp \Delta (\xi )}{2a}}$$

$${\mathbf{\varphi }}_{L/U}=\left\{\begin{array}{c}\frac{-{\omega }_{L/H}^2{m}_2+J_2+G_2}{J_2}\\ 1\end{array}\right\}$$

$${\omega }^2\bar{m}\frac{{\partial }^2{z}_n}{\partial {\theta }_n^2}+\bar{k}{z}_n={\mathbf{\varphi }}^H{\mathbf{f}}^1$$

$$\frac{\partial A}{\partial \tau }+\lambda \frac{\partial A}{\partial s}+i\mu {\left|A\right|}^{2}A=0$$

$$A\left(s,\tau \right)=\alpha \left(s,\tau \right)e^{-i\beta \left(s,\tau \right)}$$

$$\begin{array}{rl}\frac{\partial \alpha }{\partial \tau }+\lambda \frac{\partial \alpha }{\partial s}& =0\\ \frac{\partial \beta }{\partial \tau }+\lambda \frac{\partial \beta }{\partial s}& =\mu {\alpha }^2\end{array}$$

$$\begin{array}{rl}\alpha & ={\alpha }_0\left(s-\lambda \tau \right)\\ \beta & ={\beta }_0\left(s-\lambda \tau \right)+{\beta }^{\ast }\end{array}$$

$${\mathbf{u}}_n^0=A_0\mathbf{\varphi }{e}^{i\left[\xi n-\left(\omega +{ϵ}_1\mu A_0^2\right)t\right]}+c.c.$$

$${\mathbf{u}}_n^0=A_0\mathbf{\varphi }{e}^{i\left[\left(\xi -{ϵ}_1\mu A_0^2/\lambda \right)n-\omega t\right]}+c.c.$$

$$\frac{\partial A}{\partial \tau }+{\lambda }_0\sin \left(\xi -{ϵ}_1\frac{\beta }{s}\right)\frac{\partial A}{\partial s}+i\mu {\left|A\right|}^{2}A=0$$

$$\begin{array}{rl}\frac{\partial \alpha }{\partial \tau }+{\lambda }_0\sin \left(\xi -{ϵ}_1\frac{\beta }{s}\right)\frac{\partial \alpha }{\partial s}& =0\\ \frac{\partial \beta }{\partial \tau }+{\lambda }_0\sin \left(\xi -{ϵ}_1\frac{\beta }{s}\right)\frac{\partial \beta }{\partial s}& =\mu {\alpha }^2\end{array}$$

$${\lambda }_0{C}_1\sin \tilde{\xi }=\mu A_0^2$$

$${\mathbf{u}}_n^0=A_0\mathbf{\varphi }{e}^{i\left(\tilde{\xi }n-\omega t\right)}+c.c.$$

$$\begin{array}{r}{f}_u^2=-8i{K}_3{A}_0^2{\varphi }_u^2\sin \tilde{\xi }{\sin }^2\frac{\tilde{\xi }}{2}{e}^{2i{\tilde{\theta }}_n}+c.c.\end{array}$$

$${\mathbf{u}}_n^{2f}=B_1\mathbf{\phi }{e}^{2i{\tilde{\theta }}_n}+c.c.$$

$$B_1\mathbf{\phi }={\left[-4{\omega }^2\mathbf{M}+\mathbf{K}(\tilde{\xi })\right]}^{-1}{\bar{\mathbf{f}}}^2$$

$${\mathbf{u}}_n^{2h}=i{B}_2\mathbf{\chi }{e}^{i{\psi }_n}+c.c.$$

$$B_2\left(s,\tau \right)=a\left(s,\tau \right)e^{-ib\left(s,\tau \right)}$$

$${\hat{\mathbf{u}}}_n^0=A\mathbf{\varphi }{e}^{i{\theta }_n}+A^{\ast }{\mathbf{\varphi }}^{\ast }{e}^{-i{\theta }_n}+i{ϵ}_2{B}_2\mathbf{\chi }{e}^{i{\psi }_n}-i{ϵ}_2{B}_2^{\ast }{\mathbf{\chi }}^{\ast }{e}^{-i{\psi }_n}$$

$$\begin{array}{rl}{f}_{cubic}^1& =-G_4{\left(v_n^0\right)}^3=-G_4{\left(A{\varphi }_v{e}^{i{\theta }_n}+A^{\ast }{\varphi }_v^{\ast }{e}^{-i{\theta }_n}+i{ϵ}_2{B}_2{\chi }_v{e}^{i{\psi }_n}-i{ϵ}_2{B}_2^{\ast }{\chi }_v^{\ast }{e}^{-i{\psi }_n}\right)}^3\\ & =-3G_4\left({\left|A\right|}^2{\left|{\varphi }_v\right|}^{2}A{\varphi }_v+2{ϵ}_2^2{\left|B_2\right|}^2{\left|{\chi }_v\right|}^{2}A{\varphi }_v\right)e^{i{\theta }_n}-3i{ϵ}_2{G}_4\left(2{\left|A\right|}^2{\left|{\varphi }_v\right|}^2{B}_2{\chi }_v+{ϵ}_2^2{\left|B_2\right|}^2{\left|{\chi }_v\right|}^2{B}_2{\chi }_v\right)e^{i{\psi }_n}+c.c.\end{array}$$

$$\begin{gathered} \frac{\partial A}{\partial \tau }+{\lambda }_0\sin \left(\xi -{ϵ}_1\frac{\beta }{s}\right)\frac{\partial A}{\partial s} \\ +i\left(\mu {\left|A\right|}^2+{ϵ}_2^2\eta {\left|B_2\right|}^2\right)A=0 \end{gathered}$$

Full text

I Introduction

The most distinctive and fundamental property of phononic crystals (and periodic structures at large) is their inherent ability to attenuate waves whose spectral content falls within frequency intervals known as phononic bandgaps Hussein et al. (2014). This property stems directly from the periodicity of the material and strictly requires spatial compatibility between the wavelength of the excitation and the characteristic size of the unit cell. By incorporating internal resonators in the unit cells, bandgaps with subwavelength characteristics can be opened at low frequencies Liu et al. (2000), with onsets at the natural frequencies of the resonators. These gaps survive the relaxation of periodicity and can therefore be observed in spatially randomized systems Celli and Gonella (2015). A medium with locally-resonant mechanisms displays negative effective mass density Huang et al. (2009), which has also been exploited, in conjunction with other negative effective moduli, to achieve negative refraction Zhu et al. (2014). In addition to modifying the bandgap spectrum, the introduction of an auxiliary microstructure of resonators enriches the band structure of a crystal with new wave modes characterized by deformation mechanisms and dispersive properties that are germane to the resonators motion, thus stretching its functionality landscape.

While the bulk of the literature on phononic crystals has been concerned with the linear regime, recently, nonlinearity has received considerable attention as a mean to endow materials with unprecedented tunability. The main manifestation of nonlinearity is the amplitude dependence of the dynamic response, which offers additional opportunities to tune the spectro-spatial characteristics of a periodic structure by simply controlling the amplitude of excitation.

The weakly nonlinear regime has been studied with special focus on cubic and quadratic nonlinearity, due to their tractability and to the fact that they represent the dominant terms in many arbitrary nonlinear material behaviors. The main manifestation of cubic nonlinearity is a frequency correction of the dispersion bands, which can be predicted by perturbation methods, as reported in classical applied mathematics textbooks Holmes (2012). Narisetti et al. Narisetti et al. (2010, 2011) revisited this concept and updated the perturbation technique in the context of periodic structures to show amplitude dependence of the dispersion branches of phononic media, while experimental evidence of the frequency shifts achievable in granular chains was reported in Cabaret et al. (2012); Bonanomi et al. (2015). Other interesting cubic nonlinear phenomena, not strictly limited to weak nonlinearity, include intrinsic localized modes Sievers and Takeno (1988), chaotic bands Fang et al. (2017), and third harmonic generation (THG) Fang et al. (2018). In addition, Moore et al. Moore et al. (2018) achieved nonreciprocity in coupled oscillators featuring strong cubic nonlinearity realized using thin acrylic strips under transverse deformation.

In the case of quadratic nonlinearity, the main nonlinear effect of interest is the second harmonic generation (SHG) Nayfeh and Mook (1995); Zarembo and KrasilNikov (1971); Hamilton et al. (1998), while the correction of the dispersion relation has been shown to be negligible Ganesh and Gonella (2013). Along these lines, perturbation analysis has been used to find the explicit solution of second harmonic terms in monoatomic granular chains, and to explore the interplay of SHG and dispersion Sánchez-Morcillo et al. (2013). Recently, Mehrem et al. Mehrem et al. (2017) have conducted a series of experiments to confirm the existence of SHG in a nonlinear chain consisting of repulsive magnets. A few recent works have addressed the problem of nonlinear wave propagation in periodic structures with multiple degrees of freedom, where the availability of complex dispersive characteristics results in the nonlinear enrichment of the default linear response and in the development of new wave manipulation capabilities. Among these effects, we recall modal mixing in granular systems Ganesh and Gonella (2015); Wallen and Boechler (2017) and augmented spatial directivity in $2$D metamaterial lattices Ganesh and Gonella (2017). Jiao and Gonella Jiao and Gonella (2018a) incorporated locally resonant mechanisms in nonlinear waveguides to realize via SHG inter-modal energy tunneling from a fundamental flexural mode to a nonlinearly-activated axial mode with subwavelength characteristics, and subsequently confirmed this phenomenon experimentally Jiao and Gonella (2018b). Another example of this quadratic nonlinear tunneling is an acoustic rectifier realized by coupling a superlattice with a nonlinear microbubble suspension, in which a rectified energy flux of acoustic waves can be nonlinearly generated Liang et al. (2010). Moreover, phase matching condition (PMC) can be established when the fundamental harmonic and the second harmonic have the same phase velocity, resulting in cumulative second harmonic propagation Konotop (1996); De Lima and Hamilton (2003). This unique feature significantly increases the strength of SHG, thus providing a practical way for experimental verification Deng et al. (2005); Bermes et al. (2007) and applications of nonlinear effects Jiao and Gonella (2018b); Deng and Pei (2007).

While, in most relevant contexts, the effects of quadratic and cubic nonlinearity have been studied separately, a few works have dealt with these two types of nonlinearity jointly to study wave propagation in classical anharmonic lattices Huang et al. (1993); Askar and Lighthill (1973); Huang and Hu (1998) and in granular systems obeying Hertz contact law Cabaret et al. (2012); Wallen et al. (2017). However, little attention has been paid to the implications of the quadratic-cubic (Q-C) interplay for wave control and tunability. The scope of this article is precisely to fill this gap. To this end, we consider a medium - here a spring-mass chain on nonlinear elastic foundation - that features simultaneously both sources of nonlinearity (at different orders) embedded in the chain and in the support elastic elements. The coexistence of Q-C effects opens unprecedented opportunities for amplitude-controlled wave tuning. We are interested in how the cubic nonlinearity can be used to effectively modify the way in which the quadratic nonlinear effects manifest in the response. Specifically, we investigate whether the dispersion correction caused by cubic nonlinearity can be used as an internal switch (or mode selector) capable of switching on/off certain high-frequency features of the response that are generated by the quadratic mechanisms, or at least controlling their strength. This additional tuning capability is fully intrinsic to the system, as the selection of the functionality is solely due to the nonlinear effects, and therefore controlled by the amplitude of excitation, without resorting to any active external tuning.

In section 2, a multiple scales analysis is applied to obtain the full nonlinear response of the chain. Due to the Q-C interaction, SHG is observed with some unique modal characteristics, including magnified dispersion shift and tunable PMCs. In section 3, we proceed to verify our analytical model by performing a suite of numerical simulations with absorbing boundary layers. Finally, in section 4, two examples of tunable systems with self-adaptive switches are given to illustrate the wealth of wave manipulation capabilities that directly leverage the Q-C nonlinear interaction. The highlights of this investigation are summarized in the concluding section 5.

II Multiple scales analysis

II.1 Governing equations and linear dispersion relation

Our reference system is the spring-mass chain with a locally-resonant nonlinear foundation shown in Fig. 1.

In each cell, the primary mass $m_{1}$ is connected to an auxiliary mass $m_{2}$ through linear spring $1$, and $m_{2}$ is connected to a fixed point via spring $2$, which is assumed to feature cubic nonlinearity, thus implementing a locally-resonant nonlinear foundation. Spring $3$, connecting $m_{1}$ to its neighbors along the chain, is assumed to feature quadratic nonlinearity. Under weak nonlinearity assumptions, the force-stretch relations in the three springs can be expressed as:

[TABLE]

where $\epsilon_{1}$ and $\epsilon_{2}$ are small parameters for the cubic and quadratic nonlinearities, respectively. In this analysis, we assume that $\epsilon_{1}>\epsilon_{2}$ in order to establish the hierarchy between nonlinear terms necessary to achieve the desired tuning effects.

The equations of motion for the two masses ($m_{1}$ and $m_{2}$) in $n^{th}$ unit cell can be written in matrix form as:

[TABLE]

In the spirit of multiple scales analysis, we introduce a fast spatio-temporal variable $\theta_{n}=\xi n-\omega t$, where $\xi$ and $\omega$ are the normalized wavenumber and frequency, respectively, and $n$ is an integer mass index, as well as the two slow variables $s=\epsilon_{1}n$ (spatial), and $\tau=\epsilon_{1}t$ (temporal). Accordingly, the cell nodal response is expressed as an expansion, up to $O(\epsilon_{2})$, in successive approximations as:

[TABLE]

Substituting this into Eq. 2 and using the chain rule to resolve the time derivatives, one obtains the following system of cascading linear equations at each order of expansion:

[TABLE]

where $\mathbf{M}=\begin{bmatrix}m_{1}&0\\ 0&m_{2}\end{bmatrix}$; $\mathbf{K}_{1}=\begin{bmatrix}2K_{2}+J_{2}&-J_{2}\\ -J_{2}&J_{2}+G_{2}\end{bmatrix}$; $\mathbf{K}_{2}=\begin{bmatrix}-K_{2}&0\\ 0&0\end{bmatrix}$ and the forcing terms at $O(\epsilon_{1})$ and $O(\epsilon_{2})$ are:

[TABLE]

The general solution at $O(1)$ can be expressed as

[TABLE]

where $A$ is an amplitude term that only depends on the slow scale variables $(s,\tau)$, $\boldsymbol{\phi}=\begin{Bmatrix}\phi_{u}\\ \phi_{v}\\ \end{Bmatrix}$ is a modal vector, and $(\cdot)^{*}$ denotes the complex conjugate of a variable. Imposing Bloch conditions on the fast scale variable $\theta_{n}$ between neighboring cells, we obtain

[TABLE]

Substituting Eq. 9 and Eq. 10 into Eq. 4, the linear dispersion relation is obtained by solving the eigenvalue problem

[TABLE]

where $\mathbf{K(\xi)}=\begin{bmatrix}2K_{2}\left(1-\cos\xi\right)+J_{2}&-J_{2}\\ -J_{2}&J_{2}+G_{2}\end{bmatrix}$. Two branches, denoted as “L” for the lower branch and “U” for the upper branch, are obtained analytically as:

[TABLE]

where $\Delta(\xi)=\sqrt{b^{2}-4ac}$, with $a=m_{1}m_{2}$, $b=\left(m_{1}+m_{2}\right)J_{2}+2m_{2}K_{2}\left(1-\cos\xi\right)+m_{1}G_{2}$, and

$c=2K_{2}\left(J_{2}+G_{2}\right)\left(1-\cos\xi\right)+J_{2}G_{2}.$ The amplitude ratio for each mode is given by

[TABLE]

Fig. 2 depicts the dispersion relation of the linearized model obtained for an arbitrary choice of parameters: $m_{1}=1$, $m_{2}=0.5$, $K_{2}=J_{2}=G_{2}=1$ (standard SI units are adopted throughout the paper and omitted for simplicity). Visual inspection of the plot reveals a few peculiar dispersive characteristics of the chain. (1) The system is fully gapped at low frequencies, between [math] and the cut-on frequency for the first branch $\omega^{1}_{\mathrm{{cut-on}}}$, as expected for a chain on an elastic foundation. (2) A locally-resonant bandgap is opened in the interval between $\omega^{1}_{\mathrm{{cut-off}}}$ and $\omega^{2}_{\mathrm{{cut-on}}}$. The positive group velocity exhibited by the upper branch at $\xi=0$ is another clear signature of a locally-resonant mechanism.

II.2 Nonlinear correction of the dispersion relation

We now proceed to solve the equations at $O(\epsilon_{1})$. Following standard modal decomposition Meirovitch (2010), we introduce the modal coordinates $\mathbf{z}_{n}$, such that $\mathbf{u}^{1}_{n}=\boldsymbol{\Phi}\mathbf{z}_{n}$, where $\boldsymbol{\Phi}$ is the modal matrix (featuring the modal vectors $\boldsymbol{\phi}$ on the columns). Pre-multiplying Eq. 5 by $\boldsymbol{\phi}^{H}$, where $(\cdot)^{H}$ denotes conjugate transpose, yields

[TABLE]

where $\bar{m}=\boldsymbol{\phi}^{H}\mathbf{M}\boldsymbol{\phi}$ and $\bar{k}=\boldsymbol{\phi}^{H}\mathbf{K}\boldsymbol{\phi}$ are the modal mass and stiffness matrices. To prevent unbounded solution, the secular terms in $\boldsymbol{\phi}^{H}\mathbf{f}^{1}$ need to be eliminated Manktelow et al. (2014), which, upon some algebraic steps, leads to the required condition

[TABLE]

where $\lambda=\frac{K_{2}\mathinner{\!\left\lvert\phi_{u}\right\rvert}^{2}\sin\xi}{\omega\bar{m}}$ and $\mu=\frac{3G_{4}\mathinner{\!\left\lvert\phi_{v}\right\rvert}^{4}}{2\omega\bar{m}}$. Note that Eq. 15 is a first-order PDE that controls the evolution of $A$ with respect to slow space and time variables. In other words, eliminating secularity from the solution results in a constraint on the amplitude of the fundamental harmonic. The amplitude $A$ can be written in polar form as

[TABLE]

Substituting it into Eq. 15, yields a complex-variable algebraic equation, the soution of which requires that real and imaginary parts vanish individually, yielding the two conditions:

[TABLE]

The general solutions for $\alpha$ and $\beta$ are

[TABLE]

The full expression for $\beta$ contains a particular solution $\beta^{*}$, which can be expressed either in terms of variable $\tau$ as $\beta^{*}=\mu\alpha^{2}\tau$, or in terms of variable $s$ as $\beta^{*}=\mu\alpha^{2}s/\lambda$. Accordingly, in order to determine $\beta^{*}$, either initial conditions or boundary conditions must be supplied, based on the problem of interest. For example, by imposing the initial harmonic amplitude profile $\mathbf{u}_{n}=A_{0}\boldsymbol{\phi}\sin\xi n$ at $t=0$, it follows that $\alpha=A_{0}$, $\beta_{0}=0$ and $\beta^{*}=\mu\alpha^{2}\tau$. Thus, the fundamental solution at $O(1)$ is

[TABLE]

where $c.c.$ denotes the complex conjugate of all preceding terms. Eq. 19 can be seen as a harmonic wave with wavenumber $\xi$ and frequency shifted by the correction term $\epsilon_{1}\mu A^{2}_{0}$. In contrast, if the harmonic displacement time history $\mathbf{u}_{n}=A_{0}\boldsymbol{\phi}\sin\omega t$ is imposed at one end of the chain (e.g., at $n=0$), it follows that $\alpha=A_{0}$, $\beta_{0}=0$ and $\beta^{*}=\mu\alpha^{2}s/\lambda$. With this, the fundamental solution reads

[TABLE]

in which the term $\epsilon_{1}\mu A^{2}_{0}/\lambda$ is a correction shift for the wavenumber.

From Eq. 19 and Eq. 20, it is clear that that a shift can be induced either in the frequency domain or in the wavenumber domain, depending on whether initial conditions or boundary conditions are imposed. Moreover, the shifts are proportional to $A^{2}_{0}$, which means that the dispersion relations are amplitude-dependent. While the two formulations appear interchangeable, the manifestations of the two shifts and their effects on the dispersion relation tuning are profoundly different. In Fig. 3, we plot the corrected dispersion branches predicted by our analysis with the following parameters: $m_{1}=1$, $m_{2}=0.5$, $K_{2}=J_{2}=G_{2}=G_{4}=1$, $\epsilon_{1}=0.1$, $A_{0}=1$. The curves calculated from Eq. 19 and Eq. 20 coincide over a large range of wavenumbers, but significant deviations are observed when $\xi$ approaches [math] or $\pi$. In the latter limit, the dispersion relation predicted by Eq. 20 diverges and loses physical meaning (a similar issue is reported in Chakraborty and Mallik (2001) where perturbation expansion is applied for the wavenumber), while Eq. 19 shows a robust estimation on the frequency shift over the whole Brillouin zone.

To resolve this discrepancy, we have to realize that the Bloch condition (Eq. 10), which is written in terms of wavenumber, needs to be updated once a wavenumber shift is determined, which means that the original wavenumber $\xi$ should be replaced by the modified one $\xi-\epsilon_{1}\beta/s$ (this expression can be obtained by plugging Eq. 16 in Eq. 9 and combining the phase terms). The updated Bloch condition leads to a new equation (replacing Eq. 15) for the elimination of the secular terms at $O(\epsilon_{1})$:

[TABLE]

where $\lambda_{0}=\frac{K_{2}\mathinner{\!\left\lvert\phi_{u}\right\rvert}^{2}}{\omega\bar{m}}$. Substituting Eq. 16 in Fig. 21 leads to the following equations for $\alpha$ and $\beta$

[TABLE]

If the wavenumber shift $\epsilon_{1}\beta/s$ is small when compared to $\xi$, the above system of equations can be reduced to Eq. 17 to the first order approximation. In all other cases, we need to solve Eq. 22 to find the solution of $\beta$. The particular solution of $\beta$ can be assumed as $\beta^{*}=C_{1}s=C_{1}\epsilon_{1}n$ if plane wave solution is allowed, where $C_{1}$ is a constant. Substituting it into Eq. 22 and noticing that $\alpha=A_{0}$ under boundary excitation with constant amplitude $A_{0}$, yields the transcendental equation for $C_{1}$

[TABLE]

where $\tilde{\xi}=\xi-\epsilon_{1}C_{1}$. Once the root of Eq. 23 has been numerically determined, the corrected dispersion relation can be computed using the updated wavenumber $\tilde{\xi}$, resulting in the curves of Fig. 4. It can be noticed that the spurious dispersive behavior observed in Eq. 3 (d) for large values of $\xi$ is successfully removed.

This result leads to a few interesting observations. When a correction in $\xi$ is enforced (as a result of an imposed boundary condition), the corrected dispersion curves are not defined over the entire reduced Brillouin zone. This is in contrast with the frequency correction Narisetti et al. (2010), which is defined for every $\xi$. This implies that, at certain wavenumbers (especially, approaching the [math] and $\pi$ limits), wave propagation is forbidden in the nonlinear system if a boundary condition is imposed. The reasonableness and physical implications of this result will be substantiated via numerical simulations in the following chapters.

In conclusion, regardless of whether boundary conditions or initial conditions are considered, cubic nonlinearity results in a tangible modification of the dispersion relation. Since imposing boundary conditions is vastly more practical than prescribing initial conditions in most problems of engineering relevance, we focus on this case in the remainder of the paper. Then, the fundamental solution with the updated wavenumber shift $\tilde{\xi}$ can be expressed compactly as

[TABLE]

.

II.3 Second harmonic generation and its correction shift

It has been shown that, in the presence of quadratic nonlinearity, the complete solution at $O(\epsilon_{2})$ encompasses a non-oscillatory term and an oscillatory contribution, the latter comprising the solution of the corresponding homogeneous problem and the particular (forced) solution due to the forcing term $\mathbf{f}^{2}$ Sánchez-Morcillo et al. (2013). The question we address here is whether the presence of cubic nonlinear springs in the foundation introduces additional tuning effects that are appreciated at the level of these second harmonic contributions. To this end, in the following we derive a general expression for the second harmonic solution that explicitly incorporates any spectral modulation due to such nonlinear interplay.

If we substitute the associated fundamental solution (i.e., Eq. 24) in Eq. 5, the expression for the forcing function (note that only the first component $f_{u}^{2}$ of $\mathbf{f}^{2}$ is non-zero) can be determined as (details about the derivation are given in Appendix)

[TABLE]

where $\tilde{\theta}_{n}=\tilde{\xi}n-\omega t$ and $\tilde{\xi}$ is the corrected wavenumber introduced above and, by working with $\tilde{\xi}$, we implicitly assume that the excitation is prescribed as a boundary condition. Since the coefficient in front of $e^{2i\tilde{\theta}_{n}}$ is purely imaginary, this expression suggests a phase shift in the second harmonic.

II.3.1 Forced solution and its correction

Following the standard procedure for heterogeneous PDEs, the forced solution can be found as

[TABLE]

where the coefficient can be obtained by plugging Eq. 26 in Eq. 5, yielding

[TABLE]

where $\mathbf{\bar{f}}^{2}=\begin{Bmatrix}-8iK_{3}A^{2}_{0}\phi_{u}^{2}\sin\tilde{\xi}\sin^{2}\frac{\tilde{\xi}}{2}\\ 0\ \end{Bmatrix}$ and $\bm{\varphi}$ is subjected to the same normalization followed for the fundamental modal vector $\bm{\phi}$. The phase $2i\tilde{\theta}_{n}$ indicates that $\mathbf{u}_{n}^{2f}$ represents a harmonic wave traveling with $2\tilde{\xi}$ and $2\omega$, i.e., twice the frequency and twice the (corrected) wavenumber of the fundamental harmonic $\mathbf{u}_{n}^{0}$. While the doubling of frequency and wavenumber is an expected feature of the forced response in systems with quadratic nonlinearity Sánchez-Morcillo et al. (2013); Mehrem et al. (2017), the dependence on $\tilde{\xi}$ introduces a new spectral control that is germane to the problem with additional cubic nonlinearity. In essence, cubic nonlinearity introduces an additional wavenumber shift in the signature of the quadratically-generated forced second harmonic. We conclude that the interplay between the two types of nonlinearity results in a new amplitude-dependent correction mechanism that endows traveling waves with an additional self-manipulation capability. We now proceed to check whether a similar effect is observed in the homogeneous part of the SHG.

II.3.2 Homogeneous solution and its correction

The homogeneous solution can be expressed as

[TABLE]

In contrast with the forced response, here the phase $\psi=\xi(2\omega)n-2\omega t$ involves the frequency-wavenumber pair predicted by the branch of the band diagram available at $2\omega$, thus conforming to the dispersion relation of the corresponding liner problem. Accordingly, $\bm{\chi}=\begin{Bmatrix}\chi_{u}\\ \chi_{v}\\ \end{Bmatrix}$ is the modal vector at $\left(\xi(2\omega),2\omega\right)$, and the imaginary constant $i$ is assume a priori in Eq. 28 is to simplify the analysis.

In the spirit of multiple scale analysis, the amplitude term $B_{2}$ is, in general, taken to be dependent upon the slow variables $s$ and $\tau$, and can be expressed as

[TABLE]

The full characterization of the homogeneous response relies on the non-trivial task of determining $B_{2}\left(s,\tau\right)$. To this end, we introduce a manipulation step aimed at virtually extending to the second harmonic the treatment, classically reserved for fundamental harmonic, that we followed in section 1.1 for the determination of $A(s,\tau)$. Specifically, we formally combine the homogeneous second harmonic and the fundamental harmonic into a new baseline expression, which we will treat as the “new fundamental solution” at $O(1)$:

[TABLE]

Accordingly, the equation at $O(\epsilon_{1})$ is also updated, and the cubic nonlinear term in $\mathbf{f}^{1}$ becomes

[TABLE]

where higher order terms have been omitted. Following the same logic described in section $1.1$, the removal of the secular terms for the wave modes associated with $e^{i\theta_{n}}$ and $e^{i\psi_{n}}$ requires that

[TABLE]

and

[TABLE]

where $\eta=\frac{3G_{4}\mathinner{\!\left\lvert\chi_{v}\right\rvert}^{2}\mathinner{\!\left\lvert\phi_{v}\right\rvert}^{2}}{\omega\bar{m}}$; $\lambda^{\prime}_{0}=\frac{K_{2}\mathinner{\!\left\lvert\chi_{u}\right\rvert}^{2}}{2\omega\bar{m}}$; $\mu^{\prime}=\frac{3G_{4}\mathinner{\!\left\lvert\phi_{v}\right\rvert}^{2}\mathinner{\!\left\lvert\chi_{v}\right\rvert}^{2}}{2\omega\bar{m}}$; $\eta^{\prime}=\frac{3G_{4}\mathinner{\!\left\lvert\chi_{v}\right\rvert}^{4}}{4\omega\bar{m}}$. Ignoring the high order term $O(\epsilon^{2}_{2})$, Eq. 32 reduces back to Eq. 21, whose solution yields the dispersion shift in the fundamental harmonic, as discussed before. Similarly, neglecting high order term in Eq. 33, the equation becomes

[TABLE]

which is a partial differential equation in $s$ and $\tau$ controlling the slow spatio-temporal evolution of $B_{2}$. Substituting Eq. 29 in Eq. 34 and assuming $b=C_{2}s=\epsilon_{1}C_{2}n$ (an apriori assumption which leads to a traveling wave solution if $C_{2}$ can be solved as a real number), yields a transcendental equation formally similar to Eq. 23, which can be solved numerically for $C_{2}$. With $C_{2}$, the homogeneous solution can be rewritten as

[TABLE]

where $\hat{\psi}_{n}=\hat{\xi}n-2\omega t$ and $\hat{\xi}=\xi(2\omega)-\epsilon_{1}C_{2}$. This result shows that the Q-C interplay leads to a new dispersion shift that is germane to the homogeneous second harmonic.

II.3.3 Phase matching conditions

A careful inspection on Eq. 27 suggests that the above procedures for finding the solutions at $O(\epsilon_{2})$ works only when $-4\omega^{2}\mathbf{M}+\mathbf{K}\equiv\mathbf{Q}$ is invertible. It is possible, however, that $\mathbf{Q}$ is singular ($det(\mathbf{Q})=0$). If that is the case, the term containing the variable $2\tilde{\theta}$ in the forcing function $\mathbf{f}^{2}$ satisfies the linear eigenvalue problem of Eq. 11. In other words, the forcing function resonates with the homogeneous solution (i.e., $2\tilde{\theta}_{n}=\psi_{n}$) and the solution becomes secular. Such scenario is referred to as phase matching condition (PMC) and in general has been shown to be conducive to stronger nonlinear response signatures Jiao and Gonella (2018a, b). The solution, which is assumed to grow in direct proportion to the propagation distance, is expressed as

[TABLE]

Substituting it in Eq. 5, yields

[TABLE]

where $\mathbf{p}=-K_{2}\begin{Bmatrix}2i\chi_{u}\sin 2\tilde{\xi}\\ 0\\ \end{Bmatrix}$. The first term in the LHS of Eq. 37 automatically vanishes since it satisfies the eigenvalue problem Eq. 11, as stated above. Then, the coefficient $D$ can be determined as: $D=\frac{4K_{3}A^{2}_{0}\phi^{2}_{u}\sin\tilde{\xi}\sin^{2}\frac{\tilde{\xi}}{2}}{K_{2}\chi_{u}\sin 2\tilde{\xi}}$.

As we mentioned earlier in the manuscript, this type of nonlinear system also supports a non-oscillatory solution - a general nonlinear effect in anharmonic lattices Mehrem et al. (2017). For completeness, let us briefly re-examine the forcing term at $O(\epsilon_{2})$ (see Appendix), in which we notice that a constant force can be induced in each quadratic spring by the underlined terms in Eq. 42.

[TABLE]

This constant force results in a stationary displacement $\mathbf{u}^{2s}_{n}$ in the system, and such static problem can be easily solved once the size (i.e., the total number of degrees of freedom) of the system is specified. In conclusion, the complete solution at $O(\epsilon_{2})$ can be expressed as

If $2\tilde{\theta}_{n}\neq\psi_{n}$ (Non-PMC):

[TABLE]

If $2\tilde{\theta}_{n}=\psi_{n}$ (PMC):

[TABLE]

Imposing the absence of second harmonic and of non-oscillatory terms at the boundary ($n=0$), it yields $b_{0}=0$ and $a=iB_{1}\varphi_{v}/\chi_{v}$, and the unknowns in $\mathbf{u}^{2s}_{n}$ can be fully determined.

Equations 39 and 40 shed light over some unique characteristics of the nonlinearly-generated second harmonics, which are summarized below:

1. The shift in wavenumber introduced by the cubic nonlinearity at $O(\epsilon_{1})$, is also observed in the forced solution at the second harmonic. However, the shift at the second harmonic is magnified by a factor of $2$ (i.e., 2$\Delta\xi=2(\tilde{\xi}-\xi)$), which implies that the spectral signatures of the forced second harmonic can be modified to a larger extent than the fundamental harmonic by playing with the amplitude of excitation.

2. We have shown that the interaction between the quadratic and cubic nonlinearities produces a new dispersion shift in the second harmonics, which is conceptually distinct from the shift in the fundamental harmonic. In other words, at a given frequency, two shifts are defined: the one experienced by a fundamental harmonic excited directly at that frequency and the one experienced by the second harmonic generated by an excitation with half frequency content. We note that the dispersion shift in the homogeneous second harmonic (dashed red lines) is significantly larger (precisely twice larger) when compared to the shift in the fundamental solution (red lines) at that frequency. A visual comparison of these shifts is provided in Fig. 5 for different frequencies and strengths of nonlinearity. This unique feature enables an additional array of wave manipulation opportunities, as will be discussed in the application section below.

3. When PMC is achieved, the amplitude of the second harmonic grows with the propagation distance, and the full nonlinear response conforms to the linear dispersion relation. As discussed in Jiao and Gonella (2018a, b), this unique property has important practical implications. The energy tunneling induced by the SHG from the fundamental to second harmonic can be boosted by enabling PMC.

4. In the presented multiple scales analysis, the selection of the scale variables is deliberately chosen to capture the main effects of quadratic and cubic nonlinearities in the system. A more general choice, up to $O(\epsilon_{2})$ would include: $\theta_{n}=\xi n-\omega t$; $s_{1}=\epsilon_{1}n$; $\tau_{1}=\epsilon_{1}t$; $s_{2}=\epsilon_{2}n$; $\tau_{2}=\epsilon_{2}t$. However, the two additional variables $s_{2}$ and $\tau_{2}$ will not add any new information to the existing solutions in our context. This can be explained by recalling that the main purpose of introducing the slow variables $s$ and $\tau$ is to capture the shift in the dispersion relation. To the first order approximation (i.e., at $O(\epsilon_{1})$), this is brought about exclusively by the cubic nonlinearity, while the shift due to quadratic nonlinearity is negligible Ganesh and Gonella (2013). Another consequence of this argument is that the presented perturbation analysis also works when $\epsilon_{2}$ is comparable to $\epsilon_{1}$, i.e., $\epsilon_{1}\approx\epsilon_{2}$. Under these circumstances, Eq. 5 and Eq. 6 can be combined into one equation with the two forcing functions $\mathbf{f}^{1}$ and $\mathbf{f}^{2}$, and the rest of the procedure remains the same.

III Full-scale simulations

In this section, we perform full-scale numerical simulations to validate our analytical model and support the findings detailed in section 2. The wave response of a finite chain is obtained by integrating the system of nonlinear equations (Eq. 2) using the Verlet Algorithm Swope et al. (1982). In order to obtain clean signatures of the desired harmonics, it is imperative to work with traveling waves that are as monochromatic as possible. This can be achieved by prescribing sustained harmonic excitations at the boundary. However, this route is impractical with finite-size systems as reflections from the boundaries would soon result in standing waves (i.e., vibration states).

To prescribe purely harmonic traveling waves while avoiding boundary reflections, we implement an absorbing layer of $30$ unit cells endowed with viscous damping at the end of the main chain which consists of $40$ undamped unit cells. The coefficient profile for the $n^{th}$ unit cell in the absorbing layer is given by

[TABLE]

where $C_{0}=3$ is the maximum damping coefficient, $N_{\mathrm{PML}}=30$ denotes the number of unit cells in the layer, $\gamma$ is set to $3$ in the simulation, and $\mathbf{I}$ is the $2\times 2$ identity matrix. A continuous harmonic excitation is applied to mass $m_{2}$ in the first unit cell of the chain. The spatio-temporal response computed at each cell site over the considered time interval is transformed to the Fourier domain via two-dimensional Discrete Fourier Transform (2D-DFT). The spectral components $(\xi,\omega)$ of the response, to be compared against the band diagrams from unit cell analysis, are extracted by tracking the maximum values of the 2D-DFT data in different frequency intervals.

III.1 Numerical estimation of dispersion shifts in the fundamental harmonic

The system parameters are $m_{1}=1$, $m_{2}=0.5$, $K_{2}=J_{2}=G_{2}=1$, $\epsilon_{1}=0.1$. The cubic nonlinearity can be positive ($G_{4}=1$) or negative ($G_{4}=-1$), which corresponds to a hardening or softening nonlinear foundation, respectively. The amplitude of the excitation $A_{0}$ is set to $1$ for the lower branch correction (Fig. 6(a)) and to $0.6$ for the upper branch correction (Fig. 6(b)).

Fig. 6 shows that the numerical results match the multiple scales analysis predictions. The branches are shifted upward or downward depending on whether hardening or softening is considered. An interesting observation can be made about the regime close to the cutoff frequencies, where the result obtained by Eq. 20 (Fig. 3) diverged spuriously. In these regimes, the responses from simulations tend to be unstable and monochromatic conditions cannot be established effectively. This can be qualitatively explained by invoking the fact that the dispersion curves become flat around the cutoff frequencies, resulting in large wavenumber shifts (if boundary conditions are imposed) that violate the weak nonlinearity assumption. As a result, the dispersion shift predicted by Eq. 20 can be used as a guideline to check validity of the assumptions underpinning our perturbation analysis for a given system and excitation conditions.

III.2 Numerical estimation of dispersion shift in SHG

So far, we have numerically matched only the shifts described by Eq. 24, which only affect the fundamental harmonic. Now, we incorporate quadratic nonlinearity ($K_{3}=1$ and $\epsilon_{2}=0.05$) in the springs of the main chain to elicit second harmonic generation and estimate its spectral shift. The spectral content of the forced second harmonic is determined by the $(\tilde{\xi},2\omega)$ pair, and the homogeneous one is dictated by the $\hat{\xi}-2\omega$ curve. Both wavenumbers $\tilde{\xi}$ and $\hat{\xi}$ are subjected to (different) correction shifts and depend on the excitation amplitude $A_{0}$. Fig. 7 shows the numerical results of the dispersion shift in the second harmonics. Two scenarios can be identified. In a certain frequency range, the forced and homogeneous solutions co-exist, as shown in Fig. 7(a). In the frequency range where no modes are available, only the forced solutions is observed, as shown in Fig. 7(b). In both cases, the agreement between the analytical and numerical results is remarkable. Comparing Fig. 6(a) against Fig. 7(b), and Fig. 6(b) against Fig. 7(a), we can easily conclude that the wavenumber shifts in the homogeneous and forced second harmonics are significantly larger than those exhibited by the fundamental harmonics under the same amplitudes of the excitation, consistently with the prediction of the multiple scales analysis.

IV Application to waveguides with self-switching capabilities

In this section we demonstrate that the extreme tunability of the second harmonic enabled by the cubic nonlinearity allows controlling the spectral characteristics of the nonlinear wave response due to the quadratic nonlinearity, and activate new modal functionalities at the second harmonic. One interesting application is to design a family of self-adaptive modal switches that can turn on and off some desired functionalities by simply adjusting the amplitude of excitation. Here we illustrate the richness of opportunities resulting from this idea through two illustrative examples.

IV.1 Application 1: Optical mode switch via homogeneous solution self-tuning

The first example exploits the possibility to control the manifestation of the homogeneous second harmonic (due to quadratic nonlinearity) through the cubic term. This allows activating/deactivating or tuning certain modal characteristics typical of the higher-frequency optical regime that are nonlinearly activated in the response via SHG.

This application is based on the notion that, when a wave mode (dispersion branch) is available at the second harmonic, the homogeneous solution will display the spectral modal characteristics of that mode Jiao and Gonella (2018a). If the second harmonic is subjected to a dispersion shift, e.g., through cubic nonlinear correction, the modal characteristics of the activated harmonic can be further manipulated with respect to the default case with only quadratic nonlinearity. This concept is illustrated in Fig. 8 and supported by numerical simulations performed with the same parameters used in section $2$. The excitation frequency $\omega_{0}$ (red dot in Fig. 8(a)) is carefully chosen such that the second harmonic $2\omega_{0}$ is above, but sufficiently close to, the cut-on frequency $\omega^{2}_{\mathrm{cut-on}}$ of the upper branch.

Thus, at low amplitude (or without cubic nonlinearity), the homogeneous second harmonic is located at the onset of the upper branch (around the blue dot in Fig. 8(b)), whose response involves large motion of mass $m_{2}$ (a typical feature of locally resonant bandgaps). The forced solution is marked by the green dot located on the green line, which represents the locus of the $2\xi-2\omega$ pairs. At high amplitude ($A_{0}=1$), the dispersion relation is corrected. As a result, the forced second harmonic is shifted from the location of the green dot to the location of the black dot, while the homogeneous second harmonic is in the bandgap extension produced by the upshifted upper branch (magenta dashed line).

In simulations, the prescribed excitation frequency is chosen to be $\omega_{0}=1.07$ rad/s such that the second harmonic is around $2\omega_{0}=2.14$ rad/s, which is close to the cut-on frequency $\omega^{2}_{\mathrm{cut-on}}$, as desired. A band-pass filter is used to extract the second harmonics from the nonlinear response. Fig. 9 depicts the spatial profiles and the normalized DFT of the filtered second harmonics for mass $m_{1}$ and mass $m_{2}$. With a low-amplitude excitation, as shown in Fig. 9(a), we observe that there are two peaks in the spectral domain: the one with lower wavenumber corresponds to the homogeneous solution, and the other to the forced solution. This is consistent with the spatial profile, where the wave presents an amplitude modulation - a typical feature of waves with two co-existing spectral components. At high amplitude, the first peak is eliminated, and the wavefield becomes monochromatic, as depicted in Fig. 9(b). In comparison, in Fig. 9(c) we consider a case where we retain a high amplitude of excitation but we remove the cubic nonlinearity in the system. We see that all the major wave signatures are qualitatively similar to case (a), albeit with higher displacement amplitude. This indicates that the SHG due to quadratic nonlinearity is triggered under large amplitude-excitations, but the additional correction (filtering of the homogeneous component) due to the cubic term in the foundation is lost. In essence, we have demonstrated that, by playing with the amplitude of the excitation and with cubic nonlinearity to exploit, we obtain a double tuning effect on the wave: 1) we trigger and control the strength of the second harmonic (this effect would be observed even without cubic term); 2) we selectively determine whether the homogeneous solution is observed in the second harmonic in addition to the forced one (this effect is germane to the case with double nonlinearity) and therefore we determine whether two distinct wavenumbers and modal characteristics will co-exist in the response. In summary, by controlling the amplitude, not only do we control the strength of the second harmonic; we also inherently dictate its spectral and modal content.

IV.2 Application 2: Adaptive optical mode switch via PMC tuning

The second example exploits the possibility to interfere with the establishment of PMC in the second harmonic through the cubic correction. Since the strength of the harmonic generation under PMC is largely magnified, the entire modal manifestation of the second harmonic can be dramatically changed depending on whether PMC is enabled or not.

As discussed in the theory section, if PMC is achieved, the amplitude of the solution at $O(\epsilon_{2})$ is proportional to the propagation distance. This unique feature can be utilized to boost SHG, overall emphasizing the nonlinear response. However, for a given system with only quadratic nonlinearity, PMC is strictly linked to specific frequencies, or even not available at all, based on the fixed characteristics of the band diagram. In contrast, by incorporating cubic nonlinearity in the system and by exploiting its manifestation on the correction of the dispersion branches, we can instead switch on/off PMC in a certain frequency range, as shown in Fig. 10. Based on its mathematical definition, the geometrical interpretation of PMC is the existence of an intersection between the linear dispersion curves and the curve encompassing all the forced second harmonic pair (the $2\tilde{\xi}-2\omega$ curve). In light of this, in the system with $\epsilon_{1}=0$ blue(i.e., no cubic nonlinearity) and $\epsilon_{2}=0.05$ depicted in Fig. 10 PMC is located at the blue dot (zoomed details shown in the inset). If some hardening cubic nonlinearity ($\epsilon_{1}=0.1$, $G_{4}=1$) is added, PMC is shifted to the black dot. In contrast, if softening cubic nonlinearity is added ($G_{4}=-1$), PMC are shifted to the dashed black circle. In essence, by resorting to different types of cubic nonlinearity, PMC becomes highly tunable, providing another way to manipulate nonlinear wave propagation at the second harmonic.

To illustrate this capability, we consider two chains. The parameters of chain 1 (with softening foundation) are $m_{1}=m_{2}=1$, $K_{2}=J_{2}=1=1$, $G_{2}=0.92$, $G_{4}=-1$, $\epsilon_{1}=0.1$, and $\epsilon_{2}=0.01$. Chain 2 (with hardening foundation) features $G_{2}=0.84$, $G_{4}=1$, while all the other coefficients are kept identical. For amplitude $A_{0}=1$, the corrected dispersion relations for both systems are plotted in Fig. 11 and superimposed to the linear ones. In the first case, the system features one PMC (blue dot in the inset of Fig. 11(a)) close to the upper cut-off frequency, when excited with low amplitude or without cubic nonlinearity. By increasing the amplitude, the PMC can be switched off due to the activated dispersion shift (black dashed curve). In the second case, there is no PMC available at low amplitude, but a new one (i.e., black dot in Fig. 11(b)) can be created as we increase the amplitude to a certain degree and produce the correction marked by the black solid curve.

Next, we perform numerical simulations for both chains. The filtered second harmonics in the spatial and spectral domains are plotted in Fig. 12. For chain 1, the excitation frequency is carefully chosen such that the second harmonic is close to the blue dot. At low excitation amplitude ($A_{0}=0.1$), PMC is achieved as designed, as confirmed by the fact that the spectrum only contains one dominant wavenumber component, as shown in Fig. 12(a). As a result, the response at the second harmonic grows in space, as predicted from the analytical solution. Note however that, in practice, the harmonic starts losing intensity beyond a certain stage in the simulation. This reduction in amplitude is not captured by our perturbation analysis, because it occurs when the second harmonic, under PMC, is accumulated to a higher order than $O(\epsilon_{2})$, which violates the order assumption in Eq. 3. At high amplitude ($A_{0}=1$), the PMC is switched off, and the second harmonic splits into two parts, one corresponding to the homogeneous response and the other denoting the forced response, which globally result a modulation in amplitude, as shown in Fig. 12(b). In comparison, the result obtained with high amplitude of excitation but without additional cubic nonlinearity is given in Fig. 12(c). We see that PMC is again enabled and the response features wave characteristics similar to those observed for the case of Fig. 12(a).

In the second case, the excitation frequency is chosen such that the second harmonic is close to the black dot, and therefore PMC is not available at low amplitude. As a result, the response at the second harmonic has two components as shown in the DFT plot of Fig. 12(d). As we increase the amplitude to $1$, the dispersion shift leads to a nonlinearly-generated PMC, and the corresponding quasi-monochromatic numerical result is plotted in Fig. 12(e). However, if we remove the cubic nonlinearity in the system, PMC is again unavailable, as confirmed by the resulting modulated filtered second harmonic shown in Fig. 12(f). Similar to the previous example, the chain behaves like a self-adaptive switch, which can turn on and off selected modal characteristics typical of the optical mode via the interplay between cubic nonlinearity and SHG. Unlike the other case, however, here the mechanism is not the selective excitation and suppression of one component (homogeneous) of the second harmonic, but rather the selective activation of PMC and of its amplitude-reinforcing effects.

V Conclusion

In this paper, we have studied wave propagation in a quadratic nonlinear spring-mass chain attached to a locally-resonant cubic nonlinear foundation, using a multiple scales analysis and full-scale simulations. We have been able to link the interplay between quadratic and cubic nonlinearity to a number of new nonlinear phenomena, as well as to revisit some more classical tuning and manipulation effects from a different perspective.

First, we have revisited the amplitude-dependent dispersion shift of the fundamental harmonic, which has been known to be the most tangible effect of cubic nonlinearity. Here, we have shown that this effects has significantly different manifestations according to whether initial conditions or boundary conditions are prescribed. We have analyzed with special care the boundary conditions case, given its natural occurrence in practical problems of engineering interest.

We have then considered the phenomenon of second harmonic generation, which is the prominent signature of quadratic nonlinearity. We have shown that, similarly to the fundamental harmonic, both homogeneous and forced second harmonics also feature amplitude-dependent dispersion shifts due to the Q-C nonlinearity interplay. We have presented two examples of tunable chains with self-adaptive switching capabilities. The first example involves activating or deactivating new modal characteristics in the nonlinear response through the dispersion shift induced in the homogeneous second harmonic. The second example shows the capability of switching on and off phase matching conditions to activate and deactivate strong nonlinear signatures at high frequencies. In conclusion, the self-interacting effect of cubic and quadratic nonlinearities provides a new conceptual platform for the design of nonlinear phononic crystals and metamaterials with superior tunability characteristics.

Acknowledgements

The authors acknowledge the support of the National Science Foundation (CAREER Award CMMI-$1452488$).

Appendix A Derivation of $f_{u}^{2}$

The complete expression for the non-zero forcing term $f_{u}^{2}$ at $O(\epsilon_{2})$, due to the quadratic nonlinearity, is

[TABLE]