Switchable phonon diodes using nonlinear topological Maxwell lattices

TL;DR

This paper demonstrates how topological Maxwell lattices can be engineered to create tunable, non-reciprocal elastic wave transmission, functioning as phonon diodes with switchable behavior through topological phase transitions.

Contribution

It introduces a novel design for topologically protected phonon diodes using nonlinear Maxwell lattices with switchable non-reciprocal wave transmission.

Findings

Achieved non-reciprocal elastic wave transmission using topological floppy edge modes.

Demonstrated switching of diode behavior via topological phase transitions.

Proposed cellular metamaterials for tunable topologically protected phonon diodes.

Abstract

Recent progress in topological mechanics have revealed a family of Maxwell lattices that exhibit topologically protected floppy edge modes. These modes lead to a strongly asymmetric elastic wave response. In this paper, we show how topological Maxwell lattices can be used to realize non-reciprocal transmission of elastic waves. Our design leverages the asymmetry associated with the availability of topological floppy edge modes and the geometric nonlinearity built in the mechanical systems response to achieve the desired non-reciprocal behavior, which can be further turned into strongly one-way phonon transport via the addition of on-site pinning potentials. Moreover, we show that the non-reciprocal wave transmission can be switched on and off via topological phase transitions, paving the way to the design of cellular metamaterials that can serve as tunable topologically protected phonon diodes.

Equations

$$m{\ddot{\text{u}}}_g^{(1)}={\text{F}}_g^{\mathrm{ext}}-\text{D}{\text{u}}_g^{(1)}-\eta {\dot{\text{u}}}_g^{(1)},$$

$$m{\ddot{\text{u}}}_g^{(2)}={\text{f}}^{(2)}({\text{u}}_g^{(1)})-\text{D}{\text{u}}_g^{(2)}-\eta {\dot{\text{u}}}_g^{(2)},$$

$$m{\ddot{\vec{u}}}_n=-\eta {\dot{\vec{u}}}_n-{\vec{\nabla }}_{n}V+{\vec{F}}_n^{\mathrm{ext}}(t),$$

$$\text{u}(\omega )=\text{G}(\omega ){\text{F}}^{\mathrm{ext}}(\omega ),$$

$$u_{\alpha }(\omega )=G_{\alpha }{F}_{\alpha }^{\mathrm{ext}}(\omega ),$$

$$m{\ddot{u}}_{n,g}$$

$${\vec{T}}_{-1}={\vec{T}}_N=0,$$

$$f_n^{(1)}=K[c_1{c}_2(u_{n+1}+u_{n-1})-(c_1^2+c_2^2)u_n],$$

$$c_{1,2}=\frac{(a\pm 2r\sin \bar{\theta })\cos \bar{\theta }}{\sqrt{a^2+4r^2{\cos }^2\bar{\theta }}},$$

$$f_n^{(2)}(u_n)$$

$$2C_1=C_5=\frac{K}{4m}(-1+\sqrt{2})$$

$$2C_4=C_2=\frac{K}{4m}(-1-\sqrt{2})$$

$$C_3=\frac{K}{24m}(5-\sqrt{2})C_6=\frac{K}{24m}(5+\sqrt{2}).$$

$$m{\ddot{\text{u}}}_g^{(1)}+m{\ddot{\text{u}}}_g^{(2)}+\mathcal{O}({\lambda }^3)=\left({\text{F}}^{\mathrm{ext}}-\text{D}{\text{u}}_g^{(1)}-\eta {\dot{\text{u}}}_g^{(1)}\right)$$

$$+\left({\text{f}}_g^{(2)}({\text{u}}_g^{(1)})-\text{D}{\text{u}}_g^{(2)}-\eta {\dot{\text{u}}}_g^{(2)}\right)+\mathcal{O}({\lambda }^3),$$

$${\text{u}}_g^{(1)}(\omega )=\text{G}(\omega ){\text{F}}^{\mathrm{ext}}(\omega ),$$

$${\text{u}}_g^{(2)}(2\omega )=\text{G}(2\omega ){\text{f}}_g^{(2)}({\text{u}}_g^{(1)}),$$

$$\text{G}(\omega )=[\text{D}+(-m{\omega }^2+i\eta \omega )\text{I}{]}^{-1}.$$

$$u_{n,g}^{(1)}=a_g{\lambda }_1^n+A_g{\lambda }_2^{n}0\le n\le g$$

$$u_{n,g}^{(1)}=b_g{\lambda }_1^n+B_g{\lambda }_2^{n}g\le n\le N,$$

$${\lambda }_{1,2}=\frac{1}{2}\left(-\gamma \pm \sqrt{{\gamma }^2-4}\right),$$

$${\chi }_{\mathrm{in}+}(\omega )=(m{c}_1|c_1-c_2{\lambda }_1|{)}^{-1}$$

$${\chi }_{\mathrm{in}-}(\omega )=(m{c}_2|c_2-c_1{\lambda }_1|{)}^{-1},$$

$$\displaystyle\begin{array}[]{ccc}{\chi_{{\rm in}+}}/{\chi_{{\rm in}-}}>1,&\qquad\omega<\omega^{*},\\ {\chi_{{\rm in}+}}/{\chi_{{\rm in}-}}<1,&\qquad\omega^{*}<\omega<c_{1}-c_{2},\end{array}\qquad$$

$${\chi }_{\mathrm{out}-}^{(1)}(\omega )={\chi }_{\mathrm{out}+}^{(1)}(\omega )=\frac{1}{m}\frac{{\lambda }_1^{-N-1}({\lambda }_2-{\lambda }_1)}{{\omega }^2-i\eta \omega }$$

$${\text{u}}_{g^{\prime }g}^{(2)}(2\omega )=\text{G}(2\omega ){\text{f}}_{g^{\prime }g}^{(2)}(2\omega ),$$

$$u_{n,g^{\prime }g}^{(2)}=c_{g^{\prime }g}{\mu }_1^n+C_{g^{\prime }g}{\mu }_2^{n}0\le n\le g^{\prime }$$

$$u_{n,g^{\prime }g}^{(2)}=d_{g^{\prime }g}{\mu }_1^n+D_{g^{\prime }g}{\mu }_2^n{g}^{\prime }\le n\le N,$$

$${\mu }_{1,2}=\frac{1}{2}\left(-\nu \pm \sqrt{{\nu }^2-4}\right),$$

$${\chi }_{\mathrm{out}+}^{(2)}$$

$${\chi }_{\mathrm{out}-}^{(2)}$$

$${\text{D}}_i=K{\hat{n}}_i{\hat{n}}_{i}i=1,2,...,6$$

$$\displaystyle\textbf{C}_{1}=\left(\begin{array}[]{cccccc}\cos\theta_{1}&\sin\theta_{1}&-\cos\theta_{1}&-\sin\theta_{1}&0&0\\ 0&0&\cos\theta_{2}&\sin\theta_{2}&-\cos\theta_{2}&-\sin\theta_{2}\\ -\cos\theta_{3}&-\sin\theta_{3}&0&0&\cos\theta_{3}&\sin\theta_{3}\\ -e^{ik_{1}}\cos\theta_{4}&-e^{ik_{1}}\sin\theta_{4}&0&0&0&0\\ e^{ik_{1}}\cos\theta_{5}&e^{ik_{1}}\sin\theta_{5}&-\cos\theta_{5}&-\sin\theta_{5}&0&0\\ 0&0&\cos\theta_{6}&\sin\theta_{6}&0&0\\ \end{array}\right),\textbf{C}_{2}=\left(\begin{array}[]{cccccc}0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&\cos\theta_{4}&\sin\theta_{4}\\ 0&0&0&0&0&0\\ 0&0&0&0&-\cos\theta_{6}&-\sin\theta_{6}\\ \end{array}\right).\qquad$$

$${\vec{F}}_i=K(\hat{n}\hat{n})\Delta {\vec{u}}_{ij}$$

Full text

I Introduction

Over the past decade, significant progress on the development of mechanical analogs of topological states of matter has fueled the new field of “topological mechanics” Süsstrunk and Huber (2015); Kane and Lubensky (2014); Pal et al. (2016); Mousavi et al. (2015); Prodan and Prodan (2009). A plethora of applications of topological mechanical metamaterials have been proposed, such as uni-directional transport using chiral edge modes Nash et al. (2015); Mitchell et al. (2018a, b); Mousavi et al. (2015); Zhang et al. (2017); Dai et al. (2018); Guo et al. (2017); Süsstrunk and Huber (2015); Yang et al. (2015), transformable topological materialsRocklin et al. (2017), structures with programmed buckling or fracturing patternsZhang and Mao (2018); Paulose et al. (2015a), mechanical laserAbbaszadeh et al. (2017), aperiodic topological metamaterialsZhou et al. (2018, 2019), geared topological metamaterialsMeeussen et al. (2016), and dislocation-localized softnessPaulose et al. (2015b).

A particularly interesting potential application of topological metamaterials is to obtain a phonon diode, i.e., a device that only allows sound transmission in one direction. The main requirement to achieve this goal is to break reciprocity. Within linear elasticity, systems with time-reversal symmetry exhibit reciprocityMaznev et al. (2013). According to Maxwell-Betti’s theoremMaxwell (1864); Betti (1872); Charlton (1960), reciprocity implies that $u_{B}^{j,(1)}/F_{A}^{i}=u_{A}^{i,(1)}/F_{B}^{j}$, where $i,j=x,y,z$ are Cartesian components, $F_{A}^{i}$ is the $i$-th component of the external force exerted at input point $A$ and $u_{B}^{j,(1)}$ is the $j$-th component of the linear elastic response probed at output point $B$. In the remainder of this paper, we define the quantity $\chi_{{\rm out},A}^{(1)}=u_{B}^{j,(1)}/F_{A}^{i}$ as the linear transmission susceptibility. To achieve non-reciprocal transmission one needs to 1) break spatial inversion symmetry and 2) either break time-reversal symmetry, or include nonlinear effects. Major efforts have been devoted to the development of strategies to violate reciprocity by breaking time reversal symmetry. For example, several active metamaterial configurations have been proposed for uni-directional edge wave propagation, such as systems of coupled gyroscopesWang et al. (2015); Nash et al. (2015); Mitchell et al. (2018a, c), chiral active fluids and plasmaEngheta et al. (2005); Souslov et al. (2019), dynamic phononic latticesWang et al. (2018), spatio-temporally modulated metamaterialsZanjani et al. (2014); Trainiti and Ruzzene (2016); Nassar et al. (2017); Swinteck et al. (2015); Fleury et al. (2016); Salerno et al. (2016); Chaunsali et al. (2016); Peano et al. (2015) and active-liquid metamaterialsFleury et al. (2014); Souslov et al. (2017).

An alternative route to break reciprocity in mechanical systems consists of leveraging the intrinsic nonlinearity of their elastic response. Recent implementations include nonlinear self-demodulation processes obtained by coupling elastically distinct layers of metamaterialsDevaux et al. (2015); Merkel et al. (2014); Cebrecos et al. (2016), unidirectional guiding of strongly nonlinear transition waves in a bistable latticeNadkarni et al. (2016), static non-reciprocal elasticsCoulais et al. (2017), acoustic switching and rectificationBoechler et al. (2011); Liang et al. (2010); Liu et al. (2015) and broadband acoustic diodesGu et al. (2016). In this paper, we present an approach for non-reciprocal wave transmission in lattice systems, in which the task of breaking space inversion symmetry is accomplished through the activation of topological floppy edge modes, and the nonlinear response requirements are fulfilled by the geometric nonlinearity of the lattice deformation. The main advantage of the proposed design stems from the topological protection of the edge modes, which endows the non-reciprocal phenomena with robustness against potential defects and disorder.

Maxwell lattices are central-force lattices with average coordination number $\langle z\rangle=2d$ ($d$ is the spatial dimension), which puts them on the verge of mechanical instabilityZhou et al. (2018); Kane and Lubensky (2014); Rocklin et al. (2016); Lubensky et al. (2015). They host topologically protected edge modes at zero frequency (floppy modes) which are governed by the topology of the equilibrium and compatibility matrices and therefore ultimately depend on the lattice geometryKane and Lubensky (2014). The topological edge modes lead to strongly asymmetric edge stiffness, which has been shown to result in asymmetric wave propagation characteristics, whereby certain edges allow waves to propagate into the bulk, and others localize energy at the boundariesMa et al. (2018). Despite this asymmetry, the transmission of linear elastic waves is still reciprocal, meaning that the linear transmission susceptibilities $\chi_{{\rm out},A}^{(1)}$ and $\chi_{{\rm out},B}^{(1)}$ between points A and B in space are equal, in accordance with Maxwell-Betti’s theorem.

To achieve non-reciprocity, we need to operate the lattice in the nonlinear regime. Here, we pursue this by propagating second harmonics whose generation and intensity are controlled by the amplitude of the excitation. The generation of higher harmonics in mechanical systems with multi-modal dispersive behavior and the resulting opportunities for unconventional wave manipulation and functionality enrichment in elastic metamaterials have been the object of a number of recent studiesDeng et al. (2005); Matlack et al. (2011); Ganesh and Gonella (2015, 2017); Jiao and Gonella (2019); Pal et al. (2018); Jiao and Gonella (2018a, b).

In this paper we show that second harmonic modes in topological Maxwell lattices are strongly nonreciprocal, due to the contrast in stiffness between floppy and non-floppy edges, which is a topologically protected property. Further, we demonstrate that, by blocking the linear (first harmonic) modes via on-site pinning potentials (which can be realized by placing the lattice on a soft substrate), the system works as a phonon diode, in which transmission (with frequency doubled) is effectively observed in only one direction. Finally, we revisit the notion that topological kagome lattices can be reversibly transformed between different topological states with contrasting edge state landscapes through a transformation, known as the “Guest mode”, which involves a soft strain of the whole lattice. As a result, these lattices can be switched between strongly non-reciprocal and nearly reciprocal states through simple reversible operations.

II Non-reciprocity in 1D topological mechanical chain

We start our discussion by revisiting the 1D topological mechanical chain introduced inKane and Lubensky (2014), as shown in fig.1(a). This is the simplest lattice with topologically protected floppy edge modes that give rise to contrasting boundary rigidity. The chain consists of rigid rotors connected to fixed pivot points separated by lattice constant $a$. The pivot points as well as the rotors are labeled from 0 to $N$. The other ends of the rotors feature particles of mass $m$, and neighboring particles are connected with harmonic springs. The chain is subjected to open boundary conditions (OBC) at rotors [math] and $N$. The equilibrium configuration is such that rotors form an angle $\bar{\theta}$ relative to the upward and downward normals. The angular displacements are denoted as $\textbf{u}=(r\delta\theta_{0},r\delta\theta_{1},...,r\delta\theta_{N})$, where $\delta\theta_{n}=\theta_{n}-\bar{\theta}$. The system consists of $N+1$ degrees of freedom and $N$ constraints, leaving only one topological floppy mode localized on the right boundary.

Now we imagine driving the chain by a monochromatic harmonic force $F_{g}^{\rm ext}(t)=Fe^{i\omega t}$ applied at the left (right) end, i.e., on rotor $g=0$ ($g=N$), while $F_{n}^{\rm ext}(t)=0$ elsewhere. $F$ is assumed small enough that $\delta\theta_{n}\ll 1$, $\forall n=0,...,N$, which validates perturbation theory. We denote $\textbf{F}_{g}^{\rm ext}=(F_{0},F_{1},...,F_{N})$ as the array of external forces and, as we mentioned previously, $\textbf{F}_{0}^{\rm ext}=(F,0,...,0)$ and $\textbf{F}_{N}^{\rm ext}=(0,0,...,F)$. By expanding $\textbf{u}=\textbf{u}^{(1)}+\textbf{u}^{(2)}+\mathcal{O}(F^{3})$ in orders of $F$, we can solve for the linear elastic mode $\textbf{u}^{(1)}$ and for second harmonic mode $\textbf{u}^{(2)}$, respectively.

We define the input linear response function as $\chi_{\rm in}=|u_{\rm in}^{(1)}|/F_{\rm in}$, where $u^{(1)}_{\rm in}$ is the linear displacement of the rotor that is being driven. We also define $\chi_{\rm out}^{(1)}=|u_{\rm out}^{(1)}|/F_{\rm in}$ ($\chi_{\rm out}^{(2)}=|u_{\rm out}^{(2)}|/F_{\rm in}$), where $u^{(1)}_{\rm out}$ ($u^{(2)}_{\rm out}$) is the linear displacement (second harmonic displacement) at the boundary rotor opposite to the driven side.

To the linear order, Newton’s equation of motion is

[TABLE]

where $\eta$ is the damping coefficient, $m$ is the particle mass, and the lower index $g$ indicates that the force is applied at the left end if $g=0$ (right end if $g=N$). The dynamical matrix is $\textbf{D}=K\textbf{C}^{T}\textbf{C}$, where $\textbf{C}_{ij}=c_{1}\delta_{ij}-c_{2}\delta_{i,j+1}$ is the compatibility matrix, and $c_{1}$, $c_{2}$ are coefficients determined by lattice geometry (see SI.II). In the static limit, the linear elasticity is captured by the compatibility matrix C: the floppy mode is localized on the right (left) end if $|c_{1}|>|c_{2}|$ ($|c_{1}|<|c_{2}|$). The topological protection of this floppy mode arises from the winding number of the compatibility matrix in the complex plane, and is therefore invariant against continuous deformations to the geometry of this 1D chain unless the gap closesKane and Lubensky (2014). Without losing generality, in the rest of this paper we always let $c_{1}>c_{2}>0$ by allowing rotors to tilt rightwards with $\bar{\theta}>0$, so the edge floppy mode is localized on the right end of the chain. Following this convention, it is convenient to denote $\chi_{\rm in+}=|u_{N}^{(1)}|/F_{N}$, $\chi_{\rm in-}=|u_{0}^{(1)}|/F_{0}$, and $\chi_{\rm out+}^{(1,2)}=|u_{0}^{(1,2)}|/F_{N}$, $\chi_{\rm out-}^{(1,2)}=|u_{N}^{(1,2)}|/F_{0}$, where $+$ ($-$) indicates that external driving is applied at the soft (rigid) end.

Above zero frequency, the dispersion relation of the bulk phonon mode reads $\omega=[(c_{1}-c_{2})^{2}+4c_{1}c_{2}\sin^{2}\left(\frac{1}{2}ak\right)]^{1/2}$, where $k$ is the wave number. Linear elastic waves driven by external excitations with $c_{1}-c_{2}<\omega<c_{1}+c_{2}$ are bulk modes. $\chi_{\rm in\,+}$ and $\chi_{\rm in\,-}$ are at the same order of magnitude. As $\omega$ falls below $c_{1}-c_{2}$, linear modes localize on lattice boundaries. $\chi_{\rm in\,+}/\chi_{\rm in\,-}$ monotonically increases to infinity as $\omega$ approaches the static limit (see SI.II for details). Although the stiffness differs dramatically (by orders of magnitude) on opposite boundaries, the linear elastic transmission is still reciprocal, meaning that $\chi_{\rm out+}^{(1)}=\chi_{\rm out-}^{(1)}$ as a manifestation of Maxwell-Betti’s theorem. We verified this equality both analytically and numerically, as shown in fig.1(d).

Interestingly, higher order harmonics with $\omega^{(n)}=n\omega$ that are nonlinearly generated by the edge modes are bulk modes as long as $c_{1}-c_{2}<\omega^{(n)}<c_{1}+c_{2}$. In what follows, we study whether these second harmonic modes carry non-reciprocal characteristics. The Newton’s equation of motion for second harmonic modes is

[TABLE]

where $\textbf{f}^{(2)}(\textbf{u}_{g}^{(1)})$ is the second harmonic effective driving force generated by the linear displacement $\textbf{u}_{g}^{(1)}$, as defined in Eq.(1) (see SI.II for details). Since the effective driving is quadratic in $\textbf{u}^{(1)}_{g}$, it triggers second harmonic modes with amplitude $|u_{n}^{(2)}|\propto|u_{\rm in}^{(1)}|^{2}$ and frequency $2\omega$. External excitations with $\frac{1}{2}(c_{1}-c_{2})<\omega<{\rm min\,}(\frac{1}{2}(c_{1}+c_{2}),c_{1}-c_{2})$ excite linear edge modes as well as second harmonic bulk waves. For a given magnitude of excitation $F$, the input-end linear response measured at the right edge is far greater than its counterpart measured at the left edge ($\chi_{\rm in+}\gg\chi_{\rm in-}$). As a result, the global wave amplitude experienced by the chain is much greater when the chain is driven from the right, leading also, in return, to significantly stronger second harmonics generation, i.e., $|\textbf{u}_{g=N}^{(2)}|\gg|\textbf{u}_{g=0}^{(2)}|$. Through analytical and numerical calculations we can show that $\chi_{\rm out+}^{(2)}\gg\chi_{\rm out-}^{(2)}$, meaning that the transmission of second harmonics is non-reciprocal, as reported in fig.1(e). This non-reciprocity result can be generalized to the $n$-th harmonic mode: we obtain that $\chi_{\rm out+}^{(n)}\gg\chi_{\rm out-}^{(n)}$ if $\frac{1}{n}(c_{1}-c_{2})<\omega<{\rm min\,}(\frac{1}{n}(c_{1}+c_{2}),c_{1}-c_{2})$. We note that, besides the low-frequency regime $\omega<c_{1}-c_{2}$, linear modes with high frequencies $\omega>c_{1}+c_{2}$ can also localize on edges. However, they are not of interest in this paper, because the associated nonlinear harmonics are also edge excitations which cannot propagate across the lattice and therefore cannot contribute to transmission.

It is interesting to ask whether non-reciprocity still holds if the monochromatic harmonic excitation is replaced by a tone burst excitation with carrier frequency $\omega$ and Gaussian amplitude modulation, having the form $F(t)\sim Fe^{i\omega t-(t-t_{0})^{2}/\tau^{2}}$, where the parameter $\tau$ controls the spread of the Gaussian and $t_{0}$ denotes the trigger time of the packet. Since, in Fourier space, the input signal is a Gaussian function with full width at half maximum $\Delta\omega=2\sqrt{\ln 2}/\omega\tau$, we expect that the transmission of nonlinear modes is still non-reciprocal. This conjecture is verified by numerical analysis as shown in fig.1(g).

It is important to note that the key ingredient to achieve non-reciprocity is the contrast in rigidity between opposite edges, and not the topological protection of the edge modes. In principle, any system with asymmetric boundary stiffness, whether this is topological or not, can realize non-reciprocity if such asymmetry is used in conjunction with nonlinear elasticityCoulais et al. (2017). However, topologically protected floppy modes enjoy the additional attribute of being robust against disorder, noise, and stochastic damage. More interestingly, topological lattices are switchable, meaning that the topological polarization can be changed via simple, reversible operations that modify their geometry. For example, the 1D chain discussed above can be flipped to have the opposite topological polarization by propagating a soliton through the chainChen et al. (2014). As we shall discuss in section IV, 2D topological kagome lattices can undergo a geometric change through a soft strain of the whole lattice, called the “Guest mode”Guest and Hutchinson (2003), to manipulate topological phases, control floppy mode localizationRocklin et al. (2017) and thus boundary stiffness. Consequently, the transmission of nonlinear waves can be switched from non-reciprocal to reciprocal by reconfiguring the lattices from their topological to their non-topological form.

III Non-reciprocity in topological kagome lattice

Having established non-reciprocity for a 1D topological chain, we now ask if the same is true for a 2D topological latticeKane and Lubensky (2014). To this end, we consider the topological kagome lattice shown in fig.2(a). The lattice is ideal, i.e., it consists of point masses connected by nearest-neighbor linear springs. The unit cell contains one equilateral and one isosceles triangles, which are constructed from 6 bonds and 3 nodes marked by A, B and C. The side length of the equilateral triangle and the longer edge of the isosceles triangle are $l_{0}$, while the shorter edge is $l_{0}/\sqrt{3}$. The twist angle of the isosceles triangle is $5^{\circ}$ counterclockwise, which makes the longer edge of isosceles triangle inclined by $5^{\circ}$ relative to the bottom edge of equilateral triangle (marked $\overline{\rm AB}$ in fig.2(a)). $\vec{a}_{1}$ and $\vec{a}_{2}$ are the lattice primitive vectors. The lattice, spanning the area $|N_{1}\vec{a}_{1}\times(N_{2}-1)\vec{a}_{2}|$ in real space, is composed of $N_{1}\times(N_{2}-1)$ unit cells and an additional layer of $C$-sites at the bottom edge of the lattice to complete the triangles. It is subjected to periodic boundary condition in $\vec{a}_{1}$ and open boundary conditions at the top and bottom edges.

We start by introducing a supercell analysis of this lattice. For convenience we denote $k_{1}=\vec{k}\cdot\vec{a}_{1}$ and $k_{2}=\vec{k}\cdot\vec{a}_{2}$ as the wavenumbers along the primitive vectors. We further decompose the lattice in supercell strips directed along $\vec{a}_{2}$, as marked in cyan in fig.2(a), and we apply Bloch’s conditions along $\vec{a}_{1}$. Within a supercell, the unit cells are labeled from 1 to $N_{2}-1$ going from top to bottom, and the C-sites on the bottom layer are labeled as $N_{2}$. The internal nodal displacements of unit cell $n_{2}$ (with $1\leq n_{2}\leq N_{2}-1$) are denoted as $\textbf{u}_{n_{2}}=(u_{n_{2}A}^{x},u_{n_{2}A}^{y},u_{n_{2}B}^{x},u_{n_{2}B}^{y},u_{n_{2}C}^{x},u_{n_{2}C}^{y})$. The displacement field of supercell strip is therefore denoted as $\textbf{u}=(\textbf{u}_{1},\textbf{u}_{2},...,\textbf{u}_{N_{2}})$. The lattice is driven by a monochromatic harmonic force acting vertically and with amplitude varying periodically in the $\vec{a}_{1}$ direction, i.e., $\vec{F}_{n_{1},g}^{\rm ext}(t)=e^{i\omega t-ik_{1}n_{1}}\vec{F}_{C}$ with $\vec{F}_{C}=(0,F)$ at the top boundary (bottom boundary) of C-sites labeled by $g=1$ ($g=N_{2}$) and $\vec{F}_{n_{1},n_{2}}^{\rm ext}(t)=0$ otherwise. $F$ is assumed to be small such that all $|\vec{u}_{n_{1}n_{2},A}|,|\vec{u}_{n_{1}n_{2},B}|$ and $|\vec{u}_{n_{1}n_{2},C}|\ll l_{0}$, validating perturbation theory. By expanding $\textbf{u}=\textbf{u}^{(1)}+\textbf{u}^{(2)}+\mathcal{O}(F^{3})$, we can solve Newton’s equation of motion for the linear mode $\textbf{u}^{(1)}$ and for second harmonic mode $\textbf{u}^{(2)}$, as detailed in the SI.

The analysis of wave propagation in the 2D lattice follows the steps used for the 1D chain, albeit with the additional wavenumber $k_{1}$ describing spatial variation in the horizontal direction (with periodic boundary conditions). Specifically, the dynamical matrix of this super cell strip is $\textbf{D}=K\textbf{C}^{\dagger}(k_{1})\textbf{C}(k_{1})$. $\textbf{C}(k_{1})$ is the compatibility matrix given by $\textbf{C}_{ij}(k_{1})=\textbf{C}_{1}(k_{1})\delta_{ij}+\textbf{C}_{2}(k_{1})\delta_{i+1,j}$, where $\textbf{C}_{1}(k_{1})$ and $\textbf{C}_{2}(k_{1})$ are intra-cell and inter-cell compatibility matrices, respectively (see SI.III for details). To the linear order of displacement, Newton’s equation of motion is the same as Eq.(1), where $\eta$ and $m$ are damping coefficient and particle mass, respectively, and $g$ indicates that the input force is applied at the layer of C-sites indexed $g=1$ at the top (layer of C-sites indexed $g=N_{2}$ at the bottom) of the lattice. The static system is characterized by the polarization vector $\vec{R}_{T}$, which is a topological invariant. Mechanical lattices with a well-defined polarization exhibit topological floppy edge modes exponentially localized at the boundary towards which $\vec{R}_{T}$ points. The configuration of fig.2(a) has a polarization vector $\vec{R}_{T}=\vec{a}_{1}-\vec{a}_{2}$. The floppy modes are therefore localized on the top edge, making this edge much softer than the bottom one. We use lower index $+$ ($-$) to indicate that the external signal is applied where the floppy modes are localized (opposite to the floppy mode localization). It is therefore convenient to denote $\chi_{\rm in+}=|u_{1}^{y,(1)}|/F_{1}^{y}$ and $\chi_{\rm in-}=|u_{N_{2}}^{y,(1)}|/F_{N_{2}}^{y}$ as the soft edge and rigid edge linear response functions, respectively. Similarly, we denote $\chi_{\rm out+}^{(1)}=|u_{N_{2}}^{y,(1)}|/F_{1}^{y}$ ($\chi_{\rm out+}^{(2)}=|\vec{u}_{N_{2}}^{(2)}|/F_{1}^{y}$) and $\chi_{\rm out-}^{(1)}=|u_{N_{2}}^{y,(1)}|/F_{1}^{y}$ ($\chi_{\rm out-}^{(2)}=|\vec{u}_{N_{2}}^{(2)}|/F_{1}^{y}$) as the linear (second harmonic) transmission susceptibilities driven by external forces applied at the soft and rigid boundaries, respectively.

Linear wave propagation is governed by the supercell band structure which stems from the eigenvalue problem $\det(\textbf{D}-m\omega^{2}\textbf{I})=0$. The band structure is gapped except for the trivial translational zero mode point at $k_{1}=0$. Given the wave number $k_{1}$ of the applied force, the linear response is a bulk mode if $\omega>\Delta(k_{1})$, where $\Delta(k_{1})$ is the lowest bulk eigenvalue in the band structure, and $\chi_{\rm in+}$ and $\chi_{\rm in-}$ are of the same order of magnitude. As $\omega$ falls below $\Delta(k_{1})$, linear modes localize on the soft boundary of the lattice. $\chi_{\rm in\,+}/\chi_{\rm in\,-}$ monotonically increases to infinity as $\omega$ approaches zero. Despite the contrasting boundary stiffness at low-frequencies, the linear elastic transmission is still reciprocal, i.e., $\chi_{\rm out+}^{(1)}(k_{1})=\chi_{\rm out-}^{(1)}(-k_{1})$, similar to what we obtained for the 1D topological chain. We validate this equality through analytical and numerical calculations, as shown in fig.2(d).

While linear modes can localize on the lattice boundaries, nonlinearly generated components with $(\omega^{(n)},k_{1}^{(n)})=(n\omega,nk_{1})$ can be bulk waves as long as $\omega^{(n)}>\Delta(k_{1}^{(n)})$. The equation of motion for the second harmonic mode is given by Eq.(2). External excitations with frequency $\frac{1}{2}\Delta(2k_{1})<\omega<\Delta(k_{1})$ generate linear boundary modes and second harmonic bulk waves. Moreover, given the same magnitude $F$ of external force, the input-end frequency response function of the floppy edge is far greater than that of the hard edge ($\chi_{\rm in+}\gg\chi_{\rm in-}$), which renders second harmonic bulk modes excited at the floppy edge much greater than those excited at the hard edge, i.e., $|\textbf{u}_{g=N_{2}}^{(2)}|\gg|\textbf{u}_{g=1}^{(2)}|$. The transmission of second harmonics is therefore strongly non-reciprocal, with $\chi_{\rm out+}^{(2)}\gg\chi_{\rm out-}^{(2)}$, which is verified numerically as shown in fig.2(e). This conclusion can be generalized to the non-reciprocal transmission of $n$-th harmonic mode with $\chi_{\rm out+}^{(n)}\gg\chi_{\rm out-}^{(n)}$ if $\frac{1}{n}\Delta(nk_{1})<\omega<\Delta(k_{1})$.

While, so far, the analysis has followed almost verbatim the same steps of the 1D problem, one important difference is that the 2D lattice phonon band depends on $k_{1}$ (the wave number in the horizontal direction imposed along the boundary). Thus, the width of the gap $\Delta$ and the resulting availability of nonreciprocal propagation depend on the choice of $k_{1}$. In particular, the lattice always has translational zero modes since $\lim_{k_{1}\to 0}\Delta(k_{1})\to 0$. As a result, if we drive the system with a point force applied at a given location on the boundary (which ostensibly excites all values of $k_{1}$), we are bound to observe weaker signatures of non-reciprocity. In other words, the differences in behavior observed by exciting the soft and hard edges will be vastly reduced, as the long wavelength components of the excited linear waves are in both cases bulk modes that do not display asymmetry. Moreover, despite the strong non-reciprocity, this kagome lattice cannot be, strictly speaking, defined as a proper phonon diode. This is because the linear mode, which is reciprocal and always transmitted both ways, is much stronger than the second harmonic mode and always dominates the total response, completely overshadowing any asymmetry in the nonlinear response. In order to mitigate the aforementioned issues, in the next section we propose an evolution of the lattice design meant to work as a proper phonon diode for all wavenumbers.

IV Non-reciprocity in kagome lattice with on-site potential

We consider a modification of the topological kagome lattice discussed in section III, where a weak on-site pinning potential is added to every mass point $i$, with $V_{i}=\frac{1}{2}K^{\prime}\vec{u}_{i}^{2}$ and $K^{\prime}\ll K$. This operation, which practically elastically connects each site to a fixed ground point, can be thought of as the equivalent model of placing the lattice on a soft substrate (or soft elastic foundation). It penalizes particles from moving away from their rest positions, and therefore eliminates the lattice trivial translational zero modes. Since the weak pinning rigidly shifts the band up by $\Delta^{\prime}=\sqrt{K^{\prime}/m}$, now the signals with frequency below $\Delta^{\prime}$ will excite edge modes. The weak pinning does not change the landscape of asymmetric boundary stiffness of the lattice, meaning that $\chi_{\rm in+}\gg\chi_{\rm in-}$ still holds as long as $\omega<\Delta^{\prime}$. Thus, linear edge modes are still preferentially localized on the soft boundary.

The fundamental consequence of having a full low-frequency gap is that, as long as $\frac{1}{2}\Delta^{\prime}<\omega<\Delta^{\prime}$, an external excitation that is periodic along the boundary will excite linear edge modes and second harmonic bulk modes for any arbitrary wavenumber $k_{1}$. Consequently, the non-reciprocal behavior will be observed in the response to a point excitation prescribed at a given location on the boundary, thus eliminating the limitation of the previous configuration. As shown in fig.3(d), the second harmonic positive transmission is indeed stronger than the linear transmission for a point excitation. We can conclude that the topological lattice with pinning potential is now a well-defined phonon diode.

These results still hold when finite bending stiffness at the hinges is includedMa et al. (2018). In fig.3(e) our numerical results indeed show that, when bending stiffness is introduced by adding the contribution of next-nearest-neighbor (NNN) interactions, the non-reciprocal transmission is still significant. Finally, the results also hold for a Gaussian tone burst excitation $F(t)\sim Fe^{i\omega t-(t-t_{0})^{2}/\tau^{2}}$. The numerical analysis results reported in fig.3(i) show that the transmission is still non-reciprocal, similar as the result of 1D topological mechanical chain.

V Lattice reconfiguration and non-reciprocity switching

An interesting feature of Maxwell lattices with spring-mass interactions is that they can undergo uniform soft deformations, in which all the unit cells are twisted in the same fashion while leaving the bond lengths unstretched. Such uniform deformation, known as the “Guest mode”, can manipulate the geometrical parameters that controls the topological phase of the kagome lattice and, consequently, its polarization and the rigidity established on opposite boundariesRocklin et al. (2017). Starting from the unit cell configuration in fig.2(a), by uniformly rotating all the isosceles triangles counterclockwise by $30^{\circ}$ relative to the hinges on the equilateral triangles, the lattice enters a non-topological phase, as shown in fig.4(a). The total number of floppy modes remains the same, but, instead of being all localized on the top edge, they localize on both lattice boundaries with nearly comparable stiffness. Thus, given an external force excitation with amplitude $F$ and frequency $\frac{1}{2}\Delta^{\prime}<\omega<\Delta^{\prime}$, the transmission of any nonlinearly generated second harmonics is reciprocal, because the linear modes driven from opposite lattice boundaries, which are ultimately responsible for second harmonic generation, have the same order of magnitude. In conclusion, through uniform soft twisting modes that allow reversible reconfiguration between topological and non-topological phases, Maxwell lattices have the ability to switch between reciprocal and non-reciprocal transmission regimes of nonlinear elastic waves without the need to physically disassemble and reassemble the system.

VI Discussion and Concluding Remarks

In this paper we have studied the connection between non-reciprocity and topology in Maxwell lattices. Here, the conditions required for the establishment of non-reciprocal behavior come from the interplay between two factors: on one hand, the availability of floppy edge modes, which yield large boundary deformations and trigger a nonlinear response; on the other hand, the topological polarization, which guarantees asymmetry across the lattice.

Different from the previous work by Coulais *et. al.*Coulais et al. (2017), who studied static non-reciprocal elasticity in both of topological and non-topological quasi 1-dimensional (1D) mechanical metamaterials, our design focuses on non-zero frequencies. The concept is developed first using a 1D topological mechanical chain and subsequently generalized to 2D topological kagome latticesKane and Lubensky (2014), paving the way to applications in realistic mechanical metamaterials. The foundational argument of the proposed concept is that topological floppy edge modes produce contrasting stiffness on opposite lattice boundaries. We found that signals with frequency $\omega<\Delta$, where $\Delta$ is the onset of a (partial or total) bandgap, excite edge modes. These modes localize asymmetrically, leading to larger deformation that promotes second harmonic generation at the floppy edge. Second harmonic contributions are bulk modes if $\omega>\frac{1}{2}\Delta$, and, in these conditions, they can propagate through the medium. Finally, because of the asymmetry mentioned above, such transmission is highly non-reciprocal.

By adding an on-site pinning potential to every particle, the band structure of the topological kagome lattice is shifted up by $\Delta^{\prime}$. In these conditions, external signals with $\frac{1}{2}\Delta^{\prime}<\omega<\Delta^{\prime}$ excite linear modes that remain localized at the edges and second harmonic bulk modes that propagate across the lattice. Hence, the second harmonic positive transmission is greater than the fundamental mode transmission and is therefore not trivially overshadowed by the linear response. With this improved configuration, this effect is found to be true for any external harmonic excitation applied at the edges, regardless of the wavenumber established along the boundary. Consequently, the result also holds for point excitations, which represent realistic conditions in practice. We have concluded that the lattice with on-site pinning potential fulfills all the requirements to be labeled a phonon diode. In addition, one can control the geometry of the Maxwell lattice through a Guest mode to switch between topological and non-topological phases. This lattice reconfiguration allows us to manipulate reciprocal and non-reciprocal transmission of elastic waves without disassembling or reassembling the structure.

The idea of nonlinear bulk waves driven by linear edge modes is not limited to second harmonics. One can observe $n$-th order harmonic bulk modes at the output end as long as the input frequency satisfies the condition $\frac{1}{n}\Delta<\omega<\frac{1}{n-1}\Delta$ ($n\geq 2$), while all lower-order harmonics are localized on the boundaries and cannot propagate. A methodological problem to be considered in performing such an extension is associated with the intrinsic limitations of perturbation theory. Since the amplitude of the output signal becomes exponentially small when the order increases, i.e., $A^{(n)}\sim A^{n}$, higher-order harmonic bulk modes become progressively more difficult to be observed. It would thus be interesting to study these phenomena in regimes of strong nonlinearity which invalidate perturbation theory. This kind of study would likely present new challenges arising from the interplay between topological states of matter and nonlinear chaos dynamical theory.

Finally, our investigation so far has been limited to 1D and 2D topological Maxwell lattices. An analogous study of non-reciprocal transmission in 3D topological lattices appears to be possible as a natural extension within the proposed framework. This would open the doors to a broader range of engineering applications and will be one of the next directions in our research.

This work was supported by the National Science Foundation (Grant No. NSF-EFRI-1741618).

Author contributions: D. Z. performed the theoretical analysis and computations. All authors contributed to the formulation of the problem, analysis of results, and preparation of the manuscript.

Appendix A Reciprocity of linear elastic systems with time reversal symmetry

With time reversal symmetry, the transmission of linear elastic modes is reciprocal, meaning that the transmission susceptibilities from point $A$ to point $B$ and from point $B$ to point $A$ are equal. This is the essence of the Maxwell-Betti’s theoremMaxwell (1864); Betti (1872); Charlton (1960). In this section we verify this theorem by considering a $d$-dimensional general lattice based on spring-mass interactions. Within linear elasticity, the Newton’s equation of motion is

[TABLE]

where $n$ denotes a lattice site, $\eta$ is the damping coefficient, $V$ is the lattice potential energy, $\vec{\nabla}_{n}=\sum_{i=1}^{d}\hat{e}_{i}\partial_{u_{n}^{(i)}}$, and $\vec{F}_{n}^{\rm ext}$ is the external driving force. We rewrite the displacement field as a $Nd$-dimensional vector ${\textbf{u}}=(\vec{u}_{1},\vec{u}_{2},...,\vec{u}_{N})$, and rewrite the external driving as a $Nd$-dimensional vector ${\textbf{F}}^{\rm ext}$. The linear elastic mode can be calculated as

[TABLE]

where ${\textbf{G}}(\omega)=\left[{\textbf{D}}+(-m\omega^{2}+i\eta\omega){\textbf{I}}\right]^{-1}$ is the frequency response function, and D is the dynamical matrix. By using an orthogonal transformation S that relates $u_{n}$ to the normal modes $u_{\alpha}$ through $u_{\alpha}=\sum_{n}S_{\alpha n}u_{n}$, we can express the normal modes as follows,

[TABLE]

where $G_{\alpha}(\omega)=[\epsilon_{\alpha}+(-m\omega^{2}+i\eta\omega)]^{-1}$, and $\epsilon_{\alpha}$ is the $\alpha^{th}$ eigenvalue of the dynamical matrix D. We plug in the driving force at point $A$ to calculate the displacement at $B$, with $u_{B}(\omega)=\sum_{\alpha}S_{\alpha B}G_{\alpha}(\omega)S_{\alpha A}F_{A}^{\rm ext}(\omega)$. Similarly, we plug in the driving force at $B$ to calculate displacement at $A$ with $u_{A}(\omega)=\sum_{\alpha}S_{\alpha A}G_{\alpha}(\omega)S_{\alpha B}F_{B}^{\rm ext}(\omega)$. It is evident that $u_{A}(\omega)/F_{B}^{\rm ext}(\omega)=u_{B}(\omega)/F_{A}^{\rm ext}(\omega)$, meaning that the transmission of linear modes is reciprocal in real space in any elastic system with time reversal symmetry.

Appendix B Analytical calculation of linear and second harmonic modes in 1D topological mechanical chain

As shown in fig.1(a), the 1D topological mechanical chain consists of rigid bars of length $r$, free to rotate about hinges separated by the distance $a$, creating repeated 2-site unit cell of length $2a$. A mass point $m$ is attached to the end of each bar, and the neighboring ends are connected by harmonic springs with spring constant $K$. The equilibrium configuration is such that each rotor makes an angle $\bar{\theta}$ relative to the upward or downward normals. The angular displacement of rotor $n$ is $u_{n}=r\delta\theta_{n}$. We label the rotors from 0 to $N$, with open boundary conditions at rotors 0 and $N$. We apply an external angular driving along the tangential direction of rotor $g$ with $F_{n}^{\rm ext}(t)=Fe^{i\omega t}\delta_{ng}$. $F$ is not large, making angular displacements of all rotors $u_{n}\ll r$, which further validates perturbation theory. It is convenient to rewrite the external force as $F\to\lambda F$, where $\lambda\ll 1$. The full Newton’s equation of motion is

[TABLE]

where the lower index $g$ indicates that the elastic mode is in response to the external driving at rotor $g$. The open boundary conditions at rotor 0 and rotor $N$ are given by

[TABLE]

where $\vec{T}_{n}=T_{n}\hat{n}_{n}$ is the tension in bond $n$ connecting sites $n$ and $n+1$. $\hat{n}_{n}$ is the unit vector of bond $n$, and $\hat{t}_{n}$ is the tangential unit vector of rotor $n$. We expand the tangential component of bond tension in orders of $u_{n}/r$, denoted by $(\vec{T}_{n-1}-\vec{T}_{n})\cdot\hat{t}_{n}=f_{n}^{(1)}+f_{n}^{(2)}+\mathcal{O}(u^{3})$. The leading order is

[TABLE]

with

[TABLE]

where $\bar{\theta}>0$ is assumed in this paper, leading to $c_{1}>c_{2}$. The second order term is

[TABLE]

with $C_{1,2,3,4,5,6}$ being constants determined by the geometric parameters of the chain. In our calculations, we choose $a=2r$ and $\theta=\pi/4$. The coefficients are given by

[TABLE]

We denote $\textbf{u}=(u_{0},u_{1},...,u_{N})$ as the angular displacement of the particles, and expand it in orders of $\lambda$, with $\textbf{u}=\textbf{u}^{(1)}+\textbf{u}^{(2)}+\mathcal{O}(\lambda^{3})$, where $\textbf{u}^{(1)}$ and $\textbf{u}^{(2)}$ are linear and second harmonic modes, respectively. We further denote $\textbf{F}^{\rm ext}=(\lambda F,0,...,0)^{T}$ ($\textbf{F}^{\rm ext}=(0,0,...,\lambda F)$) as the external driving force driven at rotor $g=0$ ($g=N$), and denote $\textbf{f}^{(2)}(\textbf{u}^{(1)})=(f_{0}^{(2)},f_{1}^{(2)},...,f_{N}^{(2)})$ as the second harmonic effective feedback force generated by the linear elastice modes. By expanding Eq.(6) up to the second order of $\lambda$, one obtains

[TABLE]

where $\textbf{D}=K[(c_{1}^{2}+c_{2}^{2})\delta_{ij}-c_{1}c_{2}(\delta_{i,j+1}+\delta_{i,j-1})]$ is the dynamical matrix. By matching the equations of motion in orders of $\lambda$, we solve for the linear mode

[TABLE]

and for the second harmonic mode

[TABLE]

subjected to the open boundary conditions at rotors [math] and $N$. The frequency response function of the 1D chain is

[TABLE]

The dispersion relation of the only bulk phonon mode is $\omega(k)=[(c_{1}-c_{2})^{2}+4c_{1}c_{2}\sin^{2}\left(ak/2\right)]^{1/2}$, where $k$ is the wave number. The band has lower limit $\Delta=|c_{1}-c_{2}|$ and upper limit $\Delta^{\prime}=|c_{1}+c_{2}|$. Thus, the linear mode is a bulk mode if $\Delta<\omega<\Delta^{\prime}$, while it is an edge mode if $\omega<\Delta$. According to Eq.(13), the generic solution for a linear mode is

[TABLE]

where $\lambda_{1,2}$ are given by

[TABLE]

with $\gamma=\frac{\omega^{2}}{c_{1}c_{2}}-\frac{c_{1}}{c_{2}}-\frac{c_{2}}{c_{1}}-\frac{i\eta\omega}{c_{1}c_{2}}$. Through Eq.(13) and Eq.(B) we can solve for $a_{g},A_{g},b_{g},B_{g}$. We let $g=0$ and $g=N$ to obtain the linear modes when the external driving is applied in at rotors [math] and $N$. Given the definition of local response function, $\chi_{\rm in}(\omega)=|u_{\rm in}^{(1)}(\omega)|/F_{\rm in}(\omega)$ at rotor [math] (rigid end, $\chi_{{\rm in}-}$) and rotor $N$ (soft end, $\chi_{{\rm in}+}$), the local response functions are given by

[TABLE]

where we have used $\lim_{N\to\infty}|\lambda_{2}/\lambda_{1}|^{N}=0$ when $\omega<c_{1}-c_{2}$. The ratio $\chi_{{\rm in}+}/\chi_{{\rm in}-}$ tells which end of the 1D chain has greater displacement in response to external loading. In the limit of $\eta\to 0$, analytical calculations reveal that

[TABLE]

where $\omega^{*}=(c_{1}^{2}-c_{2}^{2})/\sqrt{2(c_{1}^{2}+c_{2}^{2})}$. As long as $\omega<\omega^{*}$, the response linear edge mode of the right end is greater than that of the left side. The linear transmission susceptibility, defined by $\chi_{\rm out}^{(1)}(\omega)=|u_{\rm out}^{(1)}(\omega)|/F(\omega)$, is given by

[TABLE]

The fact that the linear transmission susceptibilities are equal is testament to the reciprocity of linear waves.

We further study the second harmonic modes based on Eq.(14). To calculate the second harmonic mode displacement $\textbf{u}_{g}^{(2)}$, we notice that $\textbf{u}_{g}^{(2)}=\sum_{g^{\prime}=0}^{N}\textbf{u}_{g^{\prime}g}^{(2)}$, where $\textbf{u}_{g^{\prime}g}^{(2)}$ is the displacement field of the chain in response to external force $f^{(2)}_{g^{\prime}g}$ applied at a single rotor $g^{\prime}$. The displacement field is given by

[TABLE]

where we denote the external driving at rotor $g^{\prime}$ as a $N\times 1$ vector $\textbf{f}^{(2)}_{g^{\prime}g}=(0,0,...,f^{(2)}_{g^{\prime}g},...,0)$. The generic solution of Eq.(23) is of the following form

[TABLE]

where $\mu_{1,2}$ satisfy

[TABLE]

with $\nu=\frac{4\omega^{2}}{c_{1}c_{2}}-\frac{c_{1}}{c_{2}}-\frac{c_{2}}{c_{1}}-\frac{2i\eta\omega}{c_{1}c_{2}}$. Through Eq.(23) and Eq.(B) we solve for $c_{g^{\prime}g},C_{g^{\prime}g},d_{g^{\prime}g},D_{g^{\prime}g}$ to determine $\textbf{u}^{(2)}_{g^{\prime}g}(2\omega)$. Finally, given the definition of second harmonic transmission susceptibility, $\chi_{\rm out}^{(2)}=|u_{\rm out}^{(2)}(2\omega)|/F(\omega)$, we obtain

[TABLE]

The second harmonic transmission susceptibility is determined by $\textbf{f}_{g^{\prime}g}^{(2)}$, which in turn is proportional to the square of linear elastic waves amplitude. This consideration is essential in explaining how the asymmetric local response function $\chi_{\rm in+}^{(1)}\gg\chi_{\rm in-}^{(1)}$ results in the non-reciprocal transmission of second harmonic modes ($\chi_{\rm out+}^{(2)}\gg\chi_{\rm out-}^{(2)}$).

Appendix C Analytical calculation of linear and second harmonic modes in a 2D topological kagome lattice

In this section we calculate linear and second harmonic modes in a 2D topological kagome lattice. The unit cell is shown in fig.2(a) and consists 6 bonds, with rest lengths $l_{i}$ and unit vector directions $\hat{n}_{i}=(\cos\theta_{i},\sin\theta_{i})$, $i=1,2,...,6$. We define the $2\times 2$ “dynamical matrix” of bond $i$, as

[TABLE]

$\vec{a}_{1}$ and $\vec{a}_{2}$ are primitive vectors. The lattice has periodic boundary condition in $\vec{a}_{1}$ and open boundary condition in $\vec{a}_{2}$, leaving the top and bottom boundaries open.

The compatibility matrix of the quasi-1D strip of deformed kagome lattice is $\textbf{C}_{ij}(k_{1})=\textbf{C}_{1}(k_{1})\delta_{ij}+\textbf{C}_{2}(k_{1})\delta_{i+1,j}$, where $\textbf{C}_{1}(k_{1})$ and $\textbf{C}_{2}(k_{1})$ are intra-cell and inter-cell compatibility matrices, respectively:

[TABLE]

We simplify the lattice to a quasi-1D strip by applying Bloch condition along $\vec{a}_{1}$ to obtain a finite supercell strip with $N_{2}$ unit cells. Thus, the elastic wave $\textbf{u}_{n_{1}n_{2}}=e^{ik_{1}n_{1}}\textbf{u}_{n_{2}}$, with $k_{1}=\vec{k}\cdot\vec{a}_{1}$. We denote the displacement of cell $n_{2}$ as $\textbf{u}_{n_{2}}=(u_{n_{2}A}^{x},u_{n_{2}A}^{y},u_{n_{2}B}^{x},u_{n_{2}B}^{y},u_{n_{2}C}^{x},u_{n_{2}C}^{y})$, and further denote the displacement of the entire strip as $\textbf{u}=(\textbf{u}_{1},\textbf{u}_{2},...,\textbf{u}_{N_{2}})$. To fully gap the band structure of the lattice, we introduce an on-site potential $\frac{1}{2}K^{\prime}\textbf{u}\textbf{u}^{T}$ by embedding the lattice on a soft substrate.

Before calculating the elastic waves, we first derive the tension of a bond connecting two sites $i$ and $j$. The bond rest length is $l_{0}$ and it’s unit vector is $\hat{n}=(\cos\theta,\sin\theta)$. We denote the relative displacement of the bond as $\Delta\vec{u}_{ij}=\vec{u}_{j}-\vec{u}_{i}$, and expand the tension in orders of it, with $\vec{T}_{i}=\vec{F}_{i}+\vec{f}_{i}+\mathcal{O}(\Delta u_{ij}^{3})$. The leading order term is

[TABLE]

and the second order ones are

[TABLE]

[TABLE]

Eqs.(C, C) are useful to derive the effective second harmonic feedback force.

We now apply an external driving at cell $g$ of the strip, denoted as ${F}_{n_{1}n_{2}}^{\rm ext}(t)=e^{i\omega t-ik_{1}n_{1}}{F}_{C}^{y}\delta_{n_{2}g}$. Thus, we further denote the external force as a $6N_{2}$-dimensional vector $\textbf{F}^{\rm ext}_{g}=(\textbf{F}^{{\rm ext}}_{1},\textbf{F}^{{\rm ext}}_{2},...,\textbf{F}^{{\rm ext}}_{N_{2}})$. We assume $F_{C}^{y}$ is not too large, rendering all displacements $|\vec{u}_{A,B,C}|\ll l_{i=1,2,...,6}$, which further validates perturbation theory. Thus, we expand u in orders of $F_{C}^{y}$, as $\textbf{u}=\textbf{u}^{(1)}+\textbf{u}^{(2)}+\mathcal{O}(F^{3})$, where $\textbf{u}^{(1)}$ and $\textbf{u}^{(2)}$ are linear and 2nd harmonic modes. We define the frequency response function as

[TABLE]

where $\textbf{D}(k_{1})=K\textbf{C}^{\dagger}(k_{1})\textbf{C}(k_{1})$ is the dynamical matrix of the supercell strip. Thus, the linear mode is

[TABLE]

The 2nd harmonic mode is

[TABLE]

where $\textbf{f}^{(2)}(\textbf{u}_{g}^{(1)})=(\textbf{f}_{1},\textbf{f}_{2},...,\textbf{f}_{N_{2}})$ is the second harmonic effective feedback force generated by the linear mode $\textbf{u}_{g}^{(1)}$, as shown in Eqs.(C, C).

The generic form of linear mode $\textbf{u}_{g}^{(1)}$ is given as follows:

[TABLE]

where $\lambda_{\alpha}$, $\alpha=1,2,3,4$ are the eigenvalues of $[\textbf{D}(k_{1},\lambda)+(-m\omega^{2}+K^{\prime}+i\eta\omega)\textbf{I}]$, and $\bm{\phi}_{\alpha}$ are the corresponding $6\times 1$ eigenvectors. $\textbf{D}(k_{1},\lambda)$ is the following $6\times 6$ matrix,

[TABLE]

$a_{g\alpha}$ and $b_{g\alpha}$, $\alpha=1,2,3,4$ are constants determined by the open boundary condition at cell $1$,

[TABLE]

and the open boundary condition at cell $N_{2}$,

[TABLE]

Together with Eq.(45), we solve $a_{g\alpha}$ and $b_{g\alpha}$ to obtain the linear mode $\textbf{u}_{g}^{(1)}$.

We then calculate the second harmonic mode based on Eq.(46). However, it is not easy to solve $\textbf{u}_{g}^{(2)}$, because the effective second harmonic feedback force, $\textbf{f}^{(2)}_{g}(2\omega,2k_{1})$ is applied at every cell. In order to simplify this problem, we consider the mode $\textbf{u}_{g^{\prime}g}^{(2)}(2\omega,2k_{1})$ in response to the second harmonic effective feedback force $\textbf{f}_{g^{\prime}g}^{(2)}$ applied at a single cell $g^{\prime}$,

[TABLE]

The generic form of $\textbf{u}^{(2)}_{g^{\prime}g}(2\omega,2k_{1})$ is given by

[TABLE]

where $\mu_{\beta}$, $\beta=1,2,3,4$ are the eigenvalues of $[\textbf{D}(2k_{1},\mu)+(-4m\omega^{2}+K^{\prime}+2i\eta\omega)\textbf{I}]$, and $\bm{\varphi}_{\beta}$ are the corresponding $6\times 1$ eigenvectors. $c_{g^{\prime}g\beta}$ and $d_{g^{\prime}g\beta}$ are constants determined by the boundary condition at cell $1$

[TABLE]

and open boundary condition at cell $N_{2}$

[TABLE]

Together with Eq.(54), we solve $c_{g^{\prime}g\beta}$ and $d_{g^{\prime}g\beta}$ to obtain $\textbf{u}^{(2)}_{g^{\prime}g}$. Finally, the second harmonic displacement is $\textbf{u}^{(2)}_{g}=\sum_{g^{\prime}=1}^{N_{2}}\textbf{u}^{(2)}_{g^{\prime}g}$.

Appendix D Numerical simulation of input local response function and output transmission susceptibility

In the simulation, a finite kagome lattice which spans $40\vec{a}_{1}\times(40-1)\vec{a}_{2}$ area in real space, is considered. The lattice is made of $40\times(40-1)$ unit cells with an additional layer of C-sites at the bottom. We connect the leftmost particles to the rightmost ones with harmonic springs to provide periodic boundary condition in the $\vec{a}_{1}$ direction, and we leave the C-sites of top and bottom boundaries free to realize open boundary conditions. By applying a vertical harmonic force $\vec{F}=(0,F)$ either on the top edge or bottom edge C-sites, we drive the lattice and we compute the displacement of mass points using a molecular dynamics scheme with damping. In order to make sure steady-state conditions are established, we wait $400\times(2\pi/\omega)$ before we make any displacement reading. We collect displacements $\vec{u}_{1C}(t)$ and $\vec{u}_{N_{2}C}(t)$ on the two edges. By applying Fast Fourier Transformation (FFT), we convert displacement time histories into their frequency spectra, $\vec{u}_{1C}(\omega)$ and $\vec{u}_{N_{2}C}(\omega)$. The elastic response is obtained via summation of multiple modes,

[TABLE]

where $\vec{u}_{n_{2}C}^{(1)}(\omega)=(u_{n_{2}C}^{x,(1)}(\omega),u_{n_{2}C}^{y,(1)}(\omega))$ and $\vec{u}_{n_{2}C}^{(2)}(2\omega)=(u_{n_{2}C}^{x,(2)}(2\omega),u_{n_{2}C}^{y,(2)}(2\omega))$. $u_{n_{2}C}^{x,(1)}(\omega)$, $u_{n_{2}C}^{y,(1)}(\omega)$, $u_{n_{2}C}^{x,(2)}(2\omega)$ and $u_{n_{2}C}^{y,(2)}(2\omega)$ are the amplitudes of $x$ and $y$ components of linear and second harmonic modes. The input linear response function is defined as $\chi_{\rm in+}=|u_{1C}^{y,(1)}(\omega)|/F$, and the output linear transmission susceptibility is defined as $\chi_{\rm out+}^{(1)}=|u_{N_{2}C}^{y,(1)}(\omega)|/F$. The output second harmonic transmission susceptibility is calculated through $\chi_{\rm out+}^{(2)}=|\vec{u}_{N_{2}C}^{(2)}(2\omega)|/F$.