Non-Reciprocal MEMS Periodic Structure

Abstract

In recent years, active periodic structures with in-time modulated parameters have drawn ever-increasing attention due to their peculiar (and sometimes exotic) wave propagation properties. Although many experimental works have shown the efficacy of time-modulation strategies, the benchmarks proposed until now have been mostly proof-of-concept demonstrators, with little attention to the feasibility of the solution for practical purposes. In this work, we propose a micro electro-mechanical system (MEMS) periodic structure with modulated electromechanical stiffness featuring non-reciprocal band-gaps that are frequency bands where elastic waves are allowed to travel only in one direction. To this aim, we derive a simplified analytical lumped-parameter model, which is then verified through numerical simulations of both the lumped-parameter system and the high-fidelity multiphysics finite element model including electrostatic effects. We envision that this system, which can easily be manufactured through standard MEMS production processes, may be used as a directional filter in MEMS devices such as insulators and circulators.

1. Introduction

In mechanics, periodic structures are systems composed of the repetition of a unitary cell with a well-defined topology and physical properties [1]. In the last few decades, they have drawn an ever-increasing interest due to their ability to manipulate the propagation of waves in a number of different ways [2]. From a mechanical standpoint, vibration insulation and filtering are possibly their most interesting applications. Periodic structures indeed feature frequency bands, known as band-gaps, where wave propagation is hindered.

In recent years, a wealth of works dealing with periodic structures with (periodically) time-varying properties have shown that even more complex phenomena can be attained: selective wave filtering [3], topological pumping [4], and non-reciprocal wave transmission [5] are a few examples. The latter consists, loosely speaking, of the fact that wave propagation in a non-reciprocal medium is different depending on the propagation direction. In the context of periodic structures, this reveals in the presence of the so-called directional band-gaps and in a non-symmetric dispersion diagram [5,6,7,8,9]. In other words, depending on the propagation direction, one may find band-gaps that are not present if traveling in a different direction. In the case of 1D periodically modulated waveguides, directional band-gaps for right and left-propagating waves will appear in different frequency intervals, thus realizing a directional filter. Many benchmarks have already been proposed in the literature, as in [8], where an array of masses coupled through piloted magnets have been studied, and in [7], where an elastic beam with periodic shunted piezoelectric patches was experimentally tested. These solutions, however, are mainly proof-of-concept demonstrators.

In this work, we propose a micro electro-mechanical system (MEMS) periodic structure that can be used as a directional filter. In Radio Frequency (RF) applications, indeed, non-reciprocal components have long been used for filters, insulators, and circulators [10,11]. Given the high demand for these components in modern telecommunication systems, MEMS filters [12] and non-reciprocal filters [13,14] have also been proposed to reduce the footprint area of these devices and lower the production costs. To the best of the authors’ knowledge, however, all these implementations rely upon different kinds of electric circuit resonators and typically operate in the range of MHz or even GHz. The solution we propose in this work, instead, uses a mechanical array of resonators, whose equivalent mechanical properties are modulated by means of the electrostatic forces at the interface of properly placed parallel plate electrodes. For this reason, the layout we propose can even be designed to work at frequencies in the range of a few tens of kHz.

The paper is organized as follows. In Section 2, a simplified electromechanic model for the unitary cell of the system is derived; in Section 3, with the aid of a continuum model, we find analytic expressions for the key parameters defining the non-symmetric dispersion diagram of the structure; in Section 4, we numerically verify our model by means of a high-fidelity multiphysics finite element model. Finally, conclusions are drawn in Section 5.

2. The Model

The system we propose is a 1-dimensional periodic waveguide, basically consisting of a spring-mass chain with elastic suspensions and where each mass hosts a fixed electrode in a dedicated cavity. Like the majority of inertial MEMS, the system features a 2D geometry, which is imposed by the etching process of the silicon wafer. In this context, springs are typically realized using folded beams [15]. The unitary cell of the periodic structure is in turn composed of $R=3$ sub-cells, with three being the minimum number of sub-cells to obtain a wave-like traveling stiffness profile required to observe a non-reciprocal behavior [16]. The stiffness modulation, as already mentioned, is achieved by regulating the electrostatic forces acting between the fixed electrodes and the masses. The system is qualitatively represented in Figure 1.

Let us now consider the sub-cell depicted in Figure 2, where the mass m is connected to the following/preceding mass and to the ground by springs of elastic constant k and $k_g$, respectively. For the i-th mass, the equation of motion writes

$$m{\ddot{x}}_i+(2k+k_g)x_i-k{x}_{i-1}-k{x}_{i+1}-f_{el}=0,$$

(1)

where $f_{el}$ is the net electrostatic force acting on the mass. Considering only the right side of the electrode, we can write its capacity as

$$C_R\left(x\right)=\frac{\epsilon A}{x+x_0},$$

(2)

where $\epsilon$ is the vacuum permittivity constant, A the area, and $x_0$ the gap at rest. The electromechanical attractive force can then be obtained from the capacitor potential energy U as

$$f_{el,R}=\frac{\partial U}{\partial x}=\frac{1}{2}\frac{\partial C_R}{\partial x}\Delta V^2=-\frac{1}{2}\frac{\epsilon A}{{(x+x_0)}^2}\Delta V^2,$$

(3)

with $\Delta V={\tilde{V}}_i-V_0$, ${\tilde{V}}_i$ is the time-varying voltage applied to the electrode, and $V_0$ the constant voltage applied to the structure. Using a Taylor expansion about $x=0$ and neglecting $\mathcal{O}\left(x^4\right)$ terms, we can write

$$f_{el,R}\approx -\frac{\epsilon A}{2}\left(\frac{1}{x_0^2}-\frac{2}{x_0^3}x+\frac{3}{x_0^4}{x}^2-\frac{4}{x_0^5}{x}^3\right)\Delta V^2.$$

(4)

In the same way, for the left side, we have

$$f_{el,L}=\frac{1}{2}\frac{\epsilon A}{{(x-x_0)}^2}\Delta V^2\approx \frac{\epsilon A}{2}\left(\frac{1}{x_0^2}+\frac{2}{x_0^3}x+\frac{3}{x_0^4}{x}^2+\frac{4}{x_0^5}{x}^3\right)\Delta V^2,$$

(5)

so that, summing the two contributions, we have that the net electrostatic force is

$$f_{el}=f_{el,L}-f_{el,R}\approx \epsilon A\left(\frac{2}{x_0^3}x+\frac{4}{x_0^5}{x}^3\right)\Delta V^2.$$

(6)

As it can be observed, this configuration is self-equilibrated (electrostatic forces are null for $x=0$) and extends the linearity range of the parallel plate capacitor since the quadratic terms are canceled out. Additionally, it can be seen how $f_{el}$ introduces a linear and cubic voltage-dependant stiffness into the system. Setting the electrode tension to

$${\tilde{V}}_i\left(t\right)=V_{a}cos\left({\omega }_{m}t-\frac{2i}{R}\pi \right)=V_{a}cos{\theta }_i,$$

(7)

where $i=\{1,\dots ,R\}$ is the sub-cell ordinal, $R=3$, and ${\omega }_m$ is the modulation frequency, we can develop the squared voltage difference as

$$\Delta V^2=\left(\frac{V_a^2}{2}+V_0^2\right)-2V_a{V}_{0}cos{\theta }_i+\frac{V_a^2}{2}cos\left(2{\theta }_i\right).$$

(8)

Plugging (8) into (6), we obtain

$$f_{el}=\left(k_{E0}-k_{E1}cos{\theta }_i+k_{E2}cos\left(2{\theta }_i\right)\right)\left(1+\frac{2}{x_0^2}{x}^2\right)x,$$

(9)

where

$$\begin{array}{c}\hfill k_{E0}=\frac{2\epsilon A}{x_0^3}\left(\frac{V_a^2}{2}+V_0^2\right),k_{E1}=\frac{4\epsilon A}{x_0^3}{V}_a{V}_0,k_{E2}=\frac{\epsilon A}{x_0^3}{V}_a^2.\end{array}$$

(10)

According to (9), the electrostatic stiffness features a constant term and two harmonics at ${\omega }_m$ and $2{\omega }_m$ (for both the linear and cubic terms). However, considering typical values $V_0\sim 10$ V and $V_a\sim 1$ V, we have that

$$\frac{k_{E1}}{k_{E2}}=\frac{4V_0}{V_1}\sim 40,$$

(11)

and we can neglect the second harmonic.

Under the hypotheses of small displacements, we can also neglect the cubic stiffness, and, substituting the electrostatic forces in the equation of motion, we can finally write

$$m{\ddot{x}}_i+(2k+k_g-k_{E0}+k_{E1}cos{\theta }_i)x_i-k{x}_{i-1}-k{x}_{i+1}=0$$

(12)

and recast the equations of motion for the i-th cell in matrix form as

$$\mathbf{M}{\ddot{\mathbf{x}}}_i+\mathbf{K}\left(t\right){\mathbf{x}}_i-{\mathbf{K}}_L{\mathbf{x}}_{i-1}-{\mathbf{K}}_R{\mathbf{x}}_{i+1}=\mathbf{0},$$

(13)

where ${\mathbf{x}}_i$ is the displacement vector of the i-th cell, ${\mathbf{x}}_{i-1}$ and ${\mathbf{x}}_{i+1}$ are the neighboring cell displacement vectors, while ${\mathbf{K}}_{L/R}$ are the corresponding left/right stiffness matrices, so that Equation (13) writes:

$$\begin{array}{c}\left[\begin{array}{ccc}m& 0& 0\\ 0& m& 0\\ 0& 0& m\end{array}\right]\left\{\begin{array}{c}{\ddot{x}}_1\\ {\ddot{x}}_2\\ {\ddot{x}}_3\end{array}\right\}-\left[\begin{array}{ccc}0& 0& k\\ 0& 0& 0\\ 0& 0& 0\end{array}\right]\left\{\begin{array}{c}{x}_{-2}\\ x_{-1}\\ x_0\end{array}\right\}-\left[\begin{array}{ccc}0& 0& 0\\ 0& 0& 0\\ k& 0& 0\end{array}\right]\left\{\begin{array}{c}{x}_4\\ x_5\\ x_6\end{array}\right\}+\hfill \\ \hfill \left[\begin{array}{ccc}{k}^{*}+k_{E1}cos{\theta }_1\left(t\right)& 0& 0\\ 0& k^{*}+k_{E1}cos{\theta }_2\left(t\right)& 0\\ 0& 0& k^{*}+k_{E1}cos{\theta }_3\left(t\right)\end{array}\right]\left\{\begin{array}{c}{x}_1\\ x_2\\ x_3\end{array}\right\}=0\end{array}$$

(14)

where $k^{*}=2k+k_g-k_{E0}$. Using this formalism, the asymmetric dispersion relation for the lumped-parameter spatiotemporal cell can be computed with the Bloch Based Procedure (BBP) illustrated in [6].

3. Parametric Analysis for Harmonic Modulation

In this section, we derive the modulation parameters that can be used to design a structure satisfying a user-tailored dispersion relation. This can be performed, in theory, by analyzing the quadratic eigenvalue problem (QEVP) stemming from the application of the BBP to (13). However, even considering the simple case (described until now) of three masses and truncating the solution to the first harmonic, the resulting equations feature high-order polynomials with no closed-form solution. To circumvent this problem, we turn to the analysis of a rod with elastic suspensions, which is the continuum equivalent of a spring-mass chain system, and we seek a solution using the well-known Plane Wave Expansion Method (PWEM) [17].

3.1. Dispersion Relation for Rod on Elastic Foundations

Considering a rod with Young’s modulus E, density $\rho$, and cross-section area A, the governing equation writes

$$EA\frac{{\partial }^{2}u}{\partial x^2}-\gamma u=\rho A\frac{{\partial }^{2}u}{\partial t^2},$$

(15)

where $\gamma$ is the distributed ground stiffness and is modulated as

$$\gamma (x,t)={\gamma }_0+{\gamma }_{m}cos({\omega }_{m}t-{\kappa }_{m}x),$$

(16)

where ${\gamma }_0$ and ${\gamma }_m$ are constant, and ${\kappa }_m=2\pi /{\lambda }_m$ is the modulation wavenumber, with ${\lambda }_m$ as the spatial periodicity of the spatiotemporal cell. To apply PWEM, both Floquet solution u and stiffness $\gamma$ are expanded in a Fourier series as

$$u=\sum _{n=-\infty }^{+\infty }{\widehat{u}}_n{e}^{i((\omega +n{\omega }_m)t-(\kappa +n{\kappa }_{m}x))},\gamma =\sum _{p=-\infty }^{+\infty }{\widehat{\gamma }}_p{e}^{ip({\omega }_{m}t-{\kappa }_{m}x)},$$

(17)

which are first truncated to order N and then plugged into (15). Exploiting the orthogonality property of the Fourier series, one can rid of the double summation generated by the $\gamma u$ term, and finally, write a set of equations

$$EA{(\kappa +q{\kappa }_m)}^2{\widehat{u}}_q+\sum _{n=-N}^N{\widehat{\gamma }}_{q-n}{\widehat{u}}_n=\rho A{(\omega +q{\omega }_m)}^2{\widehat{u}}_q,$$

(18)

where $q\in \mathbb{Z}$ is an additional index. Writing a number $2N+1$ of the above equation for $q\in \{-N,\dots ,N\}$, a square system in the coefficients ${\widehat{u}}_q$ is obtained and, in matrix form, writes as

$$\left({\mathbf{L}}_0\left(\kappa \right)+\omega {\mathbf{L}}_1+{\omega }^2{\mathbf{L}}_2\right)\widehat{\mathbf{u}}=\mathbf{L}(\kappa ,\omega )\widehat{\mathbf{u}}=0,$$

(19)

which is a QEVP that can be solved by imposing $\kappa$, and where $\widehat{\mathbf{u}}={[{\widehat{u}}_N,\dots ,{\widehat{u}}_{-N}]}^T$. For each $\kappa$, then $2(2N+1)$ eigenfrequencies are obtained and the dispersion diagram can be drawn.

3.2. Dispersion Analysis

To derive simplified formulas for the dispersion key characteristics, namely, directional bandgap positions and amplitude, with respect to the modulation parameters, we loosely retrace the steps presented in [5] for a simple rod with no elastic suspensions.

First, setting $N=1$, we retain only the first harmonic in the solution. We also notice that for the harmonic modulation in (16), the Fourier coefficients are ${\widehat{\gamma }}_0={\gamma }_0$, ${\widehat{\gamma }}_{\pm 1}={\gamma }_m/2$, and ${\gamma }_p=0$ for $\left|p\right|\ge 2$. Let us define

$${\omega }_{\gamma }=\sqrt{\frac{{\widehat{\gamma }}_0}{\rho A}},{\omega }_a=\sqrt{\frac{{\widehat{\gamma }}_1}{\rho A}},$$

(20)

and the following dimensionless quantities

$${\mathrm{\Omega }}_{★}=\frac{{\omega }_{★}{\lambda }_m}{2\pi c_0},\mu =\kappa {\lambda }_m,$$

(21)

where ★ denotes a generic subscript and where $c_0=\sqrt{E/\rho }$. Using the expressions above, we can write a $3\times 3$ dimensionless QEVP as

$$\mathbf{L}(\mu ,\mathrm{\Omega })=0,$$

(22)

where the components of L are

$$L_{11}={\left(\frac{\mu }{2\pi }+1\right)}^2-{(\mathrm{\Omega }+{\mathrm{\Omega }}_m)}^2+{\mathrm{\Omega }}_{\gamma }^2,$$

(23)

$$L_{22}={\left(\frac{\mu }{2\pi }\right)}^2-{\mathrm{\Omega }}^2+{\mathrm{\Omega }}_{\gamma }^2,$$

(24)

$$L_{33}={\left(\frac{\mu }{2\pi }-1\right)}^2-{(\mathrm{\Omega }-{\mathrm{\Omega }}_m)}^2+{\mathrm{\Omega }}_{\gamma }^2,$$

(25)

$$L_{12}=L_{21}=L_{23}=L_{32}={\mathrm{\Omega }}_a^2.$$

(26)

As conducted in [5], we can focus on the forward and backward directional band-gaps separately by considering the partitions of $\mathbf{L}(\mu ,\mathrm{\Omega })$ pertaining to $\{{\widehat{u}}_0,{\widehat{u}}_{-1}\}$ and $\{{\widehat{u}}_0,{\widehat{u}}_{+1}\}$, respectively. This way, the characteristic equation for the backward directional band-gap writes

$$L_{11}{L}_{22}-L_{12}{L}_{21}=0,$$

(27)

that is,

$$\left({\left(\frac{\mu }{2\pi }+1\right)}^2-{(\mathrm{\Omega }+{\mathrm{\Omega }}_m)}^2+{\mathrm{\Omega }}_{\gamma }^2\right)\left({\left(\frac{\mu }{2\pi }\right)}^2-{\mathrm{\Omega }}^2+{\mathrm{\Omega }}_{\gamma }^2\right)=-{\mathrm{\Omega }}_a^2.$$

(28)

The above equation is quartic; however, in the limit ${\mathrm{\Omega }}_a\to 0$ (no modulation), it turns into two separate quadratic algebraic equations, whose solutions are

$${\mathrm{\Omega }}_{1,2}=\pm \sqrt{{\mathrm{\Omega }}_{\gamma }^2+{\left(\frac{\mu }{2\pi }\right)}^2},$$

(29)

$${\mathrm{\Omega }}_{3,4}=-{\mathrm{\Omega }}_m\pm \sqrt{{\mathrm{\Omega }}_{\gamma }^2+{\left(\frac{\mu }{2\pi }+1\right)}^2}.$$

(30)

Similarly, from $L_{22}{L}_{33}=0$, we obtain

$${\mathrm{\Omega }}_{5,6}=+{\mathrm{\Omega }}_m\pm \sqrt{{\mathrm{\Omega }}_{\gamma }^2+{\left(\frac{\mu }{2\pi }-1\right)}^2}.$$

(31)

We can now obtain the wave number where the branches intersect by solving ${\mathrm{\Omega }}_1={\mathrm{\Omega }}_3$ and ${\mathrm{\Omega }}_1={\mathrm{\Omega }}_5$ for $\mu$:

$${\mu }_{F,B}=\pm \pi (1\pm {\mathrm{\Omega }}_m\sigma ),$$

(32)

where

$$\sigma =\sqrt{1+\frac{4{\mathrm{\Omega }}_{\gamma }^2}{1-{\mathrm{\Omega }}_m^2}}.$$

(33)

The corresponding frequencies, finally, write

$${\mathrm{\Omega }}_{F,B}={\mathrm{\Omega }}_1\left({\mu }_{F/B}\right)=\sqrt{{\mathrm{\Omega }}_{\gamma }^2+\frac{{(1\pm {\mathrm{\Omega }}_m\sigma )}^2}{4}},$$

(34)

so that the positions of the directional band-gaps openings are defined by the two sets $\{{\mu }_B,{\mathrm{\Omega }}_B\}$ and $\{{\mu }_F,{\mathrm{\Omega }}_F\}$. The frequency separation between these two points can be computed as the difference

$$\Delta {\mathrm{\Omega }}_{FB}={\mathrm{\Omega }}_F-{\mathrm{\Omega }}_B=\frac{{\mathrm{\Omega }}_m\sigma }{\sqrt{{\mathrm{\Omega }}_{\gamma }^2+{(1+{\mathrm{\Omega }}_m\sigma )}^2/4}+\sqrt{{\mathrm{\Omega }}_{\gamma }^2+{(1-{\mathrm{\Omega }}_m\sigma )}^2/4}},$$

(35)

which, considering ${\mathrm{\Omega }}_m\ll 1$, can be approximated as

$$\Delta {\mathrm{\Omega }}_{FB}\approx \frac{{\mathrm{\Omega }}_m\sigma }{\sqrt{1+4{\mathrm{\Omega }}_{\gamma }^2}}\approx {\mathrm{\Omega }}_m.$$

(36)

The analytic dispersion branches and the parameters derived above are shown in Figure 3a and compared to the numerical results obtained by exactly solving (22).

Another important parameter to consider is the band-gap amplitude. Previous works have shown that this usually has only a second-order dependence on the modulation speed; differently from [5], we thus conduct this analysis setting ${\mathrm{\Omega }}_m=0$. In this case, the dispersion is symmetric, so we can also set $\mu =\pi$. Under these hypotheses, one can again solve the characteristic equation $L_{11}{L}_{22}-L_{12}^2=0$ for $\mathrm{\Omega }$, obtaining

$${\mathrm{\Omega }}_{top}=\sqrt{\frac{1}{4}+{\mathrm{\Omega }}_{\gamma }^2+{\mathrm{\Omega }}_a^2},{\mathrm{\Omega }}_{bot}=\sqrt{\frac{1}{4}+{\mathrm{\Omega }}_{\gamma }^2-{\mathrm{\Omega }}_a^2},$$

(37)

for the top and bottom frequencies of the band-gap, so that the amplitude writes

$$\Delta {\mathrm{\Omega }}_{bg}={\mathrm{\Omega }}_{top}-{\mathrm{\Omega }}_{bot}=\frac{2{\mathrm{\Omega }}_a^2}{\sqrt{\frac{1}{4}+{\mathrm{\Omega }}_{\gamma }^2+{\mathrm{\Omega }}_a^2}+\sqrt{\frac{1}{4}+{\mathrm{\Omega }}_{\gamma }^2-{\mathrm{\Omega }}_a^2}}\approx \frac{{\mathrm{\Omega }}_a^2}{\sqrt{\frac{1}{4}+{\mathrm{\Omega }}_{\gamma }^2}},$$

(38)

where the last approximation is taken, considering that ${\mathrm{\Omega }}_{\gamma }^2\gg {\mathrm{\Omega }}_a^2$. The band-gap amplitude is reported in the diagrams of Figure 3b for ${\mathrm{\Omega }}_m=0$ and in Figure 3c for ${\mathrm{\Omega }}_m=0.05$, showing in both cases accurate predictions in spite of the approximations. Finally, it is straightforward to observe that we have complete separation of the directional band-gaps (i.e., with no overlap) for a modulation speed ${\mathrm{\Omega }}_m\ge {\mathrm{\Omega }}_{m,lim}=\Delta {\mathrm{\Omega }}_{bg}$.

3.3. From Rod to Spring-Mass Chain

As mentioned at the beginning of this section, a spring-mass chain can be seen as the discrete approximation of a rod. We can use the following relationships to switch from the rod parameters to spring-mass chain ones:

$$m=\rho Aa,k=\frac{EA}{a},k_g=\gamma a,c_0=a\sqrt{\frac{k}{m}}$$

(39)

with $a={\lambda }_m/R$ being the distance between two adjacent masses and R the number of masses per spatial period (which was assumed equal to 3 in the previous section).

In Figure 4, we show the comparison of the dispersion relations for a space-time modulated rod (obtained with PWEM) and its spring-mass chain counterpart (obtained with the BBP). These results were obtained by fixing the rod properties to $EA=\rho A=1$, ${\gamma }_0=10$, ${\gamma }_m=5$, ${\lambda }_m=1$, ${\mathrm{\Omega }}_m={\mathrm{\Omega }}_{m,lim}$, and letting R vary. As it can be observed, the spring-mass chain’s diagram converges to the rod’s one, increasing the number of masses (i.e., sub-cells). In the case of the minimum number of sub-cells $R=3$, the diagrams slightly drift apart; however, band-gap positions are very close to the ones predicted with the rod model, and their widths are almost the same. Upon these observations, we can finally conclude that the analytic formulae derived in this section for a rod on elastic suspensions can be usefully exploited to guide the design of a spring-mass chain system.

4. Numerical Study of the MEMS Device

In this section, we numerically study an example of a finite-size system, with typical MEMS dimensions, composed of the unitary cells introduced in Section 2 and qualitatively depicted in Figure 1. As anticipated, we use $R=3$ masses per cell, accepting a small drift from the continuous rod case with the aim to reduce the size and complexity of the system. With reference to the scheme in Figure 5, we can define the geometrical dimension of our system by setting the values of the lumped parameters appearing in (12) and by using the following equations

$$m=\rho t(L_x{L}_y-L_{hx}{L}_{hy}),$$

(40)

$$k_g=2\frac{E{w}_b^{3}t}{L_g^3},k=\frac{2}{N_F}\frac{E{w}_b^{3}t}{L_i^3},$$

(41)

where factor 2 in the stiffness definitions is to take into account the fact that the beams come in pairs as parallel springs, and where $N_F=2$ is the number of folds. In the same way, setting the desired electrical stiffness in (10), we can define the electrode dimensions from

$$A=L_{ey}t,x_0=(L_{hx}-L_{ex})/2.$$

(42)

Notice that the stiffness modulation parameters $k_{E0}$ and $k_{E1}$ both depend on the voltages $V_0$ and $V_a$; however, since the former depends on $V_0^2+V_a^2$ and the latter on $V_0{V}_a$, if $V_0\gg V_a$, we can tune $k_{E0}$ by setting $V_0$ and $k_{E1}$ changing $V_a$. All the geometric, mechanical, and electric parameters are collected in

Table 1. Geometric parameters (in $\mathsf{\mu }$m), as shown in Figure 5.

$L_g=150$	$L_i=150$	$L_x=100$	$L_y=200$
$L_{ex}=5$	$L_{ey}=100$	$L_{hx}=7$	$L_{hy}=102$
$x_0=1$	$w_b=2$	$w_c=8$

,

Table 2. Mechanical and electrical parameters.

$E=148$ [GPa]	silicon elastic modulus
$\rho =2330$ [kg/m${}^3$]	silicon density
$t=10$ [$\mathsf{\mu }$m]	silicon wafer thickness
${\lambda }_m=312$ [$\mathsf{\mu }$m]	cell spatial period
$\epsilon =8.854\times {10}^{-12}$ [F/m]	vacuum permittivity
$V_a=2$ [V]	AC voltage amplitude
$V_0=14$ [V]	DC voltage amplitude
${\omega }_m=12,500$ [rad/s]	voltage modulation frequency

and

Table 3. Spring-mass chain lumped parameters.

$k=3.51$ [N/m]	stiffness between masses
$k_g=7.02$ [N/m]	ground stiffness
$k_{E0}=3.51$ [N/m]	constant electrostatic stiffness
$k_{E1}=0.99$ [N/m]	time-modulated electrostatic stiffness
$m=44.81$ [ng]	mass value

.

For the numerical verification, we simulate the system with an electromechanical model in COMSOL Multiphysics to take into account in full the effect of the electrostatic forces exerted by the modulated electrodes. The system is composed of $N_c=40$ cells and is forced in the middle of the chain with a tone burst. The tone burst’s central frequency is in between the two directional band-gaps. This way, both right and left-propagating waves can be analyzed with a single run. Moreover, the amplitude of the forcing is chosen in order not to trigger nonlinear effects in the electrostatic domains (the maximum mass displacement is less than $2\%$ of the gap amplitude $x_0$– with reference to (6), the ratio between the linear and cubic terms is equal to ${(x_0/x)}^2/2$, so that the former is three orders of magnitude higher than the latter).

The results are shown in Figure 6. As it can be observed, the numerical dispersion closely matches the spring-mass chain model prediction, with only a small drift due to the flexibility of the masses and of the beam connections, which are not accounted for in the analytic formulae for the springs. The dimensionless parameters introduced in Section 3 also match: the shift of the band-gaps is equal to ${\mathrm{\Omega }}_m=0.07$ (1990 Hz), the band-gap amplitude is approximately $\Delta {\mathrm{\Omega }}_{bg}=0.05$ (1300 Hz), and the cutoff frequency of the system is given by ${\mathrm{\Omega }}_{\gamma }=0.48$ (14.1 kHz).

5. Conclusions

In this work, we analyzed a one-dimensional periodic waveguide with periodically time-varying properties. The system consists of a spring-mass chain system resting on elastic foundations, whose stiffness is modulated by means of the electrostatic forces exerted by the parallel plate electrodes collocated inside the masses. We first derived a simplified model, valid for small displacement around the equilibrium, and we fully characterized its dispersion relation through the study of its continuum equivalent (rod system). Finally, our model has been validated against the numerical results of a high-fidelity model in COMSOL Multiphysics. The proposed system can be readily produced using standard MEMS technologies and, as demonstrated by the dimensionless analysis of the dispersion, can be easily tailored to work in different frequency ranges, going (considering typical MEMS dimensions) from a few kHz onward, just by properly selecting the spring-mass chain parameters. For this reason, the presented layout can be both interesting to be used for new low-frequency applications and/or to replace already existing components in high-frequency ones.