# Analysis and Modeling of an 11.8 GHz Fin Resonant Body Transistor in a 14nm
FinFET CMOS Process

Udit Rawat, Bichoy Bahr, Dana Weinstein This work was supported in part by the
DARPA MIDAS Program. Udit Rawat is with the Department of Electrical
Engineering, Purdue University, West Lafayette, IN 47906 USA (e-mail:
rawatu@purdue.edu).Bichoy Bahr is with Kilby Labs - Texas Instruments, Dallas,
TX, USA.Dana Weinstein is with the Department of Electrical Engineering,
Purdue University, West Lafayette, IN 47906 USA (e-mail: danaw@purdue.edu).

###### Abstract

In this work, a compact model is presented for a 14 nm CMOS-based FinFET
Resonant Body Transistor (fRBT) operating at a frequency of 11.8 GHz and
targeting RF frequency generation/filtering for next generation radio
communication, clocking, and sensing applications. Analysis of the phononic
dispersion characteristics of the device, which informs the model development,
shows the presence of polarization exchange due to the periodic nature of the
back-end-of-line (BEOL) metal PnC. An eigenfrequency-based extraction process,
applicable to resonators based on electrostatic force transduction, has been
used to model the resonance cavity. Augmented forms of the BSIM-CMG (Common
Multi-Gate) model for FinFETs are used to model the drive and sense
transistors in the fRBT. This model framework allows easy integration with the
foundry-supplied process design kits (PDKs) and circuit simulators while being
flexible towards change in transduction mechanisms and device architecture.
Ultimately, the behaviour is validated against RF measured data for the
fabricated fRBT device under different operating conditions, leading to the
demonstration of the first complete model for this class of resonant device
integrated seamlessly in the CMOS stack.

###### Index Terms:

CMOS-Microelectromechanical Systems (MEMS), FinFET, compact model, Fin
Resonant Body Transistor (fRBT), phononic crystals, Radio frequency (RF) MEMS,
resonators, circuit simulation.

##  I Introduction

Increasing proliferation of cm- and mm-Wave 5G mobile communication technology
to address the demand for high data rates, better reliability and low-latency
necessitates innovation in the field of front-end electroacoustic devices for
filtering and carrier generation. AlN and LiNbO3 bulk acoustic wave (BAW)
resonators and filters utilizing higher modes of operation have recently shown
promise [1][2] with regards to applicability in the currently allocated
n257-n260 (24.5 to 40 GHz) 5G mm-Wave bands. However, these resonators have
specialized fabrication and packaging requirements making their integration
with CMOS prohibitive even with MEMS-last [3], Front-end-of-line (FEOL) [4]
and Back-end-of-line (BEOL) post-processed MEMS [5][6][7] approaches. The
typical solution of packaging the MEMS die separately from CMOS exacerbates
the problem because of interconnect bandwidth limitations across multiple
chips at cm- and mm-Wave frequencies. Moreover, at these frequencies, carrier
generation and distribution to all the channels in a multi-element phased
array IC via the conventional off-chip crystal and PLL combination, as shown
in [8], results in significant routing power dissipation. Availability of
high-Q, integrated CMOS-MEMS resonator-based oscillators with good phase noise
levels at the targeted frequencies would result in a reduction in the carrier
power since the central PLL would no longer be necessary. Considering the
aforementioned challenges and potential opportunities at the aforementioned
high frequencies, monolithic integration of RF/mmWave MEMS resonators in a
conventional CMOS process becomes an attractive proposition.

Fully-integrated, solid-state, CMOS-MEMS RF/mmWave resonators have previously
been demonstrated using different technology nodes [9][10]. These resonators
make use of acoustic waveguiding confinement based on Back-end-of-Line (BEOL)
metal phononic crystals (PnCs) and adiabatic terminations for mode
localization to form a resonance cavity concentrated at the transistor layers
of the CMOS stack. The targeted mode is excited differentially using metal-
oxide-semiconductor (MOS) capacitors, or MOSCAPs, and sensed using a pair of
transistors incorporated into the resonant cavity, biased in saturation. The
drain current in these sense FETs is modulated by the stress generated in the
channel due to vibration, resulting in a differential readout. Since these
resonators are fully integrated within a given CMOS technology, high
performance oscillators can be designed (e.g. [11]) with significant savings
in terms of carrier power generation and distribution as well as area and
cost. To design and correctly estimate the performance of such oscillators,
precise, physics-accurate compact models are required for the constituent
CMOS-MEMS resonant devices. These models are expected to capture the coupled
physics of the various transduction mechanisms as well as the mechanical
behaviour of the device without having to resort to computationally intensive
Finite Element Method (FEM) simulations.

An initial circuit model for a transistor-sensed CMOS-MEMS resonator, as
depicted in [12] has previously been devised to capture the small-signal
behaviour of a Resonant Body Transistor based on a simplified analytical model
of the drive capacitor and sense transistor. In that model, the passive
section of the device consists of a drive MOSCAP in accumulation and a cavity
whose resonant behaviour is captured using a series resistance (R), inductance
(L) and capacitance (C) equivalent circuit. The transistor sensing is
implemented using an electromechanical transconductance $g_{m,em}$ which
represents the transduction from the cavity mechanical resonance to drain
current modulation. This model explains the fundamental behaviour of the
resonator, but is not sufficiently detailed and suffers from two significant
drawbacks. First, owing to the low motional capacitance of these devices with
respect to the static drive capacitance $c_{0}$, it is very challenging to
extract the equivalent R, L, and C values for the passive section
experimentally at high frequencies. Second, oscillators typically operate in
the regime of large signal swings to achieve low phase noise, so a small
signal equivalent circuit does not capture the effect of various
nonlinearities inherent to the device.

Figure 1: (a) Complete 3D FEM Simulation model depicting two adjoining fRBT
unit cells. Mx (x=1-3) and Cy (y=4-6) represent the first 6 metal levels that
form a part of the BEOL PnC. (b) Inset showing the front-end-of-line (FEOL)
channel region with the gate stack hidden. The quantities $a$, $b$ and $L$
represent lattice constants in the $x_{1}$ and $x_{2}$ directions and the gate
length respectively.

Compact models were subsequently developed for a 1-D Unreleased [13] and
released [14] Resonant Body Transistor (RBT) which use modified MIT Virtual
Source and BSIM planar FET models for the sense transistors only. While these
model are large-signal in nature and improve upon the drawbacks of the small-
signal model in [12], they do not translate well to CMOS-integrated RBTs built
using FinFETs as in [10] and the resonator considered in this paper. The
primary reason for this being the 3D nature of the silicon fin together with
the complex mode shape requiring additional analysis for the modeling of
various effects of the stress developed in the sense transistor channel. The
model in [13] considers drain current modulation in the sense transistor only
due to mobility modulation. However, in the actual device, stress and strain
in the transistor channel cause changes to multiple other device parameters
which need to be considered. These models also do not consider electrostatic
drive using MOS capacitors. For the complete device model to be compatible
with Electronic Design Automation (EDA) tools and the foundry-supplied process
design kits (PDKs) used for simulation and design, the industry standard BSIM-
CMG model [15] is required to be augmented to capture these effects
accurately.

Figure 2: Cross-sectional SEM of the 11.8 GHz fRBT resonator depicting the
resonant cavity bound laterally by the termination regions. The BEOL metal
phononic crystal above the gate region can also be seen.

Waveguiding-based CMOS-MEMS resonators frequently exhibit complicated mode
shapes with interesting phononic dispersion behaviour owing to the intricate
geometry and range of materials in advanced-node FinFET CMOS technology. A
standard abstraction technique is required to reduce the mode shape to an
equivalent mechanical description that can be integrated into the overall
compact model for the resonator. Theoretical techniques for extracting the
equivalent mass, stiffness, and damping of a resonator are not applicable in
this scenario. Thus, a technique such as that described in [16] can used for
CMOS-MEMS RBTs.

In this paper, a large-signal compact model for a 11.8 GHz Fin Resonant Body
Transistor (fRBT) as shown in Fig. 2, fabricated using a commercial 14 nm
FinFET (GlobalFoundries® 14LPP) process has been presented. The organisation
of this paper is as follows: Section II gives a detailed description of the
unit-cell based 3-D FEM model framework as well as an analysis of the phonon
dispersion. In Section III, the equivalent mechanical parameter extraction
procedure for the waveguide cavity is described. Section IV addresses the
implementation details of the individual constituent modules in the complete
fRBT model as well as their interconnections. In Section V, the model is
benchmarked against measured resonator data to demonstrate its utility in
real-world simulation scenarios. Finally, Section VI presents a conclusion to
the study.

##  II Mechanical Resonance and Dispersion Analysis

A robust, FEM-based, mechanical simulation and analysis framework forms the
basis for the compact model of the fRBT. To analyze the mechanical resonance
characteristics of the designed fRBT device, a 3D simulation model for a unit
cell, as shown in Fig. 1, is constructed in COMSOL Multiphysics®. Simulations
are used to extract the phononic band structure of the acoustic waveguide,
which is formed by arraying the fRBT unit cells. The mode of interest and its
corresponding dispersion are then analyzed to understand the nature of the
displacement and stress fields in the FEOL region. This is relevant for
modeling the effect on the drain current as well as equivalent parameter
extraction.

Figure 3: (a) Cavity region layout including layers up to metal layer M1,
highlighting the unit cell along with the X and Y direction PBCs. Higher metal
layers and vias are not shown. PC and CA layers in the layout represent the
gate and source/drain contacts respectively. (b) Irreducible Brillouin Zone
(IBZ) region for the unit cell highlighting symmetry points. Point X is most
amenable to excitation with the designed unit cell.

###  II-A Unit Cell 3D FEM Model

Simulation of an entire 3D model of the resonator is not computationally
feasible. Therefore, a pair of repeating unit cells are simulated using
periodic boundary conditions (PBCs) along the $x_{1}$ and $x_{2}$ directions.
While the use of PBCs to model the resonator cavity introduces a certain
degree of inaccuracy because of the finite number of gate interdigitated
transducer (IDT) fingers, this is later accounted for in the compact model
through an adjustment factor. The FEOL section geometry of the unit cell is
constructed using process parameters such as fin height $HFIN$, fin thickness
$TFIN$ etc., as well as an understanding of the GF14LPP FinFET process flow. A
single gate IDT represents a transistor (number of fingers $NF=1$) with a fin
array (number of fins $NFIN$) along the $x_{2}$ direction separated by the
technology-defined parameter $FPITCH$ which represents the fin pitch as shown
in the layout in Fig. 3(a). The only design parameter available for
modification in this section of the geometry is the gate length $L$, which
determines the resonance frequency in the particular fRBT mode under
consideration in this work. The BEOL portion of the unit cell consists of Mx
and Cx level Copper metal phononic crystal (PnC) which is designed such that
the phononic bandgap (PnBG) encompasses the resonance frequency to provide
confinement along the positive $x_{3}$ direction (out of the plane of the
chip). The PnC design takes place within the limits set by the design rules
for the process. Each element of the PnC in the BEOL extends uniformly along
the gate finger direction $x_{2}$. Appropriate material assignments are done
for all regions of the structure to complete the unit cell design.

Eigenfrequency analysis is required to obtain the mode shapes and
corresponding resonance frequencies. Since we only have a single gate IDT per
unit cell, the electrical excitation couples most efficiently to the modes
corresponding to $k_{x}=\pi/a$ and $k_{y}=0$ ie. point $X$ at the edge of the
Irreducible Brillouin Zone (IBZ) as shown in Fig. 3(b). The unit cell is set
up to be excited by the wave vector $\overrightarrow{k}=k_{x}\hat{x}_{1}$ to
obtain the eigenstates and the eigenfrequencies where $\hat{x}_{1}$ represents
direction in reciprocal space.

###  II-B Theoretical Formulation for Modal and Dispersion Analysis

The Plane Wave Expansion framework as described in [17] can be applied to the
phononic waveguide comprising of the fRBT unit cells to obtain a qualitative
understanding of the nature of the mode shapes and the dispersion
characteristics of the device. The mass density $\rho$ and the elastic moduli
$c_{ijkl}$ in the waveguide vary with and are periodic functions of the
position vector $\overrightarrow{r}$. If $u_{i}(\overrightarrow{r})$ ($i=1-3$)
denotes the displacement field components along the $x_{i}$ directions and
$T_{ij}(\overrightarrow{r})$ is the stress, then the Hooke’s law can be
written in the form:

$T_{ij}(\overrightarrow{r})=\sum_{kl}c_{ijkl}(\overrightarrow{r})\frac{\partial
u_{k}(\overrightarrow{r})}{\partial x_{l}}$ (1)

where, $i$, $j$, $k$ and $l$ can be 1, 2 or 3. The equation of motion in
accordance with Newton’s second law can we written in the form:

$\rho(\overrightarrow{r})\frac{\partial^{2}u_{i}(\overrightarrow{r})}{\partial
t^{2}}=\sum_{j}\frac{\partial T_{ij}(\overrightarrow{r})}{\partial x_{j}}\\\
=\sum_{j}\frac{\partial}{\partial
x_{j}}\left[\sum_{kl}c_{ijkl}(\overrightarrow{r})\frac{\partial
u_{k}(\overrightarrow{r})}{\partial x_{l}}\right]$ (2)

The materials in the CMOS stack are either isotropic or cubic symmetric in
nature which results in the elimination of some of the elements of the
$c_{ijkl}$ tensor. Substituting equation (1) into (2) and converting to Voigt
notation we obtain three coupled equations of motion of the form (position
vector dependence of $\rho$ and $c$ has not been shown):

$\begin{aligned} \rho\frac{\partial^{2}u_{i}}{\partial
t^{2}}=\frac{\partial}{\partial x_{i}}\left(c_{11}\frac{\partial
u_{i}}{\partial x_{i}}+c_{12}\left(\frac{\partial u_{j}}{\partial
x_{j}}+\frac{\partial u_{k}}{\partial x_{k}}\right)\right)\\\
+\frac{\partial}{\partial x_{j}}\left(c_{44}\left(\frac{\partial
u_{i}}{\partial x_{j}}+\frac{\partial u_{j}}{\partial
x_{i}}\right)\right)+\frac{\partial}{\partial
x_{k}}\left(c_{44}\left(\frac{\partial u_{i}}{\partial x_{k}}+\frac{\partial
u_{k}}{\partial x_{i}}\right)\right)\end{aligned}$

(3)

where $i$,$j$ and $k$ are 1,2 and 3 respectively for the equation of motion
corresponding to the displacement $u_{1}$. Similar equations can be written
for the other two components, $u_{2}$ and $u_{3}$. The wave solutions to these
equations of motion are of the form
$\overrightarrow{u}(\overrightarrow{r})e^{-i\omega t}$ where $\omega$ is the
angular frequency. Waves inside a periodic structure such as the fRBT are
analogous to plane waves but are modulated by an envelope function. The
envelope function takes on the same symmetry and periodicity as the underlying
structure. According to the Bloch theorem:

$\overrightarrow{u}(\overrightarrow{r})=e^{i\overrightarrow{k}.\overrightarrow{r}}\overrightarrow{U}_{\overrightarrow{k}}(\overrightarrow{r})$
(4)

where $\overrightarrow{k}(k_{1},k_{2},k_{3})$ is the Bloch wave vector. The
term $e^{i\overrightarrow{k}.\overrightarrow{r}}$ in equation 4 represents a
plane wave like phase ”tilt” term and
$\overrightarrow{U}_{\overrightarrow{k}}(\overrightarrow{r})$ is the envelope
function. The envelope function has the same periodicity as the fRBT structure
such that
$\overrightarrow{U}_{\overrightarrow{k}}(\overrightarrow{r})=\overrightarrow{U}_{\overrightarrow{k}}(\overrightarrow{r}+\overrightarrow{R})$.
The quantity $\overrightarrow{R}$ belongs to the real space lattice. The
envelope function can be expressed in terms of a Fourier series as:

$\overrightarrow{U}_{\overrightarrow{k}}(\overrightarrow{r})=\sum_{\overrightarrow{G}^{\prime}}\overrightarrow{U}_{\overrightarrow{k}}(\overrightarrow{G}^{\prime})e^{i\overrightarrow{G}^{\prime}.\overrightarrow{r}}$
(5)

which leads to

$\overrightarrow{u}(\overrightarrow{r},t)=e^{-i\omega
t}\sum_{\overrightarrow{G}^{\prime}}\overrightarrow{U}_{\overrightarrow{k}}(\overrightarrow{G}^{\prime})e^{i(\overrightarrow{G}^{\prime}+\overrightarrow{k}).\overrightarrow{r}}$
(6)

Similarly to the displacement field, $\overrightarrow{u}(\overrightarrow{r})$,
$\rho(\overrightarrow{r})$ and $C(\overrightarrow{r})$ can be written in terms
of a Fourier expansion using
$\alpha(\overrightarrow{r})=(\rho(\overrightarrow{r}),C(\overrightarrow{r}))$:

$\alpha(\overrightarrow{r})=\sum_{\overrightarrow{G}^{\prime\prime}}\alpha(\overrightarrow{G}^{\prime\prime})e^{i\overrightarrow{G}^{\prime\prime}.\overrightarrow{r}}$
(7)

where $\overrightarrow{G}^{\prime}$ and $\overrightarrow{G}^{\prime\prime}$
belong to the reciprocal lattice domain. Substitution of (6) and (7) into the
equations of the type (3), gives the Fourier Transforms of the equations of
motion. After some manipulation and modal projection of the resulting
equations we get the eigenvalue problem:

$\omega^{2}\begin{pmatrix}Q^{(11)}_{\overrightarrow{G},\overrightarrow{G}^{\prime}}&0&0\\\
0&Q^{(22)}_{\overrightarrow{G},\overrightarrow{G}^{\prime}}&0\\\
0&0&Q^{(33)}_{\overrightarrow{G},\overrightarrow{G}^{\prime}}\end{pmatrix}\begin{pmatrix}U_{1,\overrightarrow{k}}(\overrightarrow{G}^{\prime})\\\
U_{2,\overrightarrow{k}}(\overrightarrow{G}^{\prime})\\\
U_{3,\overrightarrow{k}}(\overrightarrow{G}^{\prime})\end{pmatrix}=\\\
\begin{pmatrix}P^{(11)}_{\overrightarrow{G},\overrightarrow{G}^{\prime}}&P^{(12)}_{\overrightarrow{G},\overrightarrow{G}^{\prime}}&P^{(13)}_{\overrightarrow{G},\overrightarrow{G}^{\prime}}\\\
P^{(21)}_{\overrightarrow{G},\overrightarrow{G}^{\prime}}&P^{(22)}_{\overrightarrow{G},\overrightarrow{G}^{\prime}}&P^{(23)}_{\overrightarrow{G},\overrightarrow{G}^{\prime}}\\\
P^{(31)}_{\overrightarrow{G},\overrightarrow{G}^{\prime}}&P^{(32)}_{\overrightarrow{G},\overrightarrow{G}^{\prime}}&P^{(33)}_{\overrightarrow{G},\overrightarrow{G}^{\prime}}\end{pmatrix}\begin{pmatrix}U_{1,\overrightarrow{k}}(\overrightarrow{G}^{\prime})\\\
U_{2,\overrightarrow{k}}(\overrightarrow{G}^{\prime})\\\
U_{3,\overrightarrow{k}}(\overrightarrow{G}^{\prime})\end{pmatrix}$ (8)

which can also be written in the form:

$\omega^{2}\overleftrightarrow{Q}\overleftrightarrow{U}=\overleftrightarrow{P}\overleftrightarrow{U}$
(9)

solution to equation (9) gives the eigenstates (modes) and their corresponding
eigenfrequencies. The elements of the matrices $\overleftrightarrow{P}$ and
$\overleftrightarrow{Q}$ determine the polarization and coupling between the
modes, respectively. Considering that the fRBT is periodic in the $x_{1}$ and
$x_{2}$ directions and the BEOL PnC is also periodic for the metal layers Mx
in the positive $x_{3}$ direction, $G_{1}$, $G_{2}$ and $G_{3}$ are assumed to
be non-zero. Since we are looking at exciting modes at the symmetry point $X$
using the gate IDT along $x_{2}$, we get $k_{2},k_{3}=0$. The values of all
the elements in the matrices $\overleftrightarrow{P}$ and
$\overleftrightarrow{Q}$ can be calculated and it is seen that under the
assumptions for the components of $\overrightarrow{G}$ and
$\overrightarrow{k}$, none of the terms in the matrix $\overleftrightarrow{P}$
reduce to zero. This means that all of the displacement components can couple
to each other in different ways.

###  II-C 3-D FEM Modal and Dispersion Analysis

---

---

Figure 4: (a) Phononic dispersion of the fRBT unit cell obtained using FEM
simulation, depicting the targeted mode of operation. (b) Displacement
profiles for the PnC modes at $k=1$ (points B, D and E).

The results from the eigenfrequency simulation in COMSOL Multiphysics® are
shown in Fig. 4. As mentioned previously, since the most efficient drive using
gate IDTs corresponds to the symmetry point $X$ at the edge of the IBZ, only a
part of the $\Gamma-X$ path (highlighted in red in Fig. 3(b)) is traced for
the dispersion analysis. The fRBT unit cell is driven with wave vectors with
$k_{x}$ ranging from $0.6\pi/a$ to $\pi/a$ and $k_{y}=0$ and the corresponding
eigenfrequencies and eigenstates (mode shapes) are evaluated. The resulting
dispersion plot as seen in Fig. 4(a), shows the modes lying below the sound
cone supported by the structure. The shear wave sound-line in Silicon
(represented by the red line $\omega=c_{shear}k_{x}$) demarcates the FEOL-BEOL
waveguide mode region from the region where these waveguide modes can couple
to bulk modes in the Si substrate. As a general design principle, to reduce
scattering to bulk modes, the difference in $\overrightarrow{k}$-space between
the mode of interest and the sound-line should be maximized [18]. Moreover,
for the modes to be sensed and driven efficiently to maximize the
electromechanical transconductance $g_{m,em}$, the stress localization should
be the strongest in the Si fin of the transistors. Considering these criteria,
the mode highlighted in yellow in Fig. 4(a) is the targeted mode for the fRBT.
Other modes that are present in Fig. 4(a) are either localized entirely within
the BEOL Mx/Cx PnC or Rayleigh modes at the top and bottom extremities of the
structure arising due to the finite nature of the simulation model. The modes
localised within the BEOL PnC which cannot be driven or sensed efficiently are
depicted in Fig. 4(b).

##  III Equivalent Mechanical Parameter Extraction

Figure 5: Small signal equivalent circuit model of the passive component of
the fRBT. The drive MOSCAP capacitance corresponding to each of the two phases
in the differential scheme is represented by $c_{0}$, the capacitive
electromechanical transduction coefficient by $\eta_{eff}$, and the equivalent
mass-spring-damper by $m_{eff}$, $k_{eff}$ and $b_{eff}$, respectively. The
total differential current into the differential 1-port $i_{tot}$ is the sum
of the feed-through current $i_{ft}$ through $c_{0}$ and the motional current
$i_{mot}$.

Once the target mode for the fRBT has been identified and analyzed, we can
then develop its equivalent representation. The 1-D mechanical equivalent for
the targeted mode forms the core of the complete fRBT compact model. RF
measurements are typically used to extract the equivalent parameters for MEMS
resonators. However, this technique is not applicable in the case of fRBT
devices operating at frequencies greater than 10 GHz, the reason for which can
be explained using the small signal equivalent circuit of the passive section
of the device (without the sense transistors) as shown in Fig. 5. At the
targeted operation frequencies, the feed-through current $i_{ft}$ through the
static MOSCAP capacitance $C_{0}$ is significantly higher than the motional
current through the resonator body. This makes the detection of differential
1-port impedance/admittance dependence on frequency exceptionally challenging
even with extraction techniques such as that described in [19], since the
resonance is completely embedded in feed-through.

An alternative approach for extraction has been described in [16] which makes
use of eigenfrequency simulations. The main principle behind this approach is
the equivalence of the work done in a single degree-of-freedom (DOF)
equivalent circuit to that done in a system with multiple DOF such as the
fRBT. The transduction scheme in the fRBT is internal dielectric transduction
[20][21] via the gate dielectric of the drive MOSCAPs which is different from
the externally transduced resonator discussed in [16]. Because of a difference
in the capacitive transduction mechanism the extraction procedure is required
to be amended. The electrostatic work done by the MOSCAP actuator is given by:

$W_{e}=|F_{e}|u_{i,eq}=\frac{1}{A_{act}}\iint_{A_{act}}\overrightarrow{F_{e}}.\overrightarrow{u_{i}}dA$
(10)

Figure 6: (a) Fin region of the fRBT showing the actuation surfaces for the
electrostatic force. (b) Displacement plots for the actuation surfaces showing
the total displacement along with all the components along the three principal
directions.

where $\overrightarrow{F_{e}}$ is the applied electrostatic force, $u_{i,eq}$
the equivalent displacement and $\overrightarrow{u_{i}}$ the vector
displacement of the actuation surface having an area $A_{act}$. The
electrostatic force is applied between the charge on the metal gate and the
charge carriers at the channel surface of the fin. Fig. 6(a) presents the Si
fin of the transistor clad with the gate dielectric stack where surface 1
denotes the Si-SiO2 interface and surface 2 corresponds to the gate-metal-
dielectric interface. Assuming that $A_{1}$ and $A_{2}$ are the areas of
surfaces 1 and 2, respectively, the equivalent displacements of the two
surfaces after eliminating the force term $|\overrightarrow{F_{e}}|$ in
equation 10 are then given by:

$u_{1}=\frac{\iint_{A_{1}}\overrightarrow{u_{i}}.\hat{n}dA_{1}}{\iint_{A_{1}}dA_{1}}$
(11)

$u_{2}=\frac{\iint_{A_{2}}\overrightarrow{u_{i}}.\hat{n}dA_{2}}{\iint_{A_{2}}dA_{2}}$
(12)

where $\hat{n}$ is the unit vector normal to the actuation surfaces. The gate
dielectric is thin and the variation in thickness is negligible as
corroborated by the FEM simulations. Thus the equivalent displacement
$u_{i,eq}$ can be expressed as an average of the equivalent displacements of
the two actuation surfaces:

$u_{i,eq}=\left|\frac{u_{1}+u_{2}}{2}\right|$ (13)

The multi-DOF stored energy in the fRBT unit cell evaluated using a volume
integral in COMSOL, $E_{stored}$, is equated to the stored energy in the
equivalent single degree of freedom system to give the equivalent parameters:

$k_{eq}=\frac{2E_{stored}}{|u_{i,eq}^{2}|}$ (14)

$m_{eq}=\frac{k_{eq}}{|\omega_{0}^{2}|}$ (15)

$b_{eq}=\frac{\sqrt{m_{eq}k_{eq}}}{Q}$ (16)

where $\omega_{0}$ is the resonance frequency and $Q$ is the Q-factor as
obtained from the electromechanical transconductance RF measurement. From Fig.
6(b) it can be seen that $|u_{x3}|$ is the dominant component in the total
displacement $|u|$. The extracted values of the mechanical parameters for the
cavity obtained using FEM eigenfrequency simulation for the targeted mode have
been listed in Table I. The same procedure can be used to characterize
spurious modes that may be excited besides the main resonance.

TABLE I: Extracted Mechanical Parameters

PARAMETER | VALUE
---|---
$m_{eq}$ |  $1.0157$ x $10^{-16}$ kg
$k_{eq}$ |  $0.394$ x $10^{6}$ N/m
$b_{eq}$ |  $1.054$ x $10^{-8}$ kg/s

##  IV Compact Model Implementation

In this section, based upon the understanding of the fRBT device structure and
the targeted mode shape, the implementation details of the constituent
building blocks for the complete compact model are discussed. The measured
characteristics of the fabricated fRBT device are used for optimizing as well
as benchmarking the developed model. The model implementation using VerilogA
takes into account compatibility with standard circuit simulators and the PDK
for the GF14LPP process.

###  IV-A Drive MOSCAP Module

---

---

Figure 7: (a) Test-bench circuit schematic for the fbsimd transistor. The
circuit branch highlighted in blue represents the connection point to the
equivalent circuit for the structure’s mechanical quantities. (b) Transient
electrostatic force output of the fbsimd transistor for three drive voltage
levels of 10 mV, 50 mV and 100mV at two gate DC biases of 600 mV and 800 mV.

Standard-Vt transistors are used to model the array of MOSCAPs that serve as
gate-drive transducers. The BSIM-CMG model for FinFET transistors does not
account for the electrostatic force generated across the “plates” of the
MOSCAP. Like every other physical capacitor, a force exists between the gate
metal and the channel surface charge of the MOSCAP transistor whose source and
drain are shorted together. Application of a sinusoidal voltage at the gate
results in the generation of an AC force across the gate dielectric which
couples to the target mode shape resulting in vibration of the fRBT cavity.
The BSIM-CMG model is modified (with the modified model termed “$fbsimd$”) to
include the generation of this electrostatic actuation force and a
corresponding additional mechanical terminal is added alongside the original
gate, drain, source, and bulk electrical terminals.

From [22] it can be seen that the standard BSIM-CMG model calculates intrinsic
capacitances of a transistor using the terminal charges at all four terminals.
Similarly, in $fbsimd$, the calculated gate charge $Q_{G}$ and intrinsic gate
capacitance $c_{ox}$ are used to evaluate the driving voltage $v_{D}$ across
the gate dielectric film:

$v_{D}(t)=\frac{Q_{G}(t)}{c_{ox}(t)}$ (17)

The drive force generated by a drive MOSCAP IDT is then given by:

$F_{d}(t)=-\chi\frac{\epsilon_{0}\epsilon_{r}}{2}\frac{Av_{D}^{2}(t)}{t_{ox}^{2}}$
(18)

where A is the effective actuation area of the drive transducer that depends
on the number of fins $NFIN$, other fin geometry parameters as well as the
gate length $L$, $t_{ox}$ is a model parameter representing the gate
dielectric thickness and $\epsilon_{r}$ is the relative gate dielectric stack
permittivity. The model parameter $\chi$ is a drive force adjustment parameter
that is introduced to account for the effect of a finite transducer length
which modifies effective coupling to the cavity mode shape.

Transient simulation for the $fbsimd$ transistor is carried out using the
test-bench as shown in Fig. 7(a) and the resulting force waveforms are
analysed. As can be seen in each of the DC biasing cases of Fig. 7(b), for
small drive voltage levels such as 10 mV, the distortion of the output force
waveform $F_{d}(t)$ is minimal. However, for large drive voltage levels such
as 100 mV, the distortion in $F_{d}(t)$ is appreciable. The nonlinearity in
the $F_{d}(t)$ waveform can be attributed to the bias dependent nature of the
MOSCAP capacitance which fluctuates when the $fbsimd$ transistor gate is
driven by a large signal. The distortion increases as the DC biasing reaches
closer to the $V_{t}$ of the $fbsimd$ (which is same as the standard-$V_{t}$
n-channel FinFET). The drive MOSCAPs should therefore be operated in the
strong inversion or accumulation regimes.

###  IV-B Resonant Body Module

Using the equivalent mechanical parameters for the fRBT waveguide cavity
extracted in Section III, a mechanical resonant body module has been
constructed. Spurious mode resonances are added alongside the target
resonance, all of which can be modeled using the force-current formulation
[23]. The state-space representation for the $i^{th}$ mode in the module is as
follows:

$\nu_{i}(t)=\frac{dx_{i}(t)}{dt}$ (19a)
$\vspace{2mm}F_{d,i}(t)=m_{eq,i}\frac{d\nu_{i}(t)}{dt}+b_{eq,i}\nu_{i}(t)+k_{eq,i}x_{i}(t)+\sqrt{4k_{B}Tb_{eq,i}}$
(19b) $\vspace{2mm}F_{d,i}(t)=\gamma_{i}F_{d}(t)$ (19c)

where $v_{i}(t)$ is the velocity associated with the mechanical node
$x_{i}(t)$, $k_{B}$ is the Boltzmann Constant and $T$ is the temperature. The
$\sqrt{4k_{B}Tb_{eq,i}}$ term added to the model to account for the mechanical
force noise in the device. The term $\gamma_{i}$ is used to model the
efficiency with which the generated electrostatic force $F_{d}(t)$ couples to
the $i^{th}$ mode. The force-current analogy is best suited for implementation
of the resonant body since it allows the addition of multiple drive modules,
each of which contributes a drive force current. To combine the effects of all
the modes considered in the model, the velocities corresponding to each mode
are summed:

$\nu(t)=\sum_{i=1}^{N}\nu_{i}(t)$ (20)

Typically, series $RLC$ branches are added in parallel corresponding to each
mode as shown in [24]. The complement of the this circuit convention is
considered in the case of the resonant body module of the fRBT since the
force-current analogy is employed.

###  IV-C Sense Transistor Module

Although the drive MOSCAP transducers together with the resonant cavity are
sufficient to create a resonator, owing to the difficulties with RF detection
of the resonator response in an all-passive implementation as highlighted in
section III, transistor-based readout is the optimal choice for the fRBT
device. Two standard-$V_{t}$ transistors which act as vibration sensors are
embedded in the centre of the waveguide cavity. These transistors are exactly
the same as the drive MOSCAP transistors except that they are connected and
biased so as to conduct current through their channels. The time-varying
stress in the waveguide cavity causes a modulation in the drain current of the
sense transistors which results in a differential current readout. Stress
induced in the channel of the sense FinFET transistors causes changes to the
electronic band structure of the channel material (in this case Silicon). This
effect primarily manifests itself in the form of variation in carrier mobility
$\mu_{n}$ through the piezoresistive effect, threshold voltage $V_{t}$, and
saturation velocity $\nu_{sat}$.

The effect of stress on transistor characteristics has been studied
extensively in literature in the case of Layout Dependent Effects (LDEs) due
to shallow trench isolation (STI) etc. [25][26][27] and strained-silicon
transistors for mobility enhancement [28]. The parameters used for modeling
the dependence of drain current on stress in the conventional BSIM-CMG model
become time-variant in the case of the fRBT. A new modified BSIM-CMG model
called “$fbsims$” is therefore developed with an additional mechanical
terminal (similar to the $fbsimd$ model) and which can calculate the time-
varying changes to $\mu_{n}$, $V_{t}$ and $\nu_{sat}$.

---

---

Figure 8: (a) Stress profiles of the target mode shape for the $x_{2}-x_{3}$
plane in the fin region. (b) Stress component profiles for a pair of adjacent
unit cell gate transducers in the $x_{1}-x_{3}$ plane, depicting energy
localization in the FEOL region. Opposite signs for the stress are seen in the
adjacent transducers because of the differential nature of the mode shape. A
drive stress of 6 MPa has been considered for the frequency domain simulation.

It is important to analyze the nature of time varying stresses and strains in
the sense transistor fin region before their effect on transistor properties
can be modeled. In Fig. 8(a) we can see the stress profiles in the fin region
along the $x_{2}-x_{3}$ plane corresponding to the target mode shape.
Considering the current density for the ”on” state in a tapered fin as
simulated in [29], it can be noted that the distribution is uniform across the
majority of the fin cross-section in the $x_{2}-x_{3}$ plane. Moreover, as
seen in Fig. 8(b), all three stress components are also uniform in the fin
region along the $x_{1}$ direction. Therefore, we can use a volume averaged
stress formulation in the fin region to calculate stress values for transistor
parameter variations:

$\overline{T}=\frac{\iiint_{V}Tdx_{1}dx_{2}dx_{3}}{\iiint_{V}dx_{1}dx_{2}dx_{3}}$
(21)

####  IV-C1 Mobility Modulation Model

TABLE II: FINFET PIEZORESISTIVITY COEFFICIENTS [30]

|  $\pi_{11}^{\prime}$ (Pa-1) |  $\pi_{12}^{\prime}$ (Pa-1) |  $\pi_{12}$ (Pa-1)
---|---|---|---
NFET |  $-45.2$x$10^{-11}$ |  $-25.6$x$10^{-11}$ |  $57.6$x$10^{-11}$
PFET |  $45.0$x$10^{-11}$ |  $-23.8$x$10^{-11}$ |  $-10.1$x$10^{-11}$

Piezoresistivity theory for Si can be used to model the effect of channel
stress on the electron mobility $\mu_{n}$. It is assumed that the FinFET
transistor channel in the GF14LPP process is oriented along the $\langle
110\rangle$ direction so as to boost the weaker hole mobility $\mu_{p}$. Thus,
the piezoresistivity coefficients of Silicon as described in [31] and [32],
which are based on the $\langle 100\rangle$ crystallographic axis, must be
modified to account for the coordinate system rotation by 45o. Also, the
piezoresistive coefficients in Si FinFET devices differ from those in bulk Si
[30] and hence, the values (listed in Table II) for the $\langle 110\rangle$
direction are used in this work. The relative change in mobility with stress
is given by the relation [33]:

$\frac{\Delta\mu}{\mu}=\pi_{11}^{\prime}T_{x_{1}x_{1}}+\pi_{12}^{\prime}T_{x_{2}x_{2}}+\pi_{12}T_{x_{3}x_{3}}$
(22)

where

$\pi_{11}^{\prime}=\frac{\pi_{11}+\pi_{12}+\pi_{44}}{2}$ (23a)
$\pi_{12}^{\prime}=\frac{\pi_{11}+\pi_{12}-\pi_{44}}{2}$ (23b)

In equation (22), $T_{x_{1}x_{1}}$, $T_{x_{2}x_{2}}$ and $T_{x_{3}x_{3}}$
represent the time-varying, volume-averaged stresses obtained using equation
(21) from COMSOL. While $\pi_{11}$, $\pi_{12}$ and $\pi_{44}$ represent the
piezoresistivity coefficients along the $\langle 100\rangle$ crystallographic
axes, $\pi_{11}^{\prime}$ and $\pi_{12}^{\prime}$ represent piezoresistivity
coefficients for the rotated system along the $\langle 110\rangle$ directions
corresponding to the $x_{1}$, $x_{2}$ and $x_{3}$ directions.

A mobility multiplier $1+\frac{\Delta\mu}{\mu}$ is used in conjunction with
the mobility degradation factor $D_{mob}$ [22] in the drain current equation
of the $fbsims$ model to account for the variation of the transistor drain
current with the stress in the sense transistor fin.

####  IV-C2 $V_{t}$ Modulation Model

Stresses induced in the channel region cause changes to the band structure,
which result in fluctuations in the band-edge potentials, band-gap, and the
effective density of states [25]. Due to these shifts, the flatband voltage
$V_{FB}$ and channel surface potential $\psi_{s}$ change, causing a change in
the $V_{t}$. Threshold voltage is typically treated as a static parameter in
the BSIM-CMG model with variability parameters such as $DELVTRAND$ [22]
introduced to model the effect of layout etc. on the $V_{t}$. In this work,
the model used in [34] is adapted for calculating the shift in conduction
band-edge potential $\Delta E_{c}$ due to time-varying strains in the channel.
Changes to the valence band-edge potential $\Delta E_{v}$ are evaluated using
the model in [25]. Both conduction and valence band-edge potential shifts are
evaluated as:

$\Delta
E_{c}=\Xi_{d}(S_{x_{1}x_{1}}+S_{x_{2}x_{2}}+S_{x_{3}x_{3}})+\Xi_{u}S_{x_{3}x_{3}}$
(24a) $\Delta
E_{v}=u_{1}(S_{x_{1}x_{1}}+S_{x_{2}x_{2}}+S_{x_{3}x_{3}})+2u_{2}(S_{x_{3}x_{3}}-S_{x_{1}x_{1}})$
(24b)

---

---

Figure 9: (a) Plot of channel volume-averaged stress versus the equivalent
displacement as evaluated from frequency domain simulation in COMSOL (b) Plot
of channel volume-averaged strain versus the equivalent displacement.

while $S_{x_{1}x_{1}}$, $S_{x_{2}x_{2}}$ and $S_{x_{3}x_{3}}$ are time-
varying, channel-averaged strain components, $\Xi_{d}$, $\Xi_{u}$, $u_{1}$ and
$u_{2}$ are deformation potential constants (values given in Table III). The
shifts in the threshold voltage in the n- or p-type transistors are given by:

$q\Delta V_{thp}=(m-1)\Delta E_{c}-m\Delta E_{v}$ (25a) $q\Delta
V_{thn}=-m\Delta E_{c}+(m-1)\Delta E_{v}$ (25b)

where $m$ is the body effect parameter. This shift in the $V_{t}$ is
incorporated in the $fbsims$ model for the sense transistor alongside the
variability parameter $DELVTRAND$.

TABLE III: THRESHOLD VOLTAGE MODEL CONSTANTS [25][34]

PARAMETER | VALUE | UNIT
---|---|---
$\Xi_{d}$ | $1.13$ | $eV$
$\Xi_{u}$ | $9.16$ | $eV$
$u_{1}$ | $2.46$ | $eV$
$u_{2}$ | $-2.35$ | $eV$
$m$ | $1.3$ | $1$

---

---

---

Figure 10: (a) Testbench circuit schematic for evaluating the performance of
the fbsims sense transistor. The fbsims is used in conjunction with the fbsimd
and fRBT resonant body module. (b) Transient current output of the fbsims
transistor with a 10 mV drive at the fbsimd MOSCAP input. (c) Magnitude and
phase plots for the output current depicting a resonance at the frequency set
by the fRBT resonant body module.

####  IV-C3 $\nu_{vsat}$ Modulation Model

The effect of stress on the saturation velocity $\nu_{sat}$ can be modeled
using the methodology described in [25] and [35]. There is a direct
relationship between the $\nu_{sat}$ variation in highly scaled devices with
ballistic efficiency approaching 1 and variation in mobility due to stress,
through a factor $\alpha$. The $\nu_{sat}$ modulation is then given by:

Figure 11: (a) Complete model for the fRBT device comprising of the VerilogA
building blocks as well as parasitics extracted from the device layout. (b)
Testbench mimicking the RF measurement setup for testing the characteristics
of the extracted fRBT model.

$\nu_{sat}=\nu_{sat0}\left(1+\alpha\frac{\Delta\mu}{\mu}\right)$ (26)

where $\nu_{sat0}$ is the static saturation velocity and the constant $\alpha$
has a value of $0.85$. Source-drain current degradation due to $\nu_{sat}$ is
modeled in the BSIM-CMG model using the degradation parameter $D_{vsat}$ [22]
so the $\nu_{sat}$ modulation is added in conjunction to this parameter in the
$fbsims$ model.

The sense transistor module employing the $fbsims$ model uses stress and
strain values to calculate the changes in $\mu$, $V_{t}$ and $\nu_{sat}$. The
resonant behaviour of the device however, is captured using equivalent
displacement parameter through the resonant body module as described in
Section IV-B. Therefore, the equivalent displacement must be converted to
stress and strain components which can then be utilized. To evaluate the
relationship between the equivalent displacement and the stress/strain in the
channel, frequency domain simulation was carried out in COMSOL at the
resonance frequency while varying the drive stress level. The linear
relationship between the equivalent displacement and the channel stress/strain
is presented in Fig. 9. The constants of proportionality obtained from this
analysis are used to perform the equivalent displacement to stress/strain
conversion in the $fbsims$ model. We implement a testbench shown in Fig. 10 to
obtain transient and frequency sweep characteristics of the output current of
a $fbsims$ transistor. The $fbsimd$ drive MOSCAP and fRBT body modules are
added, and the mechanical terminals for all three components are connected to
the same mechanical node. When both the $fbsims$ and $fbsimd$ are biased in
strong inversion and a drive voltage of 10 mV is applied to the input, a
transient output waveform is observed as shown in Fig. 10(b). Under the same
biasing conditions, an ac simulation is performed to obtain the frequency
response of the ”fbsims” transistor. Fig. 10(c) shows the resonant
characteristic of the output current, as expected from the model when used
with the drive and body modules.

Once the individual building blocks of the fRBT model have been implemented
and tested for correct functionality, they are connected together in the same
configuration as the actual device as shown in Fig. 11(a). The differential
drive section of the model consists of $N_{t}$ $fbsimd$ drive transducers,
half of which are connected to one phase of the input drive voltage and half
to the opposite phase. The force contribution from the drive units is then
added to obtain the total drive force for the fRBT body module. To ensure that
the force from the differential drive transducers add in phase, 1:1 mechanical
transformers are implemented which invert the polarity of the mechanical
quantities being carried by the mechanical network. The mechanical node
$x_{disp}$ is common to all the modules in the fRBT model. The same mechanical
transformer is also used in conjunction with one of the $fbsims$ sense
transistors to make sure that the appropriate phase of channel stresses and
strains are generated. The model is incomplete without the inclusion of
parasitic capacitances and resistances associated with the metal traces, as
well as self and coupling capacitances of each net. The Calibre xACT™ tool is
used to perform the parasitic extraction on the layout of the fRBT device. The
extent of layout corresponds to the de-embedding plane used in the RF
measurements to extract device characteristics. Once a netlist is generated
after parasitic extraction, the standard transistors are replaced with their
$fbsim$ counterparts while keeping a track of the location of the each. The
netlist is also augmented to include the fRBT Body module and the mechanical
connections between each of the components in the model. Owing to the high
frequency of operation, EM extraction can also be performed for greater
accuracy. Since EM analysis is comparatively computationally expensive and
difficult to integrate with the rest of the model, it was excluded in this
iteration of the fRBT model development. In the testbench for the complete
extracted fRBT model depicted in Fig. 11(b), ideal baluns are used to handle
the differential output and input signals for the ease of calculating the
differential electromechanical transconductance $g_{m,em}$. The drive and
sense transistors are biased using large inductors to mimic the biasing
through bias-Ts in the measurement setup. The output is maintained at 0V DC
bias to extract the output current in the simulation.

##  V Results and Discussion

Figure 12: (a) Measurement setup for the fRBT DUT depicting the GSSG and DC
probe landing. The quantities $V_{in}^{+}$, $V_{in}^{-}$, $I_{out}^{+}$ and
$I_{out}^{-}$ represent a combination RF signals and DC levels at the
different ports of the DUT. (b) Comparison of measured data from the fRBT
device to the extracted fRBT model with variation in the sense transistor
drain voltage $V_{D,s}$ for fixed drive and sense gate bias voltage of 0.8 V.

On-chip RF measurements of the fRBT devices were conducted [36] using the
setup shown in Fig. 12(a). Prior to the measurement, Hybrid SOLR-LRRM[37]
calibration is carried out using impedance standard substrates. Once the
quality of the calibration was ascertained to be within recommended
tolerances, the DUT was connected to a Agilent® Parametric Network Analyzer
(PNA) through Cascade® GSSG probes landed on the DUT probe pads. The internal
bias-Ts of the PNA were connected to two Source-Measure Units (SMUs), one each
for the input gate drive bias (S3 and S3) and for the output sense transistor
drain bias (S1 and S2). A third SMU provided the sense transistors’ gate bias
using a DC probe. An input signal level of -10 dBm was used for all of the
measurements. Standard 4-port S-parameter measurement was performed at each
bias point and the resulting single-ended parameters were converted to mixed-
mode or differential parameters. Conversion of the differential S-parameters
to differential Y-parameters provides extraction of the electromechanical
transconductance of the fRBT DUT:

$|g_{m,em}|=|Y_{21dd}-Y_{12dd}|$ (27)

where $Y_{21dd}$ and $Y_{12dd}$ are differential $Y$-parameters. Open and
short structures present on-chip are used for de-embedding parasitic elements
associated with pads and routing up to the de-embedding plane specified in
Fig. 12(a). This includes bias-dependent capacitance of electrostatic
discharge (ESD) diodes necessary to protect the fRBTs. In practice, devices
would be routed directly to adjacent circuits in low-level metal layers,
irradicating the need for ESD diodes and extensive routing for each device.

Figure 13: Variation of the fRBT $g_{m,em}$ response with drive gate bias
$V_{G,d}$ for sense transistor drain ($V_{D,s}$) and gate ($V_{G,s}$) biased
at 0.8 V. Along with the measured data and the extracted model $g_{m,em}$, the
TM, SM1 and SM2 curves represent the individual resonant characteristics of
the target mode, spurious mode 1 and 2 respectively. The three modes
collectively combine to give the overall response matching the measured data.

In the first set of measurements, the drive MOSCAPs are biased into inversion
with $V_{G,drive}=0.8V$. The sense transistor gate bias is fixed at
$V_{G,s}=0.8V$ and the drain bias voltage is varied to assess the performance
in different regions of operation. From Fig. 12(b), we see that as $V_{d,s}$
is lowered from $0.8V$ (saturation regime) to $0.2V$ (linear regime), the peak
$g_{m,em}$ reduces in accordance with a reduction in DC drain current. In each
of the three measured curves, two spurious modes in the vicinity of the
desired resonance peak can be observed, one of which is merged in close
proximity with the target mode. This likely results from the finite nature of
the resonance cavity in the $x_{2}$ direction. Harmonics generated by
transverse modes such as these cannot be captured by unit-cell-based 3D FEM
models, and would require extensive modeling of the complete 3D structure
which is currently computationally prohibitive. The Q of the main resonance
peak is extracted and fed back to the fRBT model as an initial point for the
final model fit. As can be seen, the model captures the presence of the two
spurious modes along with the targeted mode. Some discrepancies are observed
between the measured response and the model in the region between the spurious
modes and in a transmission zero beyond the measured main resonance. These can
be attributed to limitations of the parasitic extraction which excludes some
of the coupling capacitance across the fRBT structure. The model shows
commensurate variation with $V_{d,s}$ as the measured data, validating that
the current sensing mechanism modeling is able to accurately capture the
effect of bias variation.

An important feature of electrostatically-driven electromecahnical devices is
the ability to control the strength of the drive transducer with DC bias. The
resonance is not completely attenuated since a capacitance still exists when
the bias goes down to $0$ $V$, which results in transduction. As can be seen
from the corresponding model curve for the $V_{G,d}=0$ $V$ bias, the
$g_{m,em}$ is attenuated.

---

---

Figure 14: (a) PSS and PAC testbench circuit schematic for emulating the mixer
measurement setup typically used for resonators embedded in high feedthrough.
(b) PAC analysis results showing the presence of the resonance characteristic
of the output current in the sideband corresponding to $f_{0}=f_{RF}+f_{LO}$.

The developed fRBT model is also compatible with Periodic Steady State (PSS),
Periodic AC (PAC) and Harmonic Balance simulation scenarios. Since the model
is designed to capture the nonlinear mechanisms inherent to the device, we can
use these simulation techniques to verify functionality. An RF/LO-based mixing
measurement technique is used to extract the performance of resonators
embedded in high feedthrough [38]. In this method, along with an RF signal
$v_{RF}$ that is applied at the device input, a lower frequency signal
$v_{LO}$ is superimposed on the bias voltage. The frequencies of the RF and LO
signals are set such that $f_{0}=f_{RF}+f_{LO}$ where $f_{0}$ is the resonance
frequency of the resonator. Even though frequencies different from the $f_{0}$
are applied to the resonator terminals, due to the nonlinear electromechanical
transduction based up-conversion mixing, a mechanical force is generated at
the resonance frequency $f_{0}$. To test the functionality of the fRBT model
under these nonlinear simulation conditions that are frequently used for
oscillator/filter designs, a simulation setup emulating resonator mixing
measurements is created as shown in Fig. 14(a). An RF signal $v_{RF}$ of
frequency $f_{RF}$ is applied to the gate of the drive transistor along with
the gate bias $V_{G,d}$ of 0.8 V. Instead of clamping the source-drain to
ground, a LO signal $v_{LO}$ of frequency $f_{LO}$ is applied as shown. PSS
and PAC simulation results corresponding to the output current are obtained as
shown in Fig. 14(b). It can be seen that, in the sideband corresponding to the
frequencies 11-12.5 GHz the resonance characteristic is exhibited while in the
6 to 7.5 GHz frequency range it is not. This example confirms that the device
model is able to capture the effect of MOSCAP drive nonlinearity as
ascertained by the PSS-PAC simulation results.

##  VI Conclusion

A compact model for an 11.8 GHz Fin Resonant Body Transistor fabricated in a
commercial 14nm Fin-FET process has been presented, which captures for the
first time all prominent device characteristics including electrical and
mechanical nonlinearities. Analysis has also been provided of the unique
nature of the phononic dispersion in the device owing to the presence of the
BEOL PnC, via a theoretical framework and FEM simulations. An Eigenfrequency-
based methodology amenable to dielectric transduction in the fRBT for the
extraction of equivalent mechanical circuit parameters for the main resonant
cavity has also been described, and can be readily adapted to alternate device
geometries and resonance modes. The fRBT model presented in this work is fully
parameterized and flexible allowing for the addition of more drive/sense
transistors, and also captures all major sources of device non-linearity. The
model is also completely integrated within the simulation framework alongside
the foundry-supplied PDK for the technology. Future improvements to the model
include temperature variation and more accurate parasitic extraction. This
large-signal model enables direct integration with interface circuitry,
necessary for design of RF and mm-wave oscillators and provides more accurate
predictions of overall performance during the design and simulation phase of
systems employing the fRBT for electromechanical signal processing.

## References

  * [1] G. Chen and M. Rinaldi, “Aluminum nitride combined overtone resonators for the 5g high frequency bands,” _Journal of Microelectromechanical Systems_ , vol. 29, no. 2, pp. 148–159, 2020.
  * [2] Y. Yang, R. Lu, L. Gao, and S. Gong, “10–60-ghz electromechanical resonators using thin-film lithium niobate,” _IEEE Transactions on Microwave Theory and Techniques_ , vol. 68, no. 12, pp. 5211–5220, 2020.
  * [3] Y. Liu, J. Kuo, B. Davaji, V. Gund, A. Lal, J. Sharma, and N. Singh, “Characterization of alscn on cmos,” in _2020 Joint Conference of the IEEE International Frequency Control Symposium and International Symposium on Applications of Ferroelectrics (IFCS-ISAF)_ , 2020, pp. 1–5.
  * [4] M. Riverola, G. Sobreviela, F. Torres, A. Uranga, and N. Barniol, “Single-resonator dual-frequency beol-embedded cmos-mems oscillator with low-power and ultra-compact tia core,” _IEEE Electron Device Letters_ , vol. 38, no. 2, pp. 273–276, 2017.
  * [5] S. Li, “Cmos-mems resonators and their applications,” in _2013 Joint European Frequency and Time Forum International Frequency Control Symposium (EFTF/IFC)_ , 2013, pp. 915–921.
  * [6] Huikai Xie and G. K. Fedder, “Fabrication, characterization, and analysis of a drie cmos-mems gyroscope,” _IEEE Sensors Journal_ , vol. 3, no. 5, pp. 622–631, 2003.
  * [7] C. Sun, M. Tsai, Y. Liu, and W. Fang, “Implementation of a monolithic single proof-mass tri-axis accelerometer using cmos-mems technique,” _IEEE Transactions on Electron Devices_ , vol. 57, no. 7, pp. 1670–1679, 2010\.
  * [8] B. Sadhu, Y. Tousi, J. Hallin, S. Sahl, S. K. Reynolds, . Renström, K. Sjögren, O. Haapalahti, N. Mazor, B. Bokinge, G. Weibull, H. Bengtsson, A. Carlinger, E. Westesson, J. Thillberg, L. Rexberg, M. Yeck, X. Gu, M. Ferriss, D. Liu, D. Friedman, and A. Valdes-Garcia, “A 28-ghz 32-element trx phased-array ic with concurrent dual-polarized operation and orthogonal phase and gain control for 5g communications,” _IEEE Journal of Solid-State Circuits_ , vol. 52, no. 12, pp. 3373–3391, 2017.
  * [9] B. Bahr, R. Marathe, and D. Weinstein, “Theory and Design of Phononic Crystals for Unreleased CMOS-MEMS Resonant Body Transistors,” _Journal of Microelectromechanical Systems_ , vol. 24, no. 5, pp. 1520–1533, Oct. 2015.
  * [10] B. Bahr, Y. He, Z. Krivokapic, S. Banna, and D. Weinstein, “32GHz resonant-fin transistors in 14nm FinFET technology,” in _2018 IEEE International Solid - State Circuits Conference - (ISSCC)_ , Feb. 2018, pp. 348–350.
  * [11] A. Srivastava, B. Chatterjee, U. Rawat, Y. He, D. Weinstein, and S. Sen, “Analysis and design considerations for achieving the fundamental limits of phase noise in mmwave oscillators with on-chip mems resonator,” _IEEE Transactions on Circuits and Systems II: Express Briefs_ , pp. 1–1, 2020.
  * [12] R. Marathe, W. Wang, and D. Weinstein, “Si-based unreleased hybrid mems-cmos resonators in 32nm technology,” in _2012 IEEE 25th International Conference on Micro Electro Mechanical Systems (MEMS)_ , 2012, pp. 729–732.
  * [13] B. W. Bahr, D. Weinstein, and L. Daniel, “Unreleased 1D CMOS Resonant Body Transistor with MIT Virtual Source (URBT-MVS) Model,” March 2016. [Online]. Available: https://nanohub.org/publications/132/1
  * [14] B. Bahr, D. Weinstein, and L. Daniel, “Released Resonant Body Transistor (RBT) Model,” March 2016. [Online]. Available: https://nanohub.org/publications/17/supportingdocs?v=1
  * [15] J. P. Duarte, S. Khandelwal, A. Medury, C. Hu, P. Kushwaha, H. Agarwal, A. Dasgupta, and Y. S. Chauhan, “Bsim-cmg: Standard finfet compact model for advanced circuit design,” in _ESSCIRC Conference 2015 \- 41st European Solid-State Circuits Conference (ESSCIRC)_ , 2015, pp. 196–201.
  * [16] P. Steeneken, J. Ruigrok, S. Kang, J. van Beek, J. Bontemps, and J.-J. Koning, “Parameter extraction and support-loss in mems resonators,” _Proc. Comsol conference 2007_ , p. 725, 10 2007.
  * [17] J. Vasseur, _The Plane Wave Expansion Method_. John Wiley & Sons, Ltd, 2019, ch. 5, pp. 107–141.
  * [18] B. Bahr, “Monolithically integrated mems resonators and oscillators in standard ic technology,” Ph.D. dissertation, Massachusetts Institute of Technology, Boston, June 2016.
  * [19] J. E.-Y. Lee and A. A. Seshia, “Direct parameter extraction in feedthrough-embedded capacitive mems resonators,” _Sensors and Actuators A: Physical_ , vol. 167, no. 2, pp. 237–244, 2011.
  * [20] D. Weinstein and S. A. Bhave, “Internal dielectric transduction in bulk-mode resonators,” _Journal of Microelectromechanical Systems_ , vol. 18, no. 6, pp. 1401–1408, 2009.
  * [21] D. Weinstein and S. A. Bhave, “The resonant body transistor,” _Nano Letters_ , vol. 10, no. 4, pp. 1234–1237, 2010.
  * [22] S. Khandelwal, J. Duarte, S. V., N. Paydavosi, D. D. Lu, C.-H. Lin, M. Dunga, S. Yao, T. H. Morshed, A. Niknejad, S. Salahuddin, and C. Hu, “Bsim-cmg 108.0.0: Multi-gate mosfet compact model: technical manual,” in _Berkeley: University of California_ , 2014.
  * [23] S. D. Senturia, _Microsystem Design_. USA: Kluwer Academic Publishers, 2001.
  * [24] Y. Yang, R. Lu, and S. Gong, “High Q antisymmetric mode lithium niobate mems resonators with spurious mitigation,” _Journal of Microelectromechanical Systems_ , vol. 29, no. 2, pp. 135–143, 2020.
  * [25] X. Li, Z. Ye, Y. Tan, and Y. Wang, “A two-dimensional analysis method on STI-aware layout-dependent stress effect,” _IEEE Transactions on Electron Devices_ , vol. 59, no. 11, pp. 2964–2972, 2012.
  * [26] H. Tsuno, K. Anzai, M. Matsumura, S. Minami, A. Honjo, H. Koike, Y. Hiura, A. Takeo, W. Fu, Y. Fukuzaki, M. Kanno, H. Ansai, and N. Nagashima, “Advanced analysis and modeling of mosfet characteristic fluctuation caused by layout variation,” in _2007 IEEE Symposium on VLSI Technology_ , 2007, pp. 204–205.
  * [27] D. C. Chen, G. S. Lin, T. H. Lee, R. Lee, Y. C. Liu, M. F. Wang, Y. C. Cheng, and D. Y. Wu, “Compact modeling solution of layout dependent effect for finfet technology,” in _Proceedings of the 2015 International Conference on Microelectronic Test Structures_ , 2015, pp. 110–115.
  * [28] Ji-Song Lim, S. E. Thompson, and J. G. Fossum, “Comparison of threshold-voltage shifts for uniaxial and biaxial tensile-stressed n-mosfets,” _IEEE Electron Device Letters_ , vol. 25, no. 11, pp. 731–733, 2004.
  * [29] H. Nam and C. Shin, “Impact of current flow shape in tapered (versus rectangular) finfet on threshold voltage variation induced by work-function variation,” _IEEE Transactions on Electron Devices_ , vol. 61, no. 6, pp. 2007–2011, 2014.
  * [30] Masumi Saitoh, Akio Kaneko, Kimitoshi Okano, Tomoko Kinoshita, Satoshi Inaba, Yoshiaki Toyoshima, and Ken Uchida, “Three-dimensional stress engineering in finfets for mobility/on-current enhancement and gate current reduction,” in _2008 Symposium on VLSI Technology_ , 2008, pp. 18–19.
  * [31] Y. Kanda, “A graphical representation of the piezoresistance coefficients in silicon,” _IEEE Transactions on Electron Devices_ , vol. 29, no. 1, pp. 64–70, 1982.
  * [32] R. C. Jaeger, J. C. Suhling, R. Ramani, A. T. Bradley, and Jianping Xu, “CMOS stress sensors on [100] silicon,” _IEEE Journal of Solid-State Circuits_ , vol. 35, no. 1, pp. 85–95, 2000.
  * [33] S. K. Marella, “Performance variations due to layout-dependent stress in vlsi circuits,” Ph.D. dissertation, University of Minnesota, May 2015.
  * [34] S. K. Marella, A. R. Trivedi, S. Mukhopadhyay, and S. S. Sapatnekar, “Optimization of finfet-based circuits using a dual gate pitch technique,” in _2015 IEEE/ACM International Conference on Computer-Aided Design (ICCAD)_ , 2015, pp. 758–763.
  * [35] A. Lochtefeld and D. A. Antoniadis, “Investigating the relationship between electron mobility and velocity in deeply scaled nmos via mechanical stress,” _IEEE Electron Device Letters_ , vol. 22, no. 12, pp. 591–593, 2001\.
  * [36] J. Anderson, Y. He, B. Bahr, and D. Weinstein, “X-band fin resonant body transistors in 14nm cmos technology,” 2021.
  * [37] “Hybrid Calibration for 4-Port On-wafer Probing,” 12-Aug 2020. [Online]. Available: https://www.formfactor.com/blog/2019/hybrid-calibration-for-4-port-on-wafer-probing/
  * [38] J. R. Clark, W. . Hsu, M. A. Abdelmoneum, and C. T. . Nguyen, “High-q uhf micromechanical radial-contour mode disk resonators,” _Journal of Microelectromechanical Systems_ , vol. 14, no. 6, pp. 1298–1310, 2005.

|  Udit Rawat received the B.Tech (Hons.) degree from Uttarakhand Technical
University, India in 2012 and M.S. by Research degree from Indian Institute of
Technology Madras, Chennai, India in 2015. He is currently pursuing a Ph.D.
degree with the Department of Electrical Engineering, Purdue University, West
Lafayette, IN, USA. From 2015 to 2017 he worked as a Digital Design Engineer
at Intel Corporation, Bangalore, India. He is currently a research assistant
with the Hybrid MEMS Research group, Purdue University. His research interests
include design, fabrication and modeling of RF MEMS resonators as well as
their integration in standard CMOS technology. He is also interested in the
design of unreleased CMOS integrated gyroscopes and other physical sensors
along with their interface circuit design.
---|---
|  Bichoy Bahr received the BSc and MSc degrees with honors in 2008 and 2012,
both in electrical engineering, from Ain Shams University, Cairo, Egypt. He
received the PhD degree from Massachusetts Institute of Technology (MIT) in
2016. He has been with the HybridMEMS group at MIT from 2012 to 2016, where he
developed unreleased CMOS-MEMS resonators and monolithic GaN MMIC oscillators.
Dr. Bahr joined Texas Instruments’ Kilby Labs, Dallas, TX in 2016. He has been
elected to TI’s prestigious tech Ladder as Member Group, Technical Staff
(MGTS) in 2019. His research interests include MEMS resonators, phononic
crystals, oscillators, ultrasonic transducers, and numerical optimization.
---|---
|  Dana Weinstein received her B.A. in Physics and Astrophysics from
University of California - Berkeley in 2004 and her Ph.D. in Applied Physics
in 2009 from Cornell, working on multi-GHz MEMS. Prior to joining Purdue as an
Associate Professor in Electrical and Computer Engineering in 2015, Dr.
Weinstein joined the Department of Electrical Engineering and Computer Science
at MIT as an Assistant Professor, and served as an Associate Professor there
between 2013 and 2015. Dana is now a Purdue Faculty Scholar and serves as
Associate Dean of Graduate Education in the College of Engineering. She is the
recipient of the Transducers Early Career Award, NSF CAREER Award, the DARPA
Young Faculty Award, the Intel Early Career Award, and the IEEE IEDM Roger A.
Haken Best Paper Award.
---|---


Let $u=b_{j}X_{S}^{T}\sigma_{S-1}\left(X_{S-2}^{T}\sigma_{S-2}\left(\cdots
X_{1}^{T}a_{j}\right)\right)$, then we compute the gradient of the function
$F_{j}$ (the $j$-data $(a_{j},b_{j})$) with regard to the weight matrix of
each layers:

$\displaystyle\frac{\partial F_{j}}{\partial X_{1}[i]}$
$\displaystyle=\frac{-\exp(-u)b_{j}}{1+\exp(-u)}\cdot
X_{S}^{T}\partial\sigma_{S-1}\left(\cdots\right)\cdots\left(X_{2}[:][i]\partial\sigma_{1}(X_{1}^{T}[i]a_{j})\right)a_{j}+\lambda
X_{1}[i]$ $\displaystyle\frac{\partial F_{j}}{\partial X_{2}[i]}$
$\displaystyle=\frac{-\exp(-u)b_{j}}{1+\exp(-u)}\cdot
X_{S}^{T}\partial\sigma_{S-1}(\cdots)\cdots\left(X_{3}[:][i]\partial\sigma_{2}(X_{2}[i]^{T}\sigma_{1}(X_{1}^{T}a_{j}))\right)\sigma_{1}(X_{1}^{T}a_{j})+\lambda
X_{2}[i]$ $\displaystyle\vdots$ $\displaystyle\frac{\partial F_{j}}{\partial
X_{S-1}[i]}$ $\displaystyle=\frac{\exp(-u)}{1+\exp(-u)}\cdot-
b_{j}X_{S}[i]^{T}\partial\sigma_{S-1}\left(X_{S-1}[i]^{T}\cdots\left(X_{1}^{T}a_{j}\right)\right)\sigma_{S-2}(X_{S-2}^{T}\cdots\sigma_{1}\left(X_{1}^{T}a_{j}\right))+\lambda
X_{S-1}[i]$ $\displaystyle\frac{\partial F_{j}}{\partial X_{S}}$
$\displaystyle=\frac{\exp(-u)}{1+\exp(-u)}\cdot-
b_{j}\sigma_{S-1}\left(X_{S-1}^{T}\sigma_{S-1}\cdots\left(X_{1}^{T}a_{j}\right)\right)+\lambda
X_{S}.$

The inner product of $X$ and $\nabla F_{j}(X)$ is

$\displaystyle\left\langle X,\nabla F_{j}(X)\right\rangle$
$\displaystyle=\sum_{i=1}^{m_{1}}\left\langle X_{1}[i],\frac{\partial
F_{j}}{\partial X_{1}[i]}\right\rangle+\sum_{i=1}^{m_{2}}\left\langle
X_{2}[i],\frac{\partial F_{j}}{\partial
X_{2}[i]}\right\rangle+\cdots+\sum_{i=1}^{m_{S-1}}\left\langle
X_{S-1}[i],\frac{\partial F_{j}}{\partial
X_{S-1}[i]}\right\rangle+\left\langle X_{S},\frac{\partial F_{j}}{\partial
X_{S}}\right\rangle$

Let $\sigma_{1}(x)=\sigma_{2}(x)=\cdots=\sigma_{S-1}(x)=\max(0,x)$, then their
gradients $\partial\sigma_{1}(x)=\partial\sigma_{2}(x)=\cdots=\chi_{x\geq 0}$.
Then the activation functions with their gradient satisfy that
$\sigma_{i}(x)=x\partial\sigma_{i}(x)$ for $i=1,2,\cdots,S-1$. Then

$\displaystyle\sum_{i=1}^{m_{1}}\left\langle X_{1}[i],\frac{\partial
F_{j}(X)}{\partial X_{1}[i]}\right\rangle$
$\displaystyle=\frac{-\exp(-u)b_{j}}{1+\exp(-u)}\cdot
X_{S}^{T}\partial\sigma_{S-1}\left(\cdots\right)X_{S-1}^{T}\partial\sigma_{S-1}(\cdots)\cdots\sum_{i=1}^{m_{1}}X_{2}[:][i]\partial\sigma_{1}(X_{1}[i]^{T}a_{j})X_{1}[i]^{T}a_{j}+\lambda\sum_{i=1}^{m_{1}}\left\|X_{1}[i]\right\|^{2}$
$\displaystyle=\frac{-\exp(-u)b_{j}}{1+\exp(-u)}\cdot
X_{S}^{T}\partial\sigma_{S-1}\left(\cdots\right)X_{S-1}^{T}\partial\sigma_{S-1}(\cdots)\cdots\sum_{i=1}^{m_{1}}X_{2}[:][i]\sigma_{1}(X_{1}[i]^{T}a_{j})+\lambda\left\|X_{1}\right\|^{2}$
$\displaystyle=\frac{-\exp(-u)b_{j}}{1+\exp(-u)}\cdot
X_{S}^{T}\partial\sigma_{S-1}\left(\cdots\right)X_{S-1}^{T}\partial\sigma_{S-1}(\cdots)\cdots\partial\sigma_{2}(X_{2}^{T}\sigma_{1}(X_{1}^{T}a_{j}))X_{2}^{T}\sigma_{1}(X_{1}^{T}a_{j})+\lambda\left\|X_{1}\right\|^{2}$
$\displaystyle=\frac{\exp(-u)}{1+\exp(-u)}\cdot-u+\lambda\left\|X_{1}\right\|^{2}$
(53) $\displaystyle\sum_{i=1}^{m_{2}}\left\langle X_{2}[i],\frac{\partial
F_{j}(X)}{\partial X_{2}[i]}\right\rangle$
$\displaystyle=\frac{-\exp(-u)b_{j}}{1+\exp(-u)}\cdot
X_{S}^{T}\partial\sigma_{S-1}(\cdots)\cdots\sum_{i=1}^{m_{2}}\left(X_{3}[:][i]\partial\sigma_{2}(X_{2}[i]^{T}\sigma_{1}(X_{1}^{T}a_{j}))\right)X_{2}[i]^{T}\sigma_{1}(X_{1}^{T}a_{j})+\lambda\sum_{i=1}^{m_{2}}\left\|X_{2}[i]\right\|^{2}$
$\displaystyle=\frac{-\exp(-u)b_{j}}{1+\exp(-u)}\cdot
X_{S}^{T}\partial\sigma_{S-1}(\cdots)\cdots\sum_{i=1}^{m_{2}}\left(X_{3}[:][i]\sigma_{2}(X_{2}[i]^{T}\sigma_{1}(X_{1}^{T}a_{j}))\right)+\lambda\left\|X_{2}\right\|^{2}$
$\displaystyle=\frac{-\exp(-u)b_{j}}{1+\exp(-u)}\cdot
X_{S}^{T}\partial\sigma_{S-1}(\cdots)\cdots\partial\sigma_{3}(X_{3}^{T}\sigma_{2}(X_{1}^{T}a_{j}))X_{3}^{T}\sigma_{2}(X_{2}^{T}\sigma_{1}(X_{1}^{T}a_{j}))+\lambda\left\|X_{2}\right\|^{2}$
$\displaystyle=\frac{\exp(-u)}{1+\exp(-u)}\cdot-u+\lambda\left\|X_{2}\right\|^{2}.$
(54)

$\vdots$

$\displaystyle\left\langle X_{S},\frac{\partial F_{j}}{\partial
X_{S}}\right\rangle$ $\displaystyle=\frac{\exp(-u)}{1+\exp(-u)}\cdot-
b_{j}X_{S}^{T}\sigma_{S-1}\left(X_{S-1}^{T}\sigma_{S-1}\cdots\left(X_{1}^{T}a_{j}\right)\right)+\lambda\left\|X_{S}\right\|^{2}$
$\displaystyle=\frac{\exp(-u)}{1+\exp(-u)}\cdot-u+\lambda\left\|X_{S}\right\|^{2}.$
(55)

Incorporating the above results, we can achieve that

$\displaystyle\left\langle X,\nabla
F_{j}(X)\right\rangle=\lambda\left\|X\right\|^{2}-\frac{Su\exp(-u)}{1+\exp(-u)}.$

For any $u$, we can estimate

$\displaystyle\frac{u\exp(-u)}{1+\exp(-u)}=\frac{u}{1+\exp(u)}\leq 1.$ (56)

Finally, we get that

$\displaystyle\left\langle X,\nabla
F_{j}(X)\right\rangle\geq\lambda\left\|X\right\|^{2}-S.$ (57)

###  D.2 Applications in Section 4

    * –

Logistic regression with $\ell_{1}$ regularizer:
$f(x)=\sum_{i=1}^{n}\log(1+\exp(-b_{i}a_{i}^{x}))+\lambda\left\|x\right\|_{1}$.
The gradient is

$\displaystyle\nabla
f(x)=\frac{-b_{i}a_{i}\exp(-b_{i}a_{i}^{T}x)}{1+\exp(-b_{i}a_{i}^{T}x)}+\lambda\mbox{sign}(x)$
(58)

then

$\displaystyle\left\langle\nabla
f(x),x\right\rangle=-b_{i}a_{i}^{T}x\frac{\exp(-b_{i}a_{i}^{T}x)}{1+\exp(-b_{i}a_{i}^{T}x)}+\lambda\left\|x\right\|_{1}=\frac{u\exp(u)}{1+\exp(u)}+\lambda\left\|x\right\|_{1}$
(59)

If $u<0$, then we get
$\left|\frac{u\exp(u)}{1+\exp(u)}\right|\leq\frac{-u}{2\exp(-u)}\leq\frac{1}{2}$;
else if $u>0$, then $\frac{u\exp(u)}{1+\exp(u)}>0$. Therefore, no matter $u>0$
or not, we both have that

$\displaystyle\left\langle\nabla
f(x),x\right\rangle\geq\lambda\left\|x\right\|_{1}-\frac{1}{2}\geq\lambda\left\|x\right\|-\frac{1}{2}$
(60)

    * –

$S$-Layer Neural networks with a $\ell_{1}$ regularizer: replace the
$\ell_{2}$ regularizer of Section D.1 by $\ell_{1}$, we can derive that

$\displaystyle\left\langle\nabla
f(x),x\right\rangle\geq\lambda\left\|x\right\|_{1}-S.$ (61)

##  Appendix E Objective Function is Dissipative on its Domain

######  Theorem E.1.

Under the condition of Lemma A.1, if $f$ is $L$-smooth and $R$-dissipative for
all $x\in\mathbb{R}^{d}$. Consider the SGD algorithm with the total number of
iterations $T\geq 1$, for any step-size
$\eta_{k}\leq\theta_{1}/((1+\rho)L^{2})$, we have

$\displaystyle\mathbb{E}[\left\|x_{T+1}-x^{\ast}\right\|^{2}]\leq\Pi_{k=1}^{T}\left(1-\eta_{k}\theta_{1}\right)\left\|x_{1}-x^{\ast}\right\|^{2}+\sum_{k=1}^{T}\left(\theta_{2}+\eta_{k}\sigma^{2}\right)\eta_{k}\cdot\Pi_{s>k}^{T}\left(1-\eta_{s}\theta_{1}\right).$
(62)

###### Proof.

Applying the $L$-smooth assumption, for any $x,y\in\mathbb{R}^{d}$, we have

$\displaystyle\left\|\nabla f(y)-\nabla f(x)\right\|$ $\displaystyle\leq
L\left\|y-x\right\|,$ $\displaystyle f(y)+\left\langle\nabla
f(x),y-x\right\rangle-\frac{L}{2}\left\|y-x\right\|^{2}$ $\displaystyle\leq
f(y)\leq f(x)+\left\langle\nabla
f(x),y-x\right\rangle+\frac{L}{2}\left\|y-x\right\|^{2}.$

Let $x=x^{\ast}$. By $\nabla f(x^{\ast})=0$, we get that

$\displaystyle\left\|\nabla f(y)\right\|$ $\displaystyle\leq
L\left\|x-x^{\ast}\right\|,$ (63a) $\displaystyle f(y)-f^{\ast}$
$\displaystyle\leq\frac{L}{2}\left\|y-x^{\ast}\right\|^{2}.$ (63b)

In this case, $R=0$, it means that $\left\langle x_{k}-x^{\ast},\nabla
f(x_{k})\right\rangle\geq\theta_{1}\left\|x_{k}-x^{\ast}\right\|^{2}-\theta_{2}$
for all $x\in\mathbb{R}^{d}$. Applying this condition and (63a) into Lemma
A.1, we have

$\displaystyle\mathbb{E}[\left\|x_{k+1}-x^{\ast}\right\|^{2}]$
$\displaystyle\leq\left\|x_{k}-x^{\ast}\right\|^{2}-2\eta_{k}\theta_{1}\left\|x_{k}-x^{\ast}\right\|^{2}+\theta_{2}\eta_{k}+\eta_{k}^{2}\sigma^{2}+\eta_{k}^{2}(1+\rho)L^{2}\left\|x_{k}-x^{\ast}\right\|^{2}$
$\displaystyle=\left(1-2\eta_{k}\theta_{1}+\eta_{k}^{2}(1+\rho)L^{2}\right)\left\|x_{k}-x^{\ast}\right\|^{2}+\left(\theta_{2}+\eta_{k}\sigma^{2}\right)\eta_{k}$

Letting $\eta_{k}\leq\frac{\theta_{1}}{(\rho+1)L^{2}}$ and , then
$0<1-2\eta_{k}\theta_{1}+\eta_{k}^{2}(1+\rho)L^{2}\leq 1-\eta_{k}\theta_{1}$.
Consider the above recursion from $k=1$ to $T$, we can achieve that

$\displaystyle\mathbb{E}[\left\|x_{T+1}-x^{\ast}\right\|^{2}]\leq\Pi_{k=1}^{T}\left(1-\eta_{k}\theta_{1}\right)\left\|x_{1}-x^{\ast}\right\|^{2}+\sum_{k=1}^{T}\left(\theta_{2}+\eta_{k}\sigma^{2}\right)\eta_{k}\cdot\Pi_{s>k}^{T}\left(1-\eta_{s}\theta_{1}\right).$

Therefore, the proof is complete. ∎

######  Theorem E.2.

(Constant step-size) Under the same conditions of Theorem E.1. If the step-
size is a constant $\eta_{k}=\eta\leq\theta_{1}/((1+\rho)L^{2})$, then for any
$T\geq 1$, we have

$\displaystyle\mathbb{E}[\left\|x_{T+1}-x^{\ast}\right\|^{2}]\leq(1-\eta\theta_{1})^{T}\left\|x_{1}-x^{\ast}\right\|^{2}+\frac{\left(\theta_{2}+\eta\sigma^{2}\right)}{\theta_{1}}.$

###### Proof.

$\displaystyle\mathbb{E}[\left\|x_{T+1}-x^{\ast}\right\|^{2}]$
$\displaystyle\leq\Pi_{k=1}^{T}\left(1-\eta_{k}\theta_{1}\right)\left\|x_{1}-x^{\ast}\right\|^{2}+\sum_{k=1}^{T}\left(\theta_{2}+\eta_{k}\sigma^{2}\right)\eta_{k}\cdot\Pi_{s>k}^{T}\left(1-\eta_{s}\theta_{1}\right)$
$\displaystyle=(1-\eta\theta_{1})^{T}\left\|x_{1}-x^{\ast}\right\|^{2}+\left(\theta_{2}+\eta\sigma^{2}\right)\eta\sum_{k=1}^{T}(1-\eta\theta_{1})^{T-k}$
$\displaystyle=(1-\eta\theta_{1})^{T}\left\|x_{1}-x^{\ast}\right\|^{2}+\left(\theta_{2}+\eta\sigma^{2}\right)\eta\cdot\frac{1-(1-\eta\theta_{1})^{T}}{\eta\theta_{1}}$
$\displaystyle=(1-\eta\theta_{1})^{T}\left\|x_{1}-x^{\ast}\right\|^{2}+\frac{\left(\theta_{2}+\eta\sigma^{2}\right)}{\theta_{1}}\cdot\left(1-(1-\eta\theta_{1})^{T}\right)$
$\displaystyle\leq(1-\eta\theta_{1})^{T}\left\|x_{1}-x^{\ast}\right\|^{2}+\frac{\left(\theta_{2}+\eta\sigma^{2}\right)}{\theta_{1}}.$

∎

######  Theorem E.3.

(Time-dependent step-size) Under the same conditions of Theorem E.1. If the
step-size is time dependent on the iteration where
$\eta_{k}\leq\theta_{1}/((1+\rho)L^{2})$, then for any $T\geq 1$, we have

$\displaystyle\mathbb{E}[\left\|x_{T+1}-x^{\ast}\right\|^{2}]$
$\displaystyle\leq\exp\left(-\theta_{1}\sum_{k=1}^{T}\eta_{k}\right)\left\|x_{1}-x^{\ast}\right\|^{2}+\sum_{k=1}^{T}\left(\theta_{2}+\eta_{k}\sigma^{2}\right)\eta_{k}\cdot\exp\left(-\theta_{1}\sum_{s=k+1}^{T}\eta_{s}\right)$

  * •

$\eta_{k}=\eta_{1}/k^{p}$ for $p\in(0,1]$

    * –

If $p\in(0,1)$, we have

$\displaystyle\mathbb{E}[\left\|x_{T+1}-x^{\ast}\right\|^{2}]\leq\exp\left(-\frac{\theta_{1}\eta_{1}}{1-p}\left((T+1)^{1-p}-1\right)\right)\left\|x_{1}-x^{\ast}\right\|^{2}+\frac{\left(\theta_{2}+\eta_{1}\sigma^{2}\right)2^{p}}{\theta_{1}}.$

    * –

$p=1$, then

$\displaystyle\mathbb{E}[\left\|x_{T+1}-x^{\ast}\right\|^{2}]\leq\frac{1}{(T+1)^{\theta_{1}\eta_{1}}}\left\|x_{1}-x^{\ast}\right\|^{2}+\frac{2\left(\theta_{2}+\eta_{1}\sigma^{2}\right)}{\theta_{1}}.$

  * •

Bandwidth-based step-sizes including step-decay band and polynomial band

    * –

$m/k^{p}\leq\eta_{k}\leq M/k^{p}$

      * *

$p\in(0,1)$,

$\displaystyle\mathbb{E}[\left\|x_{T+1}-x^{\ast}\right\|^{2}]\leq\exp\left(-\frac{\theta_{1}m}{1-p}\left((T+1)^{1-p}-1\right)\right)\left\|x_{1}-x^{\ast}\right\|^{2}+\frac{\left(\theta_{2}+M\sigma^{2}\right)2^{p}M}{m\theta_{1}}.$

      * *

$p=1$

$\displaystyle\mathbb{E}[\left\|x_{T+1}-x^{\ast}\right\|^{2}]\leq\frac{1}{(T+1)^{\theta_{1}m}}\left\|x_{1}-x^{\ast}\right\|^{2}+\frac{2\left(\theta_{2}+M\sigma^{2}\right)M}{m\theta_{1}}.$

    * –

$m\delta(k)\leq\eta_{i}^{k}\leq M\delta(k)$ where $\delta(k)=1/\alpha^{k-1}$
for $i\in[S]$ and $1\leq k\leq N$

###### Proof.

$\displaystyle\mathbb{E}[\left\|x_{T+1}-x^{\ast}\right\|^{2}]$
$\displaystyle\leq\Pi_{k=1}^{T}\left(1-\eta_{k}\theta_{1}\right)\left\|x_{1}-x^{\ast}\right\|^{2}+\sum_{k=1}^{T}\left(\theta_{2}+\eta_{k}\sigma^{2}\right)\eta_{k}\cdot\Pi_{s>k}^{T}\left(1-\eta_{s}\theta_{1}\right)$
$\displaystyle=\exp\left(-\theta_{1}\sum_{k=1}^{T}\eta_{k}\right)\left\|x_{1}-x^{\ast}\right\|^{2}+\sum_{k=1}^{T}\left(\theta_{2}+\eta_{k}\sigma^{2}\right)\eta_{k}\cdot\exp\left(-\theta_{1}\sum_{s=k+1}^{T}\eta_{s}\right)$
(64)

If $\eta_{k}=\eta_{1}/k^{p}$ for $p\in(0,1)$, we have

$\displaystyle\sum_{k=1}^{T}\eta_{k}$
$\displaystyle\geq\eta_{1}\int_{k=1}^{T+1}1/k^{p}dk=\frac{\eta_{1}}{1-p}\left((T+1)^{1-p}-1\right)$
(65a) $\displaystyle\sum_{s=k+1}^{T}\eta_{s}$
$\displaystyle\geq\eta_{1}\int_{s=k+1}^{T+1}1/s^{p}ds=\frac{\eta_{1}}{1-p}\left((T+1)^{1-p}-(k+1)^{1-p}\right)$
(65b)

Then we turn to estimate the last term of (E) as below

$\displaystyle\sum_{k=1}^{T}\left(\theta_{2}+\eta_{k}\sigma^{2}\right)\eta_{k}\cdot\exp\left(-\theta_{1}\sum_{s=k+1}^{T}\eta_{s}\right)$
$\displaystyle\leq\sum_{k=1}^{T}\left(\theta_{2}+\frac{\eta_{1}}{k^{p}}\sigma^{2}\right)\frac{\eta_{1}}{k^{p}}\cdot\exp\left(-\frac{\theta_{1}\eta_{1}}{1-p}\left(T^{1-p}-(k+1)^{1-p}\right)\right)$
$\displaystyle\leq\frac{\left(\theta_{2}+\eta_{1}\sigma^{2}\right)\eta_{1}}{\exp\left(\frac{\theta_{1}\eta_{1}}{1-p}T^{1-p}\right)}\sum_{k=1}^{T}\frac{1}{k^{p}}\exp\left(\frac{\theta_{1}\eta_{1}}{1-p}(k+1)^{1-p}\right)$
$\displaystyle\leq\frac{\left(\theta_{2}+\eta_{1}\sigma^{2}\right)\eta_{1}}{\exp\left(\frac{\theta_{1}\eta_{1}}{1-p}T^{1-p}\right)}\left(1+\frac{1}{k}\right)^{p}\int_{k=1}^{T+1}\frac{1}{(k+1)^{p}}\exp\left(\frac{\theta_{1}\eta_{1}}{1-p}(k+1)^{1-p}\right)dk$
$\displaystyle=\frac{\left(\theta_{2}+\eta_{1}\sigma^{2}\right)\eta_{1}}{\exp\left(\frac{\theta_{1}\eta_{1}}{1-p}T^{1-p}\right)}\frac{2^{p}}{\theta_{1}\eta_{1}}\left(\exp\left(\frac{\theta_{1}\eta_{1}}{1-p}(T+1)^{1-p}\right)-\exp\left(\frac{\theta_{1}\eta_{1}}{1-p}2^{1-p}\right)\right)$
$\displaystyle\leq\frac{\left(\theta_{2}+\eta_{1}\sigma^{2}\right)2^{p}}{\theta_{1}}.$

If $p=1$, then

$\displaystyle\sum_{k=1}^{T}\eta_{k}$
$\displaystyle\geq\eta_{1}\int_{k=1}^{T+1}\frac{dk}{k}=\eta_{1}\log(T+1)$
(66a) $\displaystyle\sum_{s=k+1}^{T}\eta_{s}$
$\displaystyle\geq\eta_{1}\int_{s=k+1}^{T+1}\frac{ds}{s}=\eta_{1}\log\left(\frac{T+1}{k+1}\right).$
(66b)

The last term of (E) can be bounded as

$\displaystyle\sum_{k=1}^{T}\left(\theta_{2}+\eta_{k}\sigma^{2}\right)\eta_{k}\cdot\exp\left(-\theta_{1}\sum_{s=k+1}^{T}\eta_{s}\right)$
$\displaystyle\leq\sum_{k=1}^{T}\left(\theta_{2}+\eta_{1}\sigma^{2}\right)\frac{\eta_{1}}{k}\cdot\exp\left(-\theta_{1}\eta_{1}\log\left(\frac{T+1}{k+1}\right)\right)$
$\displaystyle\leq\left(\theta_{2}+\eta_{1}\sigma^{2}\right)\eta_{1}\sum_{k=1}^{T}\frac{(k+1)^{\theta_{1}\eta_{1}}}{k(T+1)^{\theta_{1}\eta_{1}}}$
$\displaystyle\leq\left(\theta_{2}+\eta_{1}\sigma^{2}\right)\eta_{1}\left(1+\frac{1}{k}\right)\int_{k=1}^{T}\frac{(k+1)^{\theta_{1}\eta_{1}-1}}{(T+1)^{\theta_{1}\eta_{1}}}$
$\displaystyle\leq
2\left(\theta_{2}+\eta_{1}\sigma^{2}\right)\eta_{1}\frac{\left((T+1)^{\theta_{1}\eta_{1}}-2^{\theta_{1}\eta_{1}}\right)}{\theta_{1}\eta_{1}(T+1)^{\theta_{1}\eta_{1}}}$
$\displaystyle\leq\frac{2\left(\theta_{2}+\eta_{1}\sigma^{2}\right)}{\theta_{1}}.$

If $m\delta(k)\leq\eta_{k}\leq M\delta(k)$ and $\delta(k)=1/\alpha^{k-1}$
where $\alpha>1$. Then

$\displaystyle\sum_{k=1}^{N}\sum_{i=1}^{S}\eta_{i}^{k}$ $\displaystyle\geq
mS\sum_{k=1}^{N}\alpha^{-k+1}=\frac{mS(1-\alpha^{-N})}{1-\alpha^{-1}}>mS.$

The second term of (E) can be estimated as:

$\displaystyle\sum_{k=1}^{N}\sum_{i=1}^{S}\left(\theta_{2}+\eta_{i}^{k}\sigma^{2}\right)\eta_{i}^{k}\cdot\exp\left(-\theta_{1}\left(\sum_{s=i+1}^{S}\eta_{s}^{k}+\sum_{s=1}^{S}\sum_{l=k+1}^{N}\eta_{s}^{l}\right)\right)$
$\displaystyle\leq\sum_{k=1}^{N}\left(\theta_{2}+M\sigma^{2}\right)M\alpha^{-k+1}\exp\left(-m\theta_{1}\left(S\sum_{l=k+1}^{N}\alpha^{-l+1}\right)\right)\left(1+\sum_{i=1}^{S-1}\exp\left(-m\theta_{1}(S-i)\alpha^{-k+1}\right)\right)$
$\displaystyle\leq\left(\theta_{2}+M\sigma^{2}\right)M\sum_{k=1}^{N}\alpha^{-k+1}\exp\left(-m\theta_{1}\left(S\sum_{l=k+1}^{N}\alpha^{-l+1}\right)\right)\left(1+\frac{\exp(-m\theta_{1}\alpha^{-k+1})}{1-\exp(-m\theta_{1}\alpha^{-k+1})}\right)$
$\displaystyle\mathop{\leq}^{(a)}\frac{\left(\theta_{2}+M\sigma^{2}\right)M}{1-\exp(-m\theta_{1})}\left(1+\sum_{k=1}^{N-1}\exp\left(-m\theta_{1}\left(S\sum_{l=k+1}^{N}\alpha^{-l+1}\right)\right)\right)$
$\displaystyle\leq\frac{\left(\theta_{2}+M\sigma^{2}\right)M}{1-\exp(-m\theta_{1})}\left(1+\sum_{k=1}^{N-1}\exp\left(-m\theta_{1}S\alpha^{-N+1}\right)^{k}\right)\leq\frac{\left(\theta_{2}+M\sigma^{2}\right)M}{1-\exp(-m\theta_{1})}\left(1+\frac{\exp\left(-m\theta_{1}S\alpha^{-N+1}\right)}{1-\exp\left(-m\theta_{1}S\alpha^{-N+1}\right)}\right)$
$\displaystyle\leq\frac{\left(\theta_{2}+M\sigma^{2}\right)M}{1-\exp(-m\theta_{1})}\frac{1}{1-\exp\left(-m\theta_{1}S\alpha^{-N+1}\right)}$

where $(a)$ follows from the fact that for all $k\in[N]$, we have

$\displaystyle\frac{\alpha^{-k+1}}{1-\exp(-m\theta_{1}\alpha^{-k+1})}\leq\frac{1}{1-\exp(-m\theta_{1})}$
(67)

due to that $h(x)=\frac{x}{1-\exp(-m\theta_{1}x)}$ is increasing with $x$. By
properly choosing $S$ and $N$, for example: $N=\log_{\alpha}T/2$ and
$S=2T/\log_{\alpha}T$, we have
$S\alpha^{-N+1}=\frac{2T\alpha}{\sqrt{T}\log_{\alpha}T}=\frac{2\alpha\sqrt{T}}{\log_{\alpha}T}$.
For sufficient large $T$, we can see that
$0<\exp\left(-m\theta_{1}S\alpha^{-N+1}\right)\ll 1$. Thus,
$1-\exp\left(-m\theta_{1}S\alpha^{-N+1}\right)\approx 1$, we can achieve that

$\displaystyle\mathbb{E}[\left\|x_{T+1}-x^{\ast}\right\|^{2}]\leq\exp\left(-\theta_{1}mS\right)\left\|x_{1}-x^{\ast}\right\|^{2}+\frac{\left(\theta_{2}+M\sigma^{2}\right)M}{1-\exp(-m\theta_{1})}.$

∎