# Three-Tone Coherent Microwave Electromechanical Measurement of a Superfluid
Helmholtz Resonator

S. Spence<EMAIL_ADDRESS>Department of Physics, University of Alberta,
Edmonton, Alberta T6G 2E9, Canada  E. Varga  Faculty of Mathematics and
Physics, Charles University, Ke Karlovu 3, 121 16 Prague, Czech Republic  C.
A. Potts  Kavli Institute of NanoScience, Delft University of Technology, PO
Box 5046, 2600 GA Delft, Netherlands  J. P. Davis<EMAIL_ADDRESS>Department of Physics, University of Alberta, Edmonton, Alberta T6G 2E9,
Canada

###### Abstract

We demonstrate electromechanical coupling between a superfluid mechanical mode
and a microwave mode formed by a patterned microfluidic chip and a 3D cavity.
The electric field of the chip-cavity microwave resonator can be used to both
drive and detect the motion of a pure superflow Helmholtz mode, which is
dictated by geometric confinement. The coupling is characterized using a
coherent measurement technique developed for measuring weak couplings deep in
the sideband unresolved regime. The technique is based on two-probe
optomechanically induced transparency/amplification using amplitude
modulation. Instead of measuring two probe tones separately, they are
interfered to retain only a signal coherent with the mechanical motion. With
this method, we measure a vacuum electromechanical coupling strength of
$g_{0}=2\pi\times 23.3\;\mathrm{\mu}$Hz, three orders of magnitude larger than
previous superfluid electromechanical experiments.

††preprint: a

Cavity optomechanics, Aspelmeyer, Kippenberg, and Marquardt (2014) the
coupling of optical and mechanical resonances, is an extremely sensitive
scheme for detecting and controlling mechanical motion. Popular optomechanical
architectures include vibrating membranes, Yuan  _et al._ (2015); Noguchi  _et
al._ (2016); Pearson  _et al._ (2020); Pate  _et al._ (2020) superconducting
drums,Teufel  _et al._ (2011); Kotler  _et al._ (2021); Cattiaux  _et al._
(2021) phononic crystal cavities,Arrangoiz-Arriola  _et al._ (2018); Serra
_et al._ (2021) and magnetic spheres.Zhang  _et al._ (2016); Potts  _et al._
(2021) While decidedly more exotic, superfluid helium has been explored as a
promising mechanical element for cavity optomechanical and electromechanical
experiments – in part due to its large bandgap,Glaberson and W. Johnson
(1975); Kandula  _et al._ (2010) low dielectric loss,Hartung  _et al._ (2006)
and ultra-low acoustic loss at millikelvin temperatures.De Lorenzo and Schwab
(2014); Souris  _et al._ (2017) Superfluid helium mechanical resonators have
been used to develop gravitational wave detectors,Singh  _et al._ (2017);
Vadakkumbatt  _et al._ (2021) and have even been suggested for generating a
mechanical qubit Sfendla  _et al._ (2021) and for detecting dark matter,Manley
_et al._ (2020) while superfluid optomechanical resonators are allowing new
studies of superfluid properties, such as novel explorations of vortex
dynamics.Sachkou  _et al._ (2019) Microfabricated superfluid mechanical
resonators — to date lacking the optomechanical component — have also proven
to be particularly useful in the study of confined helium, where confinement
can dramatically change the physics of the system.Levitin  _et al._ (2013);
Perron, Kimball, and Gasparini (2019) Previously, our group has developed
microfluidic Helmholtz resonators for this purpose,Duh  _et al._ (2012); Rojas
and Davis (2015); Souris  _et al._ (2017) which have been used in 4He to
discover bi-stable turbulence in 2D superflow Varga  _et al._ (2020) and
reveal surface-dominated finite-size effects at the nanoscale.Varga,
Undershute, and Davis (2022) These devices have also been used to search for
pair density wave states in superfluid 3He.Shook  _et al._ (2020) Here, we
combine the concepts of microwave cavity electromechanics with microfluidic
confinement of superfluid helium-4 within a microfabricated Helmholtz
resonator.

Figure 1:  The chip-cavity system. The outline shows the 3D microwave cavity
and pin coupler, while the relative electric field magnitude of the
fundamental chip-cavity microwave mode, concentrated between the central basin
capacitive plates and around the antenna, is shown according to the left color
bar. Insert: the relative pressure of the superfluid Helmholtz mode within the
fluid geometry of the chip. Pressure is concentrated in the central basin,
with relative amplitude according to the right color bar.

In the present work, a helium volume is enclosed between two nanofabricated
quartz substrates, defining a fourth soundTilley and Tilley (2019) (pure
superflow) Helmholtz mode inside the microfluidic chip.Duh  _et al._ (2012);
Rojas and Davis (2015); Souris  _et al._ (2017) The pressure field of this
mode is shown in the insert of Fig. 1. Previously, Varga et al.Varga and Davis
(2021) showed the value in electromechanical drive and detection of the
Helmholtz mode using kHz-frequency carrier signals in a cavity-less system.
Here, we bring the readout into the GHz microwave regime by integrating the
chip into a hermetically-sealed 3D microwave cavity filled with superfluid
4He, with on-chip superconducting aluminum electrodes concentrating the
cavity’s electric field into the Helmholtz basin (main image Fig. 1), similar
to work with SiN membrane chips in 3D microwave cavities.Yuan  _et al._
(2015); Noguchi  _et al._ (2016) A three-tone coherent drive and measurement
technique was developed to observe small mechanical signals and calibrate the
optomechanical coupling rate, $g_{0}$, using a strong pump microwave tone
which is amplitude modulated to give two weak probe tone sidebands. The chip-
cavity system’s effect on the two sidebands is analogous to optomechanically
induced transparency and amplification (OMIT/A).Weis  _et al._ (2010)
Demodulating the signal to destructively interfere the two probe tones
recovers only the coherent microwave signal containing information about the
mechanical motion. This provides a powerful measurement technique in the
sideband unresolved regime, especially useful for detecting mechanics in
systems with weak opto/electro-mechanical couplings. We use this technique to
characterize the chip-cavity electromechanical coupling, measuring a vacuum
electromechanical coupling of $g_{0}=2\pi\times 23.3\;\mathrm{\mu}$Hz.

Figure 2:  Left: an RLC representation of the chip-cavity system, showing the
3D microwave cavity in blue, the microfluidic chip in black, and the
superfluid Helmholtz mode in red. Right: A typical reflection measurement (via
a directional coupler) of the cavity at 600 mK when filled with superfluid
4He. The magnitude is shown in grey, while the detuned phase is shown in blue.
Frequency has been detuned around the microwave cavity frequency
$\omega_{c}=2\pi\times 2.28$ GHz.

As mentioned above, the experimental system consists of a superfluid Helmholtz
resonator, similar to Refs. Souris  _et al._ , 2017; Varga  _et al._ , 2020;
Varga and Davis, 2021, but now incorporated into a hermetic 3D microwave
cavity. The etched Helmholz geometry is a wide flat circular basin, 7 mm in
diameter and 1.01 $\mathrm{\mu}$m deep, connected by two 1.6 mm wide 1.6 mm
long channels to the surrounding helium bath. Electrodes of 50 nm thick
aluminum are deposited into the etched regions: a 6 mm diameter circle at the
center of the basin connects to a meandering 100 $\mathrm{\mu}$m-wide lead (in
one channel) that terminates in a 5 mm $\times$ 2 mm square antenna. After
fabrication, the wafer is diced into identical top and bottom chips, which are
room-temperature direct wafer bonded Tong  _et al._ (1994); Plach  _et al._
(2013) such that the basins align and the antenna are on opposite sides of the
basin (as shown in Fig. 1). The combined basins form a single volume with
aluminum top and bottom electrodes, turning the combined basin into a parallel
plate capacitor, with an antenna connected to each plate – providing
capacitive coupling to the 3D cavity mode (more detailed fabrication notes can
be found in Ref. Souris  _et al._ , 2017).

The Helmholtz mode is an acoustic resonance of the superfluid moving back and
forth in the channel, driven by pressure fluctuations in the central basin
(shown in the insert of Fig. 1). The mode can be described as a mass-spring
system,Varga  _et al._ (2020) with effective mass
$m_{\text{eff}}=2wlD\rho_{\text{s}}$; where $w$, $l$ and $D$ are the channel
width, length, and depth respectively; and $\rho_{\text{s}}$ is the superfluid
density. The resulting mode frequency is then:Varga and Davis (2021)

$\Omega_{\text{m}}=\sqrt{\frac{k_{\text{eff}}}{m_{\text{eff}}}}=\sqrt{\frac{2wD}{l\rho}\frac{\rho_{\text{s}}}{\rho}\frac{k_{\text{p}}}{4A^{2}(1+\Sigma)}},$
(1)

with $k_{\text{eff}}$ the effective spring constant, $k_{\text{p}}$ the
effective stiffness of the mean deflection of the basin walls, $A=\pi R^{2}$
the area of the basin, and $\Sigma=\chi_{\text{B}}Dk_{\text{p}}/(4A)$, where
$\chi_{\text{B}}$ is the bulk compressibility of helium. Here, the effective
spring constant of the mode is significantly softened by the flexing of the
substrate when compared to the compressibility-only case. The flexing of the
substrate changes the distance between the capacitive plates above and below
the basin, the origin of the electromechanical coupling for this system, which
is slightly suppressed by electrostriction due to the compression of the
helium.De Lorenzo (2016); Spence (2022) Considering both of these factors, the
change in basin capacitance $C_{\text{bas}}$ can be written in terms of the
displacement of the fluid in the channels $y$, as:Varga and Davis (2021)

$\frac{\mathrm{d}C_{\text{bas}}}{\mathrm{d}y}=\bar{C}_{\text{bas}}\frac{2w\rho_{\text{s}}}{A\rho}\left(1-\frac{\varepsilon-1}{\varepsilon}\Sigma\right)\frac{1}{1+\Sigma},$
(2)

where $\bar{C}_{\text{bas}}$ is the undeformed basin capacitance and
$\varepsilon$ is the dielectric constant of liquid helium.

The Helmholtz chip is placed at the center of a rectangular 3D microwave
cavity,Pozar (2011) and the 3D cavity mode is capacitively coupled to the
basin capacitor via the on-chip antenna,Schuster (2007) effectively
concentrating the electric field of the fundamental chip-cavity mode into the
basin, as shown in Fig. 1. The lumped element RLC representation of the chip-
cavity system is shown in Fig. 2 (a), and (b) shows a typical $S_{21}$ single
port bi-directional (reflection) measurement of the microwave resonance at 600
mK with the cavity filled with superfluid helium-4 close to saturated vapor
pressure (SVP). Microwave circuit schematics can be found in the Supplementary
Information (SI). From this measurement, the resonant frequency is found to be
$\omega_{\text{c}}=2\pi\times 2.28$ GHz, with total cavity loss rate
$\kappa=2\pi\times 1.76$ MHz, and external coupling rate
$\kappa_{\text{in}}=2\pi\times 128$ kHz, meaning the system is significantly
undercoupled. Rieger  _et al._ (2022)

The vacuum electromechanical coupling strength $g_{0}=Gx_{\text{zpf}}$, with
$G=-\partial\omega_{\text{c}}/\partial x$, can be considered as the shift in
cavity resonance due to the zero point motion of the mechanical oscillator.
For the Helmholtz chip-cavity system, this can be written as:

$g_{0}=-\delta y_{\text{zpf}}\frac{\partial\omega_{\text{c}}}{\partial
y}=-\delta y_{\text{zpf}}\frac{\partial\omega_{\text{c}}}{\partial
C}\frac{\partial C}{\partial y},$ (3)

where $y_{\text{zpf}}$ is the zero point motion of the superfluid in the
channels, and $C$ is the total capacitance of the chip-cavity microwave mode
given by $\omega_{\text{c}}=1/\sqrt{LC}$. Using $\delta
y_{\text{zpf}}=\sqrt{\hbar/2m_{\text{eff}}\Omega_{\text{m}}}$, Eq. 2, and that
the temperature is cold enough for $\rho_{\text{s}}\approx\rho$, $g_{0}$ can
be expressed for this electromechanical coupling scheme as:

$g_{0}\simeq\alpha\omega_{\text{c}}\sqrt{\frac{\hbar}{2m_{\text{eff}}\Omega_{\text{m}}}}\frac{w}{A}\left(1-\frac{\varepsilon-1}{\varepsilon}\Sigma\right)\frac{1}{1+\Sigma}.$
(4)

Here, $\alpha$ is the proportion of the microwave mode’s total electric field
energy within the basin capacitor; for this work, $\alpha\approx 0.005$ is
calculated via finite element method (FEM) simulations. Equation 4 gives a
method for estimating the $g_{0}$ of a fabricated device. Only
$\Omega_{\text{m}}$ and $k_{\text{p}}$ are unknown, yet $\Omega_{\text{m}}$
can be estimated using Eq. 1 while $k_{\text{p}}=9.45\times 10^{6}$ N m-1 is
measured via capacitance change with a DC bias, as in Ref. Souris  _et al._ ,
2017. Using these values an estimate of $g_{0}\approx 2\pi\times 1.1\;\mu$Hz
is obtained; although lower than the value found later in this work, this
discrepancy is likely due to inaccuracy in $\alpha$.

We began by observing electromechanical coupling using a calibrated homodyne
measurement:Gorodetksy  _et al._ (2010); Kumar  _et al._ (2023); Spence (2022)
sending a tone at the microwave resonance’s maximum gradient, demodulating the
returned signal with a tone of the same frequency, measuring noise in the time
domain, and considering the power spectral density of fluctuations near the
mechanical frequency $\Omega_{\text{m}}$ (circuit shown in the Supplementary
Information). Reflection via a directional coupler was measured as the return
signal. The mechanical spectrum was found using an amplitude-sensitive scheme,
giving a Helmholtz mode frequency of $\Omega_{\text{m}}=2\pi\times 1394.7$.
The homodyne phase-lock was achieved by applying a 10 kHz phase modulation to
the signal on the cavity arm of the circuit and then minimizing the 10 kHz
signal after down-mixing by varying the LO phase. Phase-sensitive measurement
was also possible by instead sending the microwave tone on resonance and using
10 kHz amplitude modulation to phase-lock.

Despite homodyne measurement enabling detection of the mechanical spectrum,
the signal is weak, with a small signal-to-noise ratio near the mechanical
resonance frequency, mainly due to low electromechanical coupling strength and
undercoupled microwave readout, even after more than 450 averages. Moreover,
determination of $g_{0}$ was not possible due to uncertainty in the effective
temperature of the Helmholtz mode since vibrations from the cryocooler drove
the superfluid motion in the present experiment.

Figure 3:  Diagram of the three-tone coherent measurement scheme. a) Frequency
space diagram of the three tones relative to the microwave resonance
$\omega_{\text{c}}$. Showing the pump tone $\omega_{\text{d}}$ at detuning
$\Delta$ from the cavity resonance, along with the two probe tones
$\omega^{\pm}_{\text{p}}$ created via amplitude modulation of the pump tone
with frequency (separation) $\Omega_{\text{p}}$. The cavity resonance
lineshape has been exaggerated for clarity. b) Diagram of the two probe tones
before (left - blue) and after (right - red) demodulation. The left side shows
the response at the probe tones $\omega^{\pm}_{\text{p}}$ as the modulation
frequency $\Omega_{\text{p}}$ is swept across the mechanical frequency
$\Omega_{\text{m}}$, showing the EMIT/A response for the upper and lower
sidebands respectively. The flat background is due to being deep into the
sideband unresolved regime $\kappa\gg\Omega_{\text{m}}$. The right side shows
the measured signal at relative phase $\phi=0$, after demodulating the two
probe tones with the pump frequency $\omega_{\text{d}}$. This down-mixes the
probe signals to frequency $\pm\Omega_{\text{p}}$, where the two probe tones
destructively interfere, leaving only the combined EMIT/A signal with a peak
at $\Omega_{\text{m}}$.

A coherent measurement scheme can be used to precisely determine the
electromechanical coupling rate, $g_{0}$. The scheme in this work is based on
OMIT/A:Weis  _et al._ (2010); Agarwal and Huang (2010); Safavi-Naeini  _et
al._ (2011); Kashkanova  _et al._ (2017a) an interference effect where beating
between two optical tones causes coherent oscillation of the mechanics via
radiation pressure when the two tones are separated in frequency by the
mechanical resonance. The oscillation of the mechanics, in turn, scatters
photons from the stronger driving ‘pump’ tone to the weak ‘probe’ tone
frequency, where they then destructively (constructively) interfere with the
probe tone in the case of OMIT (OMIA). In the microwave regime, the equivalent
effect is often called electromechanically induced transparency/amplification
(EMIT/A).Zhou  _et al._ (2013) For the Helmholtz devices, mechanical motion is
driven via modulation of the electrostatic force between the two capacitive
basin plates, proportional to voltage squared. Typically OMIT/A is used for
sideband resolved systems ($\kappa\ll\Omega_{\text{m}}$), with stronger
optomechanical coupling strengths compared to this work. To overcome these
limitations, a ‘two-probe’ measurement scheme is used here. The two probe
tones are produced by amplitude modulation of the pump tone, inspired by the
work of Kashkanova et al.,Kashkanova  _et al._ (2017b) and adapted from optics
to the microwave regime. For this work, modulation of a microwave signal
allows straightforward resolution of a mHz mechanical linewidth on a GHz
signal and, for kHz mechanics deep into the sideband unresolved regime,
enables a unique measurement scheme where the probe tones are canceled to
recover only a coherent electromechanical signal, which would not be visible
otherwise.

Figure 3 shows a diagram of the coherent measurement scheme in this work. In
Fig. 3(a), a strong microwave pump tone with frequency $\omega_{\text{d}}$ is
sent into the cavity, at some detuning $\Delta$ from the cavity resonance
$\omega_{\text{c}}$, such that $\Delta=\omega_{\text{d}}-\omega_{\text{c}}$.
This pump tone is amplitude modulated with frequency $\Omega_{\text{p}}$ to
produce the two weak probe signals at
$\omega^{\pm}_{\text{p}}=\omega_{\text{d}}\pm\Omega_{\text{p}}$. The signal
incident at the cavity port can then be written as
$s_{\text{in}}=\bar{s}_{\text{in}}+\delta s_{\text{in}}$, where the weak
probes provide $\delta s_{\text{in}}$, which in a frame rotating at
$\omega_{\text{d}}$ can be written:

$\delta
s_{\text{in}}=s_{\text{p}}\left(\mathrm{e}^{-i\Omega_{\text{p}}t}+\mathrm{e}^{+i\Omega_{\text{p}}t}\right)/2.$
(5)

Here, $s_{\text{p}}/2$ is the amplitude of the probe signals. Sweeping
$\Omega_{\text{p}}$ across the mechanical frequency $\Omega_{\text{m}}$,
equivalent to sweeping the probe and pump separation, causes the EMIT/A effect
close to $\Omega_{\text{m}}$ for the upper and lower sidebands respectively,
shown in the insert of Fig. 3(a). Following the full derivation in the
Supplementary Information, the intracavity field amplitude at the two probe
frequencies $\omega^{\pm}_{\text{p}}$ can be written as

$A^{\pm}(\Omega_{\text{p}})=\frac{\chi_{\mbox{\scriptsize$\pm$}}\sqrt{\kappa_{\text{in}}}s_{\text{p}}}{2}\left(1\mp\frac{g(g\chi_{\mp}^{*}+g^{*}\chi_{\mbox{\scriptsize$\pm$}})}{\Gamma_{\text{m}}/2\mp
i\delta\pm i\Sigma(\pm\Omega_{\text{p}})}\right),$ (6)

where $\Gamma_{\text{m}}$ is the mechanical linewidth,
$\delta=\Omega_{\text{p}}-\Omega_{\text{m}}$ is the detuning of the pump-probe
separation from the mechanical frequency, $g$ is the multiphoton coupling
defined as

$g=g_{0}\frac{\sqrt{\kappa_{\text{in}}}\bar{s}_{\text{in}}}{-i\Delta+\kappa/2},$
(7)

$\chi_{\mbox{\scriptsize$\pm$}}$ is the cavity susceptibility at
$\omega^{\pm}_{\text{p}}$

$\chi_{\mbox{\scriptsize$\pm$}}=\frac{1}{-i(\Delta\pm\Omega_{\text{p}})+\kappa/2},$
(8)

and the electromechanical self-energy $\Sigma(\Omega_{\text{p}})$ is defined
as

$i\Sigma(\Omega_{\text{p}})=|g|^{2}(\chi_{+}-\chi_{-}^{*}).$ (9)

When measuring in a bi-directional single port scheme (reflection via a
circulator / directional coupler) $-\sqrt{\kappa_{\text{in}}}A^{\pm}$ is
effectively the modification of the reflected probe tones due to the cavity
electromechanical system. Demodulating this reflected signal at
$\omega_{\text{d}}$ and a relative phase of $\phi=0$, will give a signal at
$\Omega_{\text{p}}$ with amplitude proportional to the upper probe minus the
lower probe, canceling out the probe tone backgrounds and recovering only the
coherent effect, shown in Fig 3(b) and (c). Measuring this signal at
$\Omega_{\text{p}}$ using a lock-in amplifier will give a signal

$X+iY=\frac{-iK_{\text{T}}g_{0}^{2}\kappa_{\text{in}}^{2}s_{\text{p}}}{\sqrt{2}(\Gamma_{\text{m}}/2+i\delta)}\frac{(\kappa^{2}/4-\Delta^{2})}{(\kappa^{2}/4+\Delta^{2})^{3}}\frac{P_{\text{in}}}{\hbar\omega_{\text{d}}},$
(10)

where $K_{\text{T}}$ is the transfer function of the circuit post reflection
from the cavity and $P_{\text{in}}$ is the total power incident on the cavity
port. A full derivation of $X+iY$ can be found in the Supplementary
Information. Figure 4(b) shows a measurement of the $Y$ quadrature, fit to Eq.
10, showing an improvement of fitting accuracy and signal-to-noise near
resonance, compared to the homodyne measurement of noise amplitude spectral
density in Fig. 4 (a). From the coherent measurement fit, the mechanical
frequency and linewidth are found to be $\Omega_{\text{m}}=2\pi\times 1394.7$
Hz and $\Gamma_{\text{m}}=2\pi\times 0.106$ Hz. Although it may instead be
possible to cancel the probe tones using an additional source, with the
current method far into the sideband unresolved regime
($\kappa\gg\Omega_{\text{m}}$), the two probe tones effectively travel the
same path, canceling near perfectly and accounting for any frequency, phase,
or amplitude noise. Using the three-tone method, the background signal at
$\Omega_{\text{m}}$ was reduced by three to four orders of magnitude, allowing
a measurement of the mechanics that would not normally be possible with
EMIT/A.

Figure 4:  a) Homodyne noise measurement: plotting noise amplitude spectral
density of the Helmholtz mode, as driven by the environmental noise of the
cryocooler. Shown is the average of 453 acquisitions. b) Coherent three-tone
measurement (EMIT/A): plotting the amplitude measured at the lock-in vs.
amplitude modulation frequency, detuned around the mechanical resonance
$\Omega_{\text{m}}$. Shown is the average of 37 acquisitions.

The circuit transfer function $K_{\text{T}}$ in Eq. 10 can be calibrated to
allow determination of $g_{0}$ using a measurement of $R=|X+iY|$ with
$\phi=\pi/2$, to give a signal proportional to the combined probe tone
amplitudes $R_{\text{cal}}=K_{\text{T}}s_{\text{p}}/2\sqrt{2}$ (further
details in the Supplementary Information). Using this calibration, $g_{0}$ can
be written as

$g_{0}=\frac{\sqrt{\hbar\omega_{\text{d}}}}{2\kappa_{\text{in}}\sqrt{P_{\text{in}}}}\sqrt{\frac{\Gamma_{m}(\kappa^{2}/4+\Delta^{2})^{3}}{(\kappa^{2}/4-\Delta^{2})}\frac{R_{\delta=0}}{R_{cal}}}.$
(11)

Here, $R_{\delta=0}$ is the measured signal amplitude at
$\Omega_{\text{p}}=\Omega_{\text{m}}$ when the probe tones cancel, either
directly from the peak signal or from fitting to Eq. 10. Using this
methodology, the vacuum electromechanical coupling strength was calculated as
$g_{0}=2\pi\times 23.3^{+6.0}_{-4.8}\;\mu$Hz, with determination of
$P_{\text{in}}$ providing the dominant error. This $g_{0}$ is in the weak
coupling limit of electromechanics,Aspelmeyer, Kippenberg, and Marquardt
(2014) but a three order of magnitude increase when compared to a bulk
superfluid-microwave coupling.De Lorenzo and Schwab (2017) While the
corresponding single-photon cooperativity $C_{0}=1.2\times 10^{-14}$ is one of
the lowest recorded,Aspelmeyer, Kippenberg, and Marquardt (2014) and will need
to be improved in subsequent experiments, this demonstrates the sensitivity of
the three-tone coherent measurement scheme.

In conclusion, the first superfluid microwave electromechanical coupling was
measured for a micromechanical system, using a microfluidic Helmholtz
resonator inside a 3D microwave cavity. The mechanical motion was detected
using a homodyne noise measurement scheme; however, thermomechanical
calibration was not possible due to strong cryocooler vibrations. To measure
an electromechanical system with weak coupling, deep in the sideband
unresolved regime, a coherent measurement scheme was developed, inspired by
OMITWeis  _et al._ (2010) and specifically the work of Kashkanova et
al.Kashkanova  _et al._ (2017b) The scheme uses amplitude modulation of a
strong pump tone to create probe sidebands, which experience coherent EMIT/A.
The two probe tones are then canceled to recover only signal coherent with the
mechanics, allowing an electromechanical calibration of the coupling strength
$g_{0}$.

The coherent measurement scheme provides a powerful tool for measuring
electro/optomechanical systems deep into the sideband unresolved regime.
Moreover, integrating a superfluid Helmholtz device into a microwave cavity
provides a new platform for superfluid microwave electromechanics. With three
orders of magnitude improvement to $g_{0}$ over previous superfluid
electromechanical systems,De Lorenzo (2016) the microwave Helmholtz design
shows promise as a sensitive method for measuring the properties of superfluid
helium. The sensitivity of this proof of concept can be improved by orders of
magnitude in many areas: confinement could decrease $\Gamma_{\text{m}}$ while
increasing $g_{0}$ (possibly with a phononic crystalSpence  _et al._ (2021);
Sfendla  _et al._ (2021)), larger antenna or galvanic bonding would increase
$\alpha$,Noguchi  _et al._ (2016) low mK dilution temperatures would decrease
superfluid and microwave dissipation, a superconducting microwave cavity could
also decrease dissipation, and a stronger coupling into the cavity
($\kappa_{\text{in}}\approx\kappa/2$) would greatly improve signal strength
for any given cavity field. The improved system may be useful for studying
non-equilibrium thermodynamics,Awschalom and Schwarz (1984) or even vector
dark matter.Manley  _et al._ (2021)

The authors acknowledge that the land on which this work was performed is in
Treaty Six Territory, the traditional territories of many First Nations,
Métis, and Inuit in Alberta. They acknowledge support from the University of
Alberta; the Natural Sciences and Engineering Research Council, Canada (Grant
Nos. RGPIN-2022-03078, and CREATE-495446-17); the NSERC Alberta Innovates
Advance program (ALLRP-568609-21); and the Government of Canada through the
NRC Quantum Sensors Program.

## Author Delclarations

### Conflict of Interest

The authors have no conflicts to disclose.

### Author Contributions

Sebastian Spence: Conceptualization (supporting), Formal Analysis (lead),
Methodology (equal), Software (lead), Visualization (lead), Writing (lead).
Emil Varga: Conceptualization (lead), Formal Analysis (supporting),
Methodology (equal), Supervision (supporting), Writing (supporting). Clinton
Potts: Formal Analysis (supporting), Software (supporting), Writing
(supporting). John Davis: Conceptualization (supporting), Funding Acquisition
(lead), Supervision (lead), Writing (supporting).

## Data Availability

The data that support the findings of this study are available from the
corresponding author upon reasonable request.

## References

  * Aspelmeyer, Kippenberg, and Marquardt (2014) M. Aspelmeyer, T. J. Kippenberg, and F. Marquardt, Rev. Mod. Phys. 86, 1391 (2014).
  * Yuan  _et al._ (2015) M. Yuan, V. Singh, Y. M. Blanter, and G. A. Steele, Nat. Commun. 6, 8491 (2015).
  * Noguchi  _et al._ (2016) A. Noguchi, R. Yamazaki, M. Ataka, H. Fujita, Y. Tabuchi, T. Ishikawa, K. Usami, and Y. Nakamura, New J. Phys. 18, 103036 (2016).
  * Pearson  _et al._ (2020) A. Pearson, K. Khosla, M. Mergenthaler, G. A. D. Briggs, E. Laird, and N. Ares, Sci. Rep. 10, 1–7 (2020).
  * Pate  _et al._ (2020) J. M. Pate, M. Goryachev, R. Y. Chiao, J. E. Sharping, and M. E. Tobar, Nat. Phys. 16, 1117–1122 (2020).
  * Teufel  _et al._ (2011) J. D. Teufel, T. Donner, D. Li, J. W. Harlow, M. Allman, K. Cicak, A. J. Sirois, J. D. Whittaker, K. W. Lehnert, and R. W. Simmonds, Nature 475, 359–363 (2011).
  * Kotler  _et al._ (2021) S. Kotler, G. A. Peterson, E. Shojaee, F. Lecocq, K. Cicak, A. Kwiatkowski, S. Geller, S. Glancy, E. Knill, R. W. Simmonds, _et al._ , Science 372, 622–625 (2021).
  * Cattiaux  _et al._ (2021) D. Cattiaux, I. Golokolenov, S. Kumar, M. Sillanpää, L. Mercier de Lépinay, R. Gazizulin, X. Zhou, A. Armour, O. Bourgeois, A. Fefferman, _et al._ , Nat. commun. 12, 6182 (2021).
  * Arrangoiz-Arriola  _et al._ (2018) P. Arrangoiz-Arriola, E. A. Wollack, M. Pechal, J. D. Witmer, J. T. Hill, and A. H. Safavi-Naeini, Phys. Rev. X 8, 031007 (2018).
  * Serra  _et al._ (2021) E. Serra, A. Borrielli, F. Marin, F. Marino, N. Malossi, B. Morana, P. Piergentili, G. A. Prodi, P. M. Sarro, P. Vezio, _et al._ , J. Appl. Phys. 130, 064503 (2021).
  * Zhang  _et al._ (2016) X. Zhang, C.-L. Zou, L. Jiang, and H. X. Tang, Sci. Adv. 2, e1501286 (2016).
  * Potts  _et al._ (2021) C. A. Potts, E. Varga, V. A. Bittencourt, S. V. Kusminskiy, and J. P. Davis, Phys. Rev. X 11, 031053 (2021).
  * Glaberson and W. Johnson (1975) W. I. Glaberson and W. W. Johnson, J. Low Temp. Phys. 20, 313–338 (1975).
  * Kandula  _et al._ (2010) D. Z. Kandula, C. Gohle, T. J. Pinkert, W. Ubachs, and K. S. Eikema, Phys. Rev. Lett. 105, 063001 (2010).
  * Hartung  _et al._ (2006) W. Hartung, J. Bierwagen, S. Bricker, C. Compton, T. Grimm, M. Johnson, D. Meidlinger, D. Pendell, J. Popielarski, L. Saxton, _et al._ , Vacuum 10, 10 (2006).
  * De Lorenzo and Schwab (2014) L. De Lorenzo and K. Schwab, New J. Phys. 16, 113020 (2014).
  * Souris  _et al._ (2017) F. Souris, X. Rojas, P. H. Kim, and J. P. Davis, Phys. Rev. Appl. 7, 044008 (2017).
  * Singh  _et al._ (2017) S. Singh, L. De Lorenzo, I. Pikovski, and K. Schwab, New J. Phys. 19, 073023 (2017).
  * Vadakkumbatt  _et al._ (2021) V. Vadakkumbatt, M. Hirschel, J. Manley, T. Clark, S. Singh, and J. Davis, Phys. Rev. D 104, 082001 (2021).
  * Sfendla  _et al._ (2021) Y. L. Sfendla, C. G. Baker, G. I. Harris, L. Tian, R. A. Harrison, and W. P. Bowen, NPJ Quantum Inf. 7, 62 (2021).
  * Manley  _et al._ (2020) J. Manley, D. J. Wilson, R. Stump, D. Grin, and S. Singh, Phys. Rev. Lett. 124, 151301 (2020).
  * Sachkou  _et al._ (2019) Y. P. Sachkou, C. G. Baker, G. I. Harris, O. R. Stockdale, S. Forstner, M. T. Reeves, X. He, D. L. McAuslan, A. S. Bradley, M. J. Davis, _et al._ , Science 366, 1480–1485 (2019).
  * Levitin  _et al._ (2013) L. Levitin, R. Bennett, A. Casey, B. Cowan, J. Saunders, D. Drung, T. Schurig, and J. Parpia, Science 340, 841–844 (2013).
  * Perron, Kimball, and Gasparini (2019) J. Perron, M. Kimball, and F. Gasparini, Rep. Prog. Phys. 82, 114501 (2019).
  * Duh  _et al._ (2012) A. Duh, A. Suhel, B. Hauer, R. Saeedi, P. Kim, T. Biswas, and J. Davis, J. Low Temp. Phys. 168, 31–39 (2012).
  * Rojas and Davis (2015) X. Rojas and J. Davis, Phys. Rev. B 91, 024503 (2015).
  * Varga  _et al._ (2020) E. Varga, V. Vadakkumbatt, A. Shook, P. Kim, and J. Davis, Phys. Rev. Lett. 125, 025301 (2020).
  * Varga, Undershute, and Davis (2022) E. Varga, C. Undershute, and J. Davis, Phys. Rev. Lett. 129, 145301 (2022).
  * Shook  _et al._ (2020) A. Shook, V. Vadakkumbatt, P. S. Yapa, C. Doolin, R. Boyack, P. Kim, G. Popowich, F. Souris, H. Christani, J. Maciejko, _et al._ , Phys. Rev. Lett. 124, 015301 (2020).
  * Tilley and Tilley (2019) D. R. Tilley and J. Tilley,  _Superfluidity and Superconductivity_ (Routledge, 2019).
  * Varga and Davis (2021) E. Varga and J. P. Davis, New J. Phys. 23, 113041 (2021).
  * Weis  _et al._ (2010) S. Weis, R. Rivière, S. Deléglise, E. Gavartin, O. Arcizet, A. Schliesser, and T. J. Kippenberg, Science 330, 1520 (2010).
  * Tong  _et al._ (1994) Q.-Y. Tong, G. Cha, R. Gafiteanu, and U. Gosele, J. Microelectromech. Syst. 3, 29–35 (1994).
  * Plach  _et al._ (2013) T. Plach, K. Hingerl, S. Tollabimazraehno, G. Hesser, V. Dragoi, and M. Wimplinger, J. Appl. Phys. 113, 094905 (2013).
  * De Lorenzo (2016) L. A. De Lorenzo,  _Optomechanics with superfluid helium-4_ , Ph.D. thesis, California Institute of Technology (2016).
  * Spence (2022) S. Spence,  _Superfluid optomechanics with nanofluidic geometries_ , Ph.D. thesis, Royal Holloway, University of London (2022).
  * Pozar (2011) D. M. Pozar,  _Microwave engineering_ (John wiley & sons, 2011).
  * Schuster (2007) D. I. Schuster,  _Circuit quantum electrodynamics_ , Ph.D. thesis, Yale University (2007).
  * Rieger  _et al._ (2022) D. Rieger, S. Günzler, M. Spiecker, A. Nambisan, W. Wernsdorfer, and I. Pop, arXiv preprint arXiv:2209.03036 (2022).
  * Gorodetksy  _et al._ (2010) M. Gorodetksy, A. Schliesser, G. Anetsberger, S. Deleglise, and T. J. Kippenberg, Opt. Express 18, 23236 (2010).
  * Kumar  _et al._ (2023) S. Kumar, S. Spence, S. Perrett, Z. Tahir, A. Singh, C. Qi, S. Perez Vizan, and X. Rojas, J. Appl. Phys. 133, 094501 (2023).
  * Agarwal and Huang (2010) G. S. Agarwal and S. Huang, Phys. Rev. A 81, 041803 (2010).
  * Safavi-Naeini  _et al._ (2011) A. H. Safavi-Naeini, T. M. Alegre, J. Chan, M. Eichenfield, M. Winger, Q. Lin, J. T. Hill, D. E. Chang, and O. Painter, Nature 472, 69–73 (2011).
  * Kashkanova  _et al._ (2017a) A. Kashkanova, A. Shkarin, C. Brown, N. Flowers-Jacobs, L. Childress, S. Hoch, L. Hohmann, K. Ott, J. Reichel, and J. Harris, Nat. Phys. 13, 74–79 (2017a).
  * Zhou  _et al._ (2013) X. Zhou, F. Hocke, A. Schliesser, A. Marx, H. Huebl, R. Gross, and T. J. Kippenberg, Nat. Phys. 9, 179–184 (2013).
  * Kashkanova  _et al._ (2017b) A. Kashkanova, A. Shkarin, C. Brown, N. Flowers-Jacobs, L. Childress, S. Hoch, L. Hohmann, K. Ott, J. Reichel, and J. Harris, J. Opt. 19, 034001 (2017b).
  * De Lorenzo and Schwab (2017) L. De Lorenzo and K. Schwab, J. Low Temp. Phys. 186, 233–240 (2017).
  * Spence  _et al._ (2021) S. Spence, Z. Koong, S. Horsley, and X. Rojas, Phys. Rev. Appl. 15, 034090 (2021).
  * Awschalom and Schwarz (1984) D. Awschalom and K. Schwarz, Phys. Rev. Lett. 52, 49 (1984).
  * Manley  _et al._ (2021) J. Manley, M. D. Chowdhury, D. Grin, S. Singh, and D. J. Wilson, Phys. Rev. Lett. 126, 061301 (2021).
  * Buters  _et al._ (2017) F. Buters, F. Luna, M. Weaver, H. Eerkens, K. Heeck, S. De Man, and D. Bouwmeester, Opt. Express 25, 12935–12943 (2017).

## Supplementary Information

##  Appendix A Three-Tone Coherent Measurement.

This section will detail the theory of three-tone coherent measurement,
inspired by the work of Kashkanova et al. Kashkanova  _et al._ (2017b),
adapted to microwave circuits and developed to include probe tone canceling.
This technique enables sensitive measurement of weak electromechanical
couplings deep into the sideband-unresolved regime. First, the theory for two-
probe EMIT/A is outlined using an amplitude modulation scheme, then developed
to include destructive interference of the probe tones, recovering only the
signal coherent with the mechanical motion. Finally, details are provided for
calibrating coherent signal strength against the constructively interfered
probe tones. Figure 5 shows the microwave measurement circuit for the three-
tone coherent measurement.

Figure 5: Three-tone coherent measurement circuit. The ‘RF Source’ microwave
signal at pump frequency $\omega_{\rm d}$ is amplitude modulated with
frequency $\Omega_{\rm p}$ and the ‘10 kHz’ reference signal, which is used
for phase-locking. The RF signal is sent to the experimental ‘Cell’ (cavity-
Helmholtz) within the 3He cryostat (at 600 mK) via a directional coupler. The
reflected microwave signal passes back through the directional coupler, is
amplified by two room temperature microwave amplifiers, and is mixed with a
signal of frequency $\omega_{\rm d}$ from the ‘LO Source’. The mixed signal is
low-pass filtered to remove the $2\omega_{\rm d}$ components and then measured
on two lock-in amplifiers. The lock-in set to $\Omega_{\rm p}$ (lower on the
diagram) measures the signal coherent with the mechanics, while the 10 kHz
lock-in effectively measures the relative phase $\phi_{\rm ref}-\phi$ at the
mixer. The 10 kHz lock-in measurement is taken as the input of a PI loop,
which maintains the relative phase at the mixer by adjusting the LO Source
output phase.

###  A.1 Two-Probe EMIT/A

This subsection will outline the theory for two-probe EMIT/A using an
amplitude modulation scheme. The method in this section is based on Ref. Weis
_et al._ (2010), modified to include two probe beams provided by amplitude
modulation of the pump tone. Refs. Kashkanova  _et al._ (2017b) and Kashkanova
_et al._ (2017a) have also been used to assist with these modifications. To
begin, the electromechanical Hamiltonian can be written as

$\hat{\mathcal{H}}=\frac{\hat{p}^{2}}{2m_{\rm eff}}+\frac{1}{2}m_{\rm
eff}\Omega_{\rm m}^{2}\hat{x}^{2}+\hbar\omega_{\rm
c}\left(\hat{a}^{\dagger}\hat{a}+\frac{1}{2}\right)+\hbar
G\hat{a}^{\dagger}\hat{a}\hat{x}+i\hbar\sqrt{\kappa_{\rm in}}\left(s_{\rm
in}(t)\hat{a}^{\dagger}-s_{\rm in}^{*}(t)\hat{a}\right)\;.$ (12)

Here, $\hat{p}$ and $\hat{x}$ are the momentum and position operators of the
mechanical mode, with resonant frequency $\Omega_{\rm m}$ and effective mass
$m_{\rm eff}$; $\hat{a}$ and $\hat{a}^{\dagger}$ are the creation and
annihilation operators of the microwave cavity mode, with resonant frequency
$\omega_{\rm c}$, external loss rate $\kappa_{\rm in}$, and incident field at
the cavity port $s_{\rm in}(t)$; and finally $G=\partial\omega_{\rm
c}/\partial x$ is the electromechanical coupling strength.

For a strong microwave pump (drive) tone at $\omega_{\rm d}$ and weak probe
tones (one or more), the incident signal can be written $s_{\rm
in}(t)=\left(\bar{s}_{\rm in}+\delta s(t)\right)\mathrm{e}^{-i\omega_{\rm
d}t}$, where $\delta s_{\rm in}(t)$ contains the probe tones and the pump tone
amplitude $\bar{s}_{\rm in}$ is both positive and real. Then, writing the
Langevin equations in a frame rotating at $\Delta=\omega_{\rm d}-\omega_{\rm
c}$ gives

$\frac{\mathrm{d}}{\mathrm{d}t}\hat{a}(t)=(+i\Delta-\kappa/2)\hat{a}(t)-iG\hat{x}(t)\hat{a}(t)+\sqrt{\kappa_{\rm
in}}s_{\rm in}(t)+\sqrt{\kappa_{\rm 0}}\delta\hat{s}_{\rm vac}(t)\;,$ (13)

$\frac{\mathrm{d}}{\mathrm{d}t}\hat{x}(t)=\frac{\hat{p}(t)}{m_{\rm eff}}\;,$
(14)

$\frac{\mathrm{d}}{\mathrm{d}t}\hat{p}(t)=-m_{\rm eff}\Omega^{2}_{\rm
m}\hat{x}(t)-\hbar G\hat{a}^{\dagger}(t)\hat{a}(t)-\Gamma_{\rm
m}\hat{p}(t)+\delta\hat{F}_{\rm th}(t)\;,$ (15)

where $\delta\hat{s}_{\rm vac}(t)$ and $\delta\hat{F}_{\rm th}(t)$ are the
quantum and thermal noise respectively, $\kappa$ is the total cavity decay
rate, $\kappa_{\rm 0}$ the internal cavity decay rate, and $\Gamma_{\rm m}$
the mechanical decay rate. The static cavity field and mechanical displacement
are then

$\bar{a}=\frac{\sqrt{\kappa_{\rm
in}}}{\kappa/2-i(\Delta-G\bar{x})}\bar{s}_{\rm in}=\frac{\sqrt{\kappa_{\rm
in}}}{\kappa/2-i\bar{\Delta}}\bar{s}_{\rm in}\;,$ (16)

$\bar{x}=\frac{-\hbar G\absolutevalue{\bar{a}}^{2}}{m_{\rm eff}\Omega_{\rm
m}^{2}}\;.$ (17)

Here, $\bar{\Delta}=\Delta-G\bar{x}$ is the ‘corrected detuning’, accounting
for the average radiation pressure force. For a single probe beam and weak or
detuned control fields, $\bar{a}$ can be assumed real Weis  _et al._ (2010);
however, for multiple probe beams $\bar{a}$ must be considered complex.
Linearizing Eqs. 13, 14 and 15 for small perturbations using the ansatz
$\hat{a}(t)=\bar{a}+\delta\hat{a}(t)$ and
$\hat{x}(t)=\bar{x}+\delta\hat{x}(t)$, while retaining only first-order terms
in the small quantities $\delta\bar{a}$, $\delta\bar{a}^{\dagger}$ and
$\delta\bar{x}$ gives

$\frac{\mathrm{d}}{\mathrm{d}t}\delta\hat{a}(t)=(+i\bar{\Delta}-\kappa/2)\delta\hat{a}(t)-iG\bar{a}\delta\hat{x}(t)+\sqrt{\kappa_{\rm
in}}\delta s_{\rm in}(t)+\sqrt{\kappa_{\rm 0}}\delta\hat{s}_{\rm vac}(t)\;,$
(18)

$\frac{\mathrm{d}^{2}}{\mathrm{d}t^{2}}\delta\hat{x}(t)+\Gamma_{\rm
m}\frac{\mathrm{d}}{\mathrm{d}t}\delta\hat{x}(t)+\Omega_{\rm
m}^{2}\delta\hat{x}(t)=-\frac{\hbar G}{m_{\rm
eff}}(\bar{a}^{\dagger}\delta\hat{a}(t)+\bar{a}\delta\hat{a}^{\dagger}(t))+\delta\hat{F}_{\rm
th}(t)\;.$ (19)

Now considering two probe tones produced by amplitude modulation of the form
$\cos(\Omega_{\rm p}t)$, where $\pm\Omega_{\rm p}=\omega_{\rm
p}^{\pm}-\omega_{\rm d}$ is the separation in frequency between the pump and
probe tones $\omega_{\rm p}^{\pm}$, the perturbation to the incident field can
be written as

$\delta s_{\rm in}=\frac{s_{\rm p}}{2}\left(\mathrm{e}^{-i\Omega_{\rm
p}t}+\mathrm{e}^{+i\Omega_{\rm p}t}\right)\;.$ (20)

Considering the microwave drives as classical coherent fields, expectation
values can be used in place of operators
$z(t)\equiv\expectationvalue{\hat{z}(t)}$, allowing both noises, which average
to zero, to be dropped. For the form of $\delta s_{\rm in}$ in Eq. 20, a
general solution can be found using the ansatz

$\delta a(t)=A^{+}\mathrm{e}^{-i\Omega_{\rm
p}t}+A^{-}\mathrm{e}^{+i\Omega_{\rm p}t}\;,$ (21)

$\delta a^{*}(t)=(A^{-})^{*}\mathrm{e}^{-i\Omega_{\rm
p}t}+(A^{+})^{*}\mathrm{e}^{+i\Omega_{\rm p}t}\;,$ (22)

$\delta x(t)=X\mathrm{e}^{-i\Omega_{\rm p}t}+X^{*}\mathrm{e}^{+i\Omega_{\rm
p}t}\;,$ (23)

where Hermitian $\delta x(t)$ has been used. Note $A^{+}$ is used for
$\exp(-i\Omega_{\rm p}t)$ to refer to the upper sideband with
$\exp(-i(\omega_{\rm d}+\Omega_{\rm p})t)$. Substituting this ansatz into Eq.
18 and Eq. 19 and sorting by frequency will yield six equations; for a single
probe tone, only the three terms at the probe frequency are required Weis  _et
al._ (2010), while for two probe tones, all six are necessary. The two $\delta
a(t)$ terms are

$(-i(\bar{\Delta}+\Omega_{\rm p})+\kappa/2)A^{+}=-iG\bar{a}X+\sqrt{\kappa_{\rm
in}}s_{\rm p}/2\;,$ (24)

$(-i(\bar{\Delta}-\Omega_{\rm
p})+\kappa/2)A^{-}=-iG\bar{a}X^{*}+\sqrt{\kappa_{\rm in}}s_{\rm p}/2\;,\;$
(25)

while the $\delta a^{*}(t)$ terms are

$(+i(\bar{\Delta}-\Omega_{\rm
p})+\kappa/2)(A^{-})^{*}=+iG\bar{a}^{*}X+\sqrt{\kappa_{\rm in}}s_{\rm p}/2\;,$
(26)

$(+i(\bar{\Delta}+\Omega_{\rm
p})+\kappa/2)(A^{+})^{*}=+iG\bar{a}^{*}X^{*}+\sqrt{\kappa_{\rm in}}s_{\rm
p}/2\;,$ (27)

and finally the $\delta x(t)$ terms are

$m_{\rm eff}(\Omega^{2}_{\rm m}-\Omega_{\rm p}^{2}-i\Gamma_{\rm m}\Omega_{\rm
p})X=-\hbar G(\bar{a}^{*}A^{+}+\bar{a}(A^{-})^{*})\;,$ (28)

$m_{\rm eff}(\Omega^{2}_{\rm m}-\Omega_{\rm p}^{2}+i\Gamma_{\rm m}\Omega_{\rm
p})X^{*}=-\hbar G(\bar{a}^{*}A^{-}+\bar{a}(A^{+})^{*})\;.$ (29)

These six equations can be used to solve for $A^{+}$ and $A^{-}$, but first,
we define the cavity and mechanical susceptibilities as

$\chi_{\mbox{\scriptsize$\pm$}}=\frac{1}{-i(\bar{\Delta}\pm\Omega_{\rm
p})+\kappa/2}\;,$ (30)

$\chi_{\rm m}=\frac{1}{m_{\rm eff}}\frac{1}{\Omega_{\rm m}^{2}-\Omega_{\rm
p}^{2}-i\Omega_{\rm p}\Gamma_{\rm m}}\;.$ (31)

Solving Eqs. 24, 26 and 28 for $A^{+}$, and substituting in the
susceptibilities, gives

$\sqrt{\kappa_{\rm in}}s_{\rm p}/2-\frac{A^{+}}{\chi_{\rm+}}=-i\hbar
G^{2}\absolutevalue{\bar{a}}^{2}\chi_{\rm
m}\left[A^{+}+\chi_{\rm-}^{*}\left(1+\frac{\bar{a}^{2}}{\absolutevalue{\bar{a}}^{2}}\right)\sqrt{\kappa_{\rm
in}}s_{\rm p}/2\;-\;\frac{\chi_{\rm-}^{*}}{\chi_{\rm+}}A^{+}\right]\;,$ (32)

which allows for the substitution near the mechanical resonance ($\Omega_{\rm
m}\approx\Omega_{\rm p}$) of

$\begin{split}\hbar G^{2}\absolutevalue{\bar{a}}^{2}\chi_{\rm m}&=\frac{\hbar
G^{2}\absolutevalue{\bar{a}}^{2}}{m_{\rm eff}}\frac{1}{\Omega_{\rm
m}^{2}-\Omega_{\rm p}^{2}-i\Omega_{\rm p}\Gamma_{\rm m}}\\\
&\approx\frac{\hbar G^{2}\absolutevalue{\bar{a}}^{2}}{2\Omega_{\rm m}m_{\rm
eff}}\frac{1}{\Omega_{\rm m}-\Omega_{\rm p}-i\Gamma/2}\\\
&=\frac{\absolutevalue{g}^{2}}{-\delta-i\Gamma_{\rm m}/2}\;,\end{split}$ (33)

using $\delta=\Omega_{\rm p}-\Omega_{\rm m}$ and $g=\bar{a}g_{\rm 0}$, the
photon-enhanced electromechanical coupling in the linearized regime, where
$g_{\rm 0}=x_{\rm zpf}G$ is the vacuum electromechanical coupling, equivalent
to the shift in cavity resonance due to the zero point motion of the
mechanical resonator. The standard definition of zero point motion $x_{\rm
zpf}=\sqrt{\hbar/2\Omega_{\rm m}m_{\rm eff}}$ is used here. Substituting Eq.
33 into Eq. 32 and rearranging gives

$A^{+}=\chi_{\rm+}\left(1+\frac{i\absolutevalue{g}^{2}\left(\frac{\bar{a}^{2}}{\absolutevalue{\bar{a}}^{2}}\chi_{\rm-}^{*}+\chi_{\rm+}\right)}{-\delta-i\Gamma_{\rm
m}/2\;+\;i\absolutevalue{g}^{2}(\chi_{\rm-}^{*}-\chi_{\rm+})}\right)\sqrt{\kappa_{\rm
in}}s_{\rm p}/2\;.$ (34)

The same method can be used to solve for $A^{-}$ using Eqs. 25, 27 and 29,
resulting in

$A^{-}=\chi_{\rm-}\left(1+\frac{i\absolutevalue{g}^{2}\left(\frac{\bar{a}^{2}}{\absolutevalue{\bar{a}}^{2}}\chi_{\rm+}^{*}+\chi_{\rm-}\right)}{-\delta+i\Gamma_{\rm
m}/2\;+\;i\absolutevalue{g}^{2}(\chi_{\rm+}^{*}-\chi_{\rm-})}\right)\sqrt{\kappa_{\rm
in}}s_{\rm p}/2\;.$ (35)

These forms are similar to other OMIT/A experiments Weis  _et al._ (2010);
Zhou  _et al._ (2013). To obtain forms equivalent to Kashkanova et al.
Kashkanova  _et al._ (2017b) define the electromechanical self-energy as

$i\Sigma(\Omega_{\rm
p})=\absolutevalue{g}^{2}(\chi_{\rm+}-\chi_{\rm-}^{*})\;.$ (36)

This also defines the electromechanical spring effect $\Delta\Omega_{\rm
m(ele)}=\mathbb{R}\left\\{\Sigma(\Omega_{\rm p})\right\\}$, and the
electromechanical damping $\Gamma_{\rm
m(ele)}=\mathbb{I}\left\\{\Sigma(\Omega_{\rm p})\right\\}$. Using the self-
energy and $g=\bar{a}/\bar{a}^{*}g^{*}$, Eq. 34 and Eq. 35 can be rewritten as

$A^{+}=\chi_{\rm+}\sqrt{\kappa_{\rm in}}s_{\rm
p}/2\;-\;\frac{g\chi_{\rm+}(g\chi_{\rm-}^{*}+g^{*}\chi_{\rm+})}{\Gamma_{\rm
m}/2-i\delta+i\Sigma(\Omega_{\rm p})}\sqrt{\kappa_{\rm in}}s_{\rm p}/2\;,$
(37)

$A^{-}=\chi_{\rm-}\sqrt{\kappa_{\rm in}}s_{\rm
p}/2\;+\;\frac{g\chi_{\rm-}(g\chi_{\rm+}^{*}+g^{*}\chi_{\rm-})}{\Gamma_{\rm
m}/2+i\delta-i\Sigma(-\Omega_{\rm p})}\sqrt{\kappa_{\rm in}}s_{\rm p}/2\;,$
(38)

which are $A^{\pm}$ in equation (6) of the main paper, effectively the cavity
field amplitudes at $\omega_{\rm p}^{\pm}$ and equivalent to
$a_{\rm+}[\Omega]$ and $a_{\rm-}[\Omega]$ in Ref. Kashkanova  _et al._
(2017b).

###  A.2 Canceling The Probe Tones

This subsection will take the more general expressions for microwave field
amplitude $A^{+}$ and $A^{-}$, and develop the specific technique of three-
tone measurement used in this work. Using the standard input-output
relationship for a single port cavity $\hat{a}_{\rm out}=\hat{a}_{\rm
in}-\sqrt{\kappa_{\rm in}}\hat{a}$ Aspelmeyer, Kippenberg, and Marquardt
(2014), the choice of $\delta s_{\rm in}$ in Eq. 20, and using $s_{\rm
in}(t)=\left(\bar{s}_{\rm in}+\delta s(t)\right)\mathrm{e}^{-i\omega_{\rm
d}t}$, the signal reflected from the cavity $s_{\rm out}(t)$ is

$\begin{split}s_{\rm out}(t)=&\;(\bar{s}_{\rm in}-\sqrt{\kappa_{\rm
in}}\bar{a})\mathrm{e}^{-i(\omega_{\rm d}t+\phi)}\\\ &+(s_{\rm
p}/2-\sqrt{\kappa_{\rm in}}A^{+})\mathrm{e}^{-i([\omega_{\rm d}+\Omega_{\rm
p}]t+\phi+\psi)}\\\ &+(s_{\rm p}/2-\sqrt{\kappa_{\rm
in}}A^{-})\mathrm{e}^{-i([\omega_{\rm d}-\Omega_{\rm
p}]t+\phi-\psi)}\;.\end{split}$ (39)

Demodulating this signal with frequency $\omega_{\rm d}$ and phase $\phi_{\rm
ref}$, using a frequency mixer that effectively multiplies by
$\cos(\omega_{\rm d}t+\phi_{\rm ref})$, then filtering out the $2\omega_{\rm
d}$ and DC terms gives

$s_{\rm demod}(t)=\frac{K_{\rm T}}{2}\left[(s_{\rm p}-\sqrt{\kappa_{\rm
in}}A^{+})\mathrm{e}^{-i(\Omega_{\rm p}t+\psi)}+(s_{\rm p}-\sqrt{\kappa_{\rm
in}}A^{-})\mathrm{e}^{+i(\Omega_{\rm p}t+\psi)}\right]\mathrm{e}^{i(\phi_{\rm
ref}-\phi)}\;,$ (40)

where $k_{\rm p}$ is the transfer function of the circuit after the cavity,
including the mixer. The additional 10 kHz amplitude modulation (far from
$\Omega_{\rm m}$) sets $\phi_{\rm ref}-\phi=\pi/2$. A double-balanced mixer is
used to minimize loss, rather than an IQ mixer as the probe tones in the
opposing quadrature would drown out any coherent signal. As a double-balanced
mixer is used only the real part of the signal is kept; therefore, expanding
the exponentials in $s_{\rm demod}(t)$ and keeping only the real terms gives

$\begin{split}s_{\rm demod}^{\mathbb{R}}(t)/K_{\rm
T}&=-\mathbb{R}\left\\{i(A^{+}+A^{-})\right\\}\frac{\sqrt{\kappa_{\rm
in}}}{2}\cos(\Omega_{\rm
p}t+\psi)+\mathbb{R}\left\\{A^{-}-A^{+}\right\\}\frac{\sqrt{\kappa_{\rm
in}}}{2}\sin(\Omega_{\rm p}t+\psi)\\\
&=\mathbb{I}\left\\{A^{+}+A^{-}\right\\}\frac{\sqrt{\kappa_{\rm
in}}}{2}\cos(\Omega_{\rm
p}t+\psi)+\mathbb{R}\left\\{A^{-}-A^{+}\right\\}\frac{\sqrt{\kappa_{\rm
in}}}{2}\sin(\Omega_{\rm p}t+\psi)\;,\end{split}$ (41)

which is the signal the lock-in amplifier then measures.

###  A.3 Unresolved-Sideband Limit

The key regime of interest here is deep into sideband-unresolved
($\kappa\gg\Omega_{\rm m}$), and while the technique is also effective for
stronger couplings this work focuses on weak an electromechanical coupling
($\Delta\Omega_{\rm m(ele)}\ll\Gamma_{\rm m}$), where standard EMIT/A cannot
provide sufficient sensitivity. A sideband-resolved two-probe measurement
scheme can be found in Ref. Buters  _et al._ (2017), where the cavity filters
out the pump and lower probe tone. To simplify the measured signal for the
sideband-unresolved regime, first approximate the cavity susceptibilities
(using $\Omega_{\rm p}^{2},\Omega_{\rm p}\Delta\ll\kappa^{2}$) to

$\chi_{\mbox{\scriptsize$\pm$}}\approx\frac{\kappa/2+i(\Delta\pm\Omega_{\rm
p})}{\kappa^{2}/4+\Delta^{2}}\;,$ (42)

where $\bar{\Delta}\rightarrow\Delta$ has been used, as the coupling is
expected to be weak. Using the approximation in Eq. 42 the terms of interest
$\mathbb{I}\left\\{A^{+}+A^{-}\right\\}$ and
$\mathbb{R}\left\\{A^{-}-A^{+}\right\\}$ can be written:

$\mathbb{I}\left\\{A^{+}+A^{-}\right\\}=\frac{\Delta\sqrt{\kappa_{\rm
in}}s_{\rm p}}{\kappa^{2}/4+\Delta^{2}}+\frac{\sqrt{\kappa_{\rm in}}s_{\rm
p}}{2}\;\mathbb{I}\left\\{\frac{g\chi_{\rm-}(g\chi_{\rm+}^{*}+g^{*}\chi_{\rm-})}{\Gamma_{\rm
m}/2+i\delta-i\Sigma(-\Omega_{\rm
p})}-\frac{g\chi_{\rm+}(g\chi_{\rm-}^{*}+g^{*}\chi_{\rm+})}{\Gamma_{\rm
m}/2-i\delta+i\Sigma(\Omega_{\rm p})}\right\\}\;,$ (43)

$\mathbb{R}\left\\{A^{-}-A^{+}\right\\}=\frac{\sqrt{\kappa_{\rm in}}s_{\rm
p}}{2}\;\mathbb{R}\left\\{\frac{g\chi_{\rm-}(g\chi_{\rm+}^{*}+g^{*}\chi_{\rm-})}{\Gamma_{\rm
m}/2+i\delta-i\Sigma(-\Omega_{\rm
p})}+\frac{g\chi_{\rm+}(g\chi_{\rm-}^{*}+g^{*}\chi_{\rm+})}{\Gamma_{\rm
m}/2-i\delta+i\Sigma(\Omega_{\rm p})}\right\\}\;.$ (44)

To expand further, use the full form of the photon-enhanced coupling

$\begin{split}g&=\frac{g_{\rm 0}\sqrt{\kappa_{\rm in}}\bar{s}_{\rm
in}}{\kappa/2-i\Delta}\\\ &=g_{\rm 0}\sqrt{\kappa_{\rm in}}\bar{s}_{\rm
in}\frac{\kappa/2+i\Delta}{\kappa^{2}/4+\Delta^{2}}\;,\end{split}$ (45)

along with Eq. 42 and the expected small self-energy $\Sigma(\Omega_{\rm
p})\ll\Gamma_{\rm m}$, to expand the factors within
$\mathbb{I}\left\\{...\right\\}$ and $\mathbb{R}\left\\{..\right\\}$ in Eq. 43
and Eq. 44, giving

$\frac{g(g\chi_{\rm-}\chi_{\rm+}^{*}+g^{*}\chi_{\rm-}^{2})}{\Gamma_{\rm
m}/2+i\delta}=g_{\rm 0}^{2}\kappa_{\rm in}\bar{s}_{\rm
in}^{2}(\kappa/2+i\Delta)\frac{\kappa^{3}/4+\kappa\Delta^{2}+\kappa\Delta\Omega_{\rm
p}+i(2\Delta^{3}+\kappa^{2}\Delta/2-\kappa^{2}\Omega_{\rm
p}-2\Delta^{2}\Omega_{\rm
p})}{\left(\kappa^{2}/4+\Delta^{2}\right)^{4}(\Gamma/2+i\delta)}\;,$ (46)

$\frac{g(g\chi_{\rm+}\chi_{\rm-}^{*}+g^{*}\chi_{\rm+}^{2})}{\Gamma_{\rm
m}/2-i\delta}=g_{\rm 0}^{2}\kappa_{\rm in}\bar{s}_{\rm
in}^{2}(\kappa/2+i\Delta)\frac{\kappa^{3}/4+\kappa\Delta^{2}-\kappa\Delta\Omega_{\rm
p}+i(2\Delta^{3}+\kappa^{2}\Delta/2+\kappa^{2}\Omega_{\rm
p}+2\Delta^{2}\Omega_{\rm
p})}{\left(\kappa^{2}/4+\Delta^{2}\right)^{4}(\Gamma/2-i\delta)}\;,$ (47)

where factors of $\Omega_{\rm p}^{2}\;(\ll\kappa^{2})$ have been discarded.
Substituting these equations into Eq. 43 and Eq. 44, then taking the imaginary
and real parts respectively gives

$\mathbb{I}\left\\{A^{+}+A^{-}\right\\}=\frac{\Delta\sqrt{\kappa_{\rm
in}}s_{\rm p}}{\kappa^{2}/4+\Delta^{2}}-\frac{\sqrt{\kappa_{\rm in}}s_{\rm
p}}{2}g_{\rm 0}^{2}\kappa_{\rm in}\bar{s}_{\rm
in}^{2}\frac{(\kappa^{3}\Omega_{\rm p})\Gamma_{\rm
m}/2+(\kappa^{4}/4-4\Delta^{4})\delta}{\left(\kappa^{2}/4+\Delta^{2}\right)^{4}(\Gamma^{2}/4+\delta^{2})}\;,$
(48)

$\mathbb{R}\left\\{A^{-}-A^{+}\right\\}=\frac{\sqrt{\kappa_{\rm in}}s_{\rm
p}}{2}g_{\rm 0}^{2}\kappa_{\rm in}\bar{s}_{\rm
in}^{2}\frac{(\kappa^{4}/4-4\Delta^{4})\Gamma_{\rm m}/2-(\kappa^{3}\Omega_{\rm
p})\delta}{\left(\kappa^{2}/4+\Delta^{2}\right)^{4}(\Gamma^{2}/4+\delta^{2})}\;.$
(49)

Ignoring the first term in $\mathbb{I}\left\\{A^{+}+A^{-}\right\\}$ as a
constant background which can be subtracted, and using $\Omega_{\rm
p}\ll\kappa$ to remove all $\Omega_{\rm p}$ terms, $s_{\rm
demod}^{\mathbb{R}}(t)$ can be written as

$s_{\rm demod}^{\mathbb{R}}(t)=K_{\rm T}g_{\rm 0}^{2}\kappa_{\rm in}^{2}s_{\rm
p}\bar{s}_{\rm
in}^{2}\frac{(\kappa^{2}/4-\Delta^{2})}{\left(\kappa^{2}/4+\Delta^{2}\right)^{3}(\Gamma^{2}/4+\delta^{2})}[\sin(\Omega_{\rm
p}t+\psi)\Gamma_{\rm m}/2-\cos(\Omega_{\rm p}t+\psi)\delta]\;,$ (50)

which is the signal reaching the measurement lock-in amplifier. The lock-in
takes a signal of form $V_{\rm s}(t)=A_{\rm s}\cos(\Omega_{\rm s}t+\theta)$
and returns

$X+iY=\frac{A_{\rm s}}{\sqrt{2}}F(\Omega_{\rm s}-\Omega_{\rm
ref})\exp{i[(\Omega_{\rm s}-\Omega_{\rm ref}t+\theta]}\;.$ (51)

So measuring $s_{\rm demod}^{\mathbb{R}}(t)$ at frequency $\Omega_{\rm p}$
will return a signal

$X+iY=\frac{-iK_{\rm T}g_{\rm 0}^{2}\kappa_{\rm in}^{2}s_{\rm
p}}{\sqrt{2}(\Gamma/2+i\delta)}\frac{(\kappa^{2}/4-\Delta^{2})}{\left(\kappa^{2}/4+\Delta^{2}\right)^{3}}\frac{P_{\rm
in}}{\hbar\omega_{\rm d}}\;,$ (52)

where $\bar{s}_{\rm in}^{2}=P_{\rm in}/\hbar\omega_{\rm d}$. This finally
recovers equation (10) in the main text, though $K_{\rm T}$ is still unknown
which prevents calibration of $g_{\rm 0}$ from a single measurement. The phase
$\mathrm{e}^{i\psi}$ has been dropped as it will be small for a low mechanical
frequency; in practice, adjustment of the lock-in phase can always set
$\psi=0$.

###  A.4 Calibration Of Vacuum Electromechanical Coupling

To calibrate the three-tone EMIT/A response a measurement of the probe tones
interfering constructively is used. This allows for a correction of equation
(10) in the main text for $K_{\rm T}$ and $s_{\rm p}$, enabling calculation of
$g_{\rm 0}$. The calibration is achieved by setting $\phi_{\rm ref}-\phi=0$ in
Eq. 40 for $s_{\rm demod}$, giving

$s_{\rm demod}^{\mathbb{R},cal}(t)/K_{\rm T}=\mathbb{R}\left\\{s_{\rm
p}-\sqrt{\kappa_{\rm in}}(A^{+}+A^{-})\right\\}\frac{1}{2}\cos(\Omega_{\rm
p}t+\psi)+\mathbb{I}\left\\{(A^{-}-A^{+})\right\\}\frac{\sqrt{\kappa_{\rm
in}}}{2}\sin(\Omega_{\rm p}t+\psi)\;.$ (53)

The coherent contribution to the signal from the electromechanical interaction
can be ignored here, either because it is small compared to the probe tones
(for weak coupling) or because the calibration can be taken detuned from the
mechanical resonance. The sine term is also small compared to the cosine as
$\Omega_{\rm m}\ll\kappa$. Therefore the signal measured at the lock-in is
just the combined probe tones reflected from the cavity (including the cavity
susceptibility), which can be written as

$X+iY=\frac{s_{\rm p}K_{\rm T}}{2\sqrt{2}}\left(1-\frac{\kappa_{\rm
in}\kappa/2}{\kappa^{2}/4+\Delta^{2}}\right)\mathrm{e}^{i\psi}\;.$ (54)

The reference value for calibration can be achieved by correcting for the
bracketed term, using fits of the microwave resonance and taking the
magnitude, giving

$R_{\rm cal}=\frac{\absolutevalue{X+iY}}{1-\frac{\kappa_{\rm
in}\kappa/2}{\kappa^{2}/4+\Delta^{2}}}=\frac{s_{\rm p}K_{\rm
T}}{2\sqrt{2}}\;.$ (55)

Normalizing $X+iY$ in Eq. 52 by $R_{\rm cal}$, taking the value at
$\Gamma_{\rm m}$ and rearranging, allows calculation of $g_{\rm 0}$ according
to

$g_{\rm
0}=\frac{\sqrt{\hbar\omega_{\text{d}}}}{2\kappa_{\text{in}}\sqrt{P_{\text{in}}}}\sqrt{\frac{\Gamma_{\rm
m}(\kappa^{2}/4+\Delta^{2})^{3}}{(\kappa^{2}/4-\Delta^{2})}\frac{R_{\rm\delta=0}}{R_{\rm
cal}}}\;,$ (56)

which is equation (11) in the main text. The magnitude of the coherent signal
at $\delta=0$ is $R_{\rm\delta=0}$, either from the maximum of the coherent
spectrum when sweeping $\Omega_{\rm p}$ or calculated via fitting of the same
spectrum. Backgrounds should be subtracted from calculations of
$R_{\rm\delta=0}$, as these normally depend on contributions from imperfect
phase-locking. To calculate $R_{\rm cal}$, the phase $\phi_{\rm ref}$ is swept
and the maximum of the magnitude is taken.

##  Appendix B Homodyne Measurement Circuit

Figure 6 shows the microwave measurement setup for amplitude-sensitive
calibrated homodyne measurement. The method here broadly follows that of Ref.
Gorodetksy  _et al._ (2010), adapted for microwave circuits in Ref. Spence
(2022); Kumar  _et al._ (2023). This method allows readout of the Helmholtz
resonator’s mechanical noise spectral density, as this noise imparts a time
dependence on the microwave resonance via electromechanical coupling. The
time-dependent cavity resonance then imparts noise on a pump signal, measured
in the time domain by down-mixing with an LO of the same frequency, and
Fourier transformed to give the mechanical noise spectral density. The figure
shows an amplitude-sensitive scheme, where $\omega_{\rm d}$ is placed at the
maximum gradient of the microwave resonance, the 10 kHz reference is applied
via phase modulation, and the calibration tone via amplitude modulation.
Phase-sensitive measurement was also used to measure the mechanics, where
$\omega_{\rm d}$ is placed at the maximum amplitude of the microwave
resonance, the 10 kHz reference is applied via amplitude modulation, and the
calibration tone via phase modulation.

Working backward from $g_{\rm 0}$, and using the calibrated homodyne technique
Gorodetksy  _et al._ (2010); Spence (2022), an estimate for effective
Helmholtz mode temperature is calculated giving $T_{\text{eff}}\approx 80,000$
K, which demonstrates the strong incoherent driving from the cryocooler. While
initially, this figure seems high, a silicon nitride membrane driven with
white noise via a piezo at 0.001 mV2/Hz had an effective temperature of 23,000
K Kumar  _et al._ (2023); meanwhile, driving the Helmholtz resonator at over
2.5 mV2/Hz in this work could not overcome the cryocooler.

Figure 6: Amplitude-sensitive calibrated homodyne measurement circuit. The RF
Source microwave signal, at frequency $\omega_{\rm d}$, is phase modulated
with the 10 kHz reference signal and amplitude modulated with the calibration
tone at $\Omega_{\rm calib}$. The RF signal is sent to the experimental ‘Cell’
via a directional coupler. The reflected signal passes back through the
directional coupler, is amplified with two microwave amplifiers, and is mixed
with the LO source signal at frequency $\omega_{\rm d}$. The mixed signal is
low-pass filtered to remove the $2\omega_{\rm d}$ components and is then
measured on two lock-in amplifiers. The lower lock-in on the diagram is set to
1.4 kHz, near the mechanical frequency $\Omega_{\rm m}$, and acts as a filter
plus amplifier. The $X$ and $Y$ quadrature outputs of this lock-in are
measured on a DAQ card, giving a time domain noise measurement of frequency
components close to $\Omega_{\rm m}$. The upper lock-in on the diagram
measures the 10 kHz reference signal, which is again used to lock the relative
phase of mixing via the ‘LO Source’. A piezo is present on the outside of the
cell to allow for mechanical driving, either with swept frequency to find
mechanical modes, or with white noise to raise the effective mode temperature.


11institutetext: School of Mathematics and Computational Science, Xiangtan
University
Xiangtan, Hunan 411105, China

# Accelerating Stochastic Recursive and Semi-stochastic Gradient Methods with
Adaptive Barzilai-Borwein Step Sizes

JiangShan Wang  YiMing Yang  Zheng Peng

###### Abstract

The mini-batch versions of StochAstic Recursive grAdient algoritHm and Semi-
Stochastic Gradient Descent method, employed the random Barzilai-Borwein step
sizes (shorted as MB-SARAH-RBB and mS2GD-RBB), have surged into prominence
through timely step size sequence. Inspired by modern adaptors and variance
reduction techniques, we propose two new variant rules in the paper, referred
to as RHBB and RHBB+, thereby leading to four algorithms MB-SARAH-RHBB, MB-
SARAH-RHBB+, mS2GD-RHBB and mS2GD-RHBB+ respectively. RHBB+ is an enhanced
version that additionally incorporates the importance sampling technique. They
are aggressive in updates, robust in performance and self-adaptive along
iterative periods. We analyze the flexible convergence structures and the
corresponding complexity bounds in strongly convex cases. Comprehensive tuning
guidance is theoretically provided for reference in practical implementations.
Experiments show that the proposed methods consistently outperform the
original and various state-of-the-art methods on frequently tested data sets.
In particular, tests on the RHBB+ verify the efficacy of applying the
importance sampling technique to the step size level. Numerous explorations
display the promising scalability of our iterative adaptors.

###### Keywords:

Variance reduction stochastic optimization random hedge Barzilai-Borwein
method importance sampling iterative adaptors

##  1 Introduction

We focus on the following stochastic optimization problem

$\min_{w\in\mathbb{R}^{d}}P(w)=\mathbb{E}_{\mathbb{P}}[f(w;\xi)_{\xi\in\Omega}],$

where $\xi$ is a random instance of an input-output pair $(x_{i},z_{i})$, with
input representation vector $x_{i}$ and target output $z_{i}$. Hence,
$f(\cdot)$ takes the form

$f(w;\xi)=f\left(w;(x_{i},z_{i})\right).$

Define $P(\cdot)$ by empirical expectation on probability space
$(\Omega,\mathbb{P})$, where $\Omega=\left\\{\xi_{1},\cdots,\xi_{n}\right\\}$
is a finite support set and $\mathbb{P}$ is the probability measure over
$\Omega$. In this case, it transforms into an unconstrained finite-sum model,
i.e.,

$\min_{w\in\mathbb{R}^{d}}P(w)=\int_{\Omega}f\left(w;(x_{i},z_{i})\right)d\mathbb{P}(x_{i},z_{i})\approx\frac{1}{n}\sum_{0\leq
i\leq|\Omega|}f(w;(x_{i},z_{i})).$ (1)

Problem (1) covers a broad range of applications. Assume that a lower bound
$P(w_{*})$ of $P(\cdot)$ exists, the goal is to iteratively update $w$ to
reduce $P(\cdot)$ steadily and swiftly. Given a sequence of $n$ labeled pairs
$\\{\left(x_{1},z_{1}\right),\ldots,\left(x_{n},z_{n}\right)\\}$ into finite
dimensional spaces $\\{(\mathbb{R}^{d},\mathbb{R})\\}$, the linear least
squares regression is of $f(w)\stackrel{{\scriptstyle\text{ def }}}{{=}}$
$\left(x_{i}^{T}w-z_{i}\right)^{2}$. In terms of the logistic regression, we
exploit $f(w)\stackrel{{\scriptstyle\text{ def
}}}{{=}}\log\left(1+\exp\left(-z_{i}x_{i}^{T}w\right)\right)$. Indeed,
efficient regularization may be taken into account for specific purposes, then
it develops into a composite model, i.e.,

$\min_{w\in\mathbb{R}^{d}}\frac{1}{n}\sum_{0\leq
i\leq|\Omega|}f(w;(x_{i},z_{i}))+\psi(w),$

where $\psi(\cdot)$ is a proper, closed and convex penalty on parameters. Its
subdifferential at $w$ defines
$\partial\psi(w)=\left\\{v\in\mathbb{R}^{d}\mid\psi(d)\geq\psi(w)+v^{T}\left(d-w\right),\forall
d\in\operatorname{dom}(\psi)\right\\}$. Throughout the paper, we mainly
utilize a smoothing regularization process, thus $\psi(\cdot)$ is
differentiable and $\partial\psi(w)=\\{\nabla\psi(w)\\}$.

Due to the productive and scalable frameworks, stochastic optimization is
prevalent among large-scale problems or complex relationships. It offers a
cost-effective alternative to deterministic schemes, stated in a universal
form as

$w^{(t)}=w^{(t-1)}-\eta_{t}g_{t}\left(w^{(t-1)},\nu_{t}\right).$

Here, $g(\cdot)$ represents the gradient estimator, $\nu$ the randomness
pointer. Vanilla stochastic gradient descent (SGD) [Robbins][shalev] specifies
$g(\cdot)$ into the basic moment estimate, where it enjoys cheap computational
cost per update and the independence with $n$ in term of complexity. However,
this straightforward estimator inevitably introduces variance and noise to the
steps. Diminishing step size $\eta_{t}=O(1/t)$ is then forced to employed in
[luo][solodov][nemirovski] for a sublinear convergence rate of $O(1/t)$
(Moulines et al.[moulines]), which should satisfy

$\sum_{t=1}^{\infty}\eta_{t}=\infty\text{ and
}\sum_{t=1}^{\infty}\eta_{t}^{2}<\infty.$

Such updates come with a side effect of halt in convergence near the eventual
limit[eon]. Batch methods of [dekel][cotter][shalev] decrease the intrinsic
variance through a bunch of samples, but at the cost of further computational
workload. Therefore, parallel processing power becomes indispensable for
explosion-scale data.

Shifts on $g(\cdot)$ have been intensively studied with a vast number of
research papers. Canonical examples include but not limit to SAG/SAGA
[Schmidt][defazio], SVRG [johnson], MSVRG[reddi], S2GD [konevcny2], mS2GD
[konevcny], MISO [mairal], SARAH[nguyen2], iSARAH [nguyen3], MB-SARAH
[nguyen2], SPIDER [fang], etc. Konecny et al. [konevcny] proposed mS2GD and
showed it reaches a pre-defined accuracy with less overall work than a method
without batching. They established a threshold for the batch size, at which
more than linear speedup can be achieved, it’s worthy to further explore. MB-
SARAH is presented by Nguyen et al. [nguyen2] for solving non-convex problems.
Recursive updates free from the storage of past gradients and avoid
oscillation of the Euclidean norm of $g(\cdot)$ in inner loops, which are
well-suited for modern complex scenarios. Researches have also extended to the
dual space of (1) by updating random dual variables or variable blocks, such
as RCDM[nesterov2], AsySPDC[liu2], SDCA[shalev2], SPDC[zhang], mSDCA[takac],
ASDCA[shalev3] and QUARTZ[zheng], Prox-SDCA [zhao]. The subsequent drawback
lies in that these algorithms rely on a tuning step size by hand, which can be
time consuming in practice.

Several methods of auto-tuning prevail among the stochastic algorithms.
Barzilai-Borwein method (BB) of the second order tuning is outstanding in the
trend, due to its simplicity and numerical efficiency. Sopyła et al. [sopyla]
employed the BB in vanilla SGD to solve linear SVM in dual space. Tan et al.
[tan] incorporated the BB into SVRG (SVRG-BB) and established the linear
convergence in strongly convex cases. To further accelerate the rates, Yang et
al. [yang3] introduced the BB to mS2GD (mS2GD-BB) for nonsmooth and strongly
convex functions. On the basis of Hessian and its eigenvalues, Ma et al. [ma]
proposed the stabilized Barzilai-Borwein (SBB) method to match SVRG (SVRG-SBB)
for the ordinal embedding problems, which avoids the denominator tending to
zero. Yang et al. [yangy] considered the inexact SARAH (iSARAH-BB) in order to
reduce the cost in deterministic steps, then showed its robustness in
implementations. Byrd et al. [byrd] utilized batch methods to approximate the
quasi-Newton property. Recently, the timely random Barzilai-Borwein method
(RBB) emerged and was primarily applied in MB-SARAH (MB-SARAH-RBB[yang]) and
mS2GD (mS2GD-RBB[yang2]) algorithmic settings. The promising performance
outperformed and matched state-of-the-art algorithms. However, they still have
flaws and aspects that can be explored and improved.

In the context of the RBB rule, when the batch size is insufficient, it can
result in a high level of noise, primarily caused by the random step sizes.
This noise leads to an increasing or oscillating trend towards divergence. As
batch size increases, the variance of stochastic curvature decreases, but
using a larger batch discouragingly slows down the convergence. In addition,
it’s insensitive to the iterative periods, and the well-worn uniform sampling
deserves to be generalized. Therefore, we introduce the random hedge Barzilai-
Borwein method (RHBB) in pursuit of improvement. We further incorporate our
RHBB with the importance sampling technique and develop another enhanced
version, RHBB+.

The key contributions in the paper are summarized as follows:

  1. 1)

We propose the random hedge Barzilai-Borwein step size rule (RHBB) for MB-
SARAH and mS2GD algorithmic settings, obtaining MB-SARAH-RHBB and mS2GD-RHBB
algorithms. The adaptive acceleration mechanism is analyzed, trade-off rules
are studied and the tuning guidance is provided.

  2. 2)

We incorporate the importance sampling technique into the RHBB to make use of
the distributed efficiency in data sets (e.g. the elementwise sparsity), which
yields the enhanced version RHBB+ and the corresponding MB-SARAH-RHBB+ and
mS2GD-RHBB+ algorithms.

  3. 3)

We establish the global convergence of all four algorithms in strongly convex
cases. Further, we analyze the ill-conditioned scenario to theoretically
validate the robustness of the new algorithms. In MB-SARAH-RHBB/RHBB+, the
square of the full gradient converges linearly in expectation. And the
expected distance of iterates to the global optimum has linear convergence in
mS2GD-RHBB/RHBB+.

  4. 4)

We conduct extensive experiments to demonstrate the exceptional performance of
our algorithms. Next, we explore a tentative, incremental scheme for the
iterative adaptor and view immense potential for its scalability.

##  2 Common Assumptions and Inequalities

We add subscripts to distinguish element functions, e.g. $f_{i}$ denotes the
$i$-th component. Unless otherwise specified, $\|\cdot\|$ denotes Euclidean
norm in this paper. Subsequently, we provide following common assumptions.
Assumption 1 (Smoothness). Each convex
$f_{i}:\mathbb{R}^{d}\rightarrow\mathbb{R}$ is $L$-smooth over any compact set
of its domain, i.e., there exists an $L>0$, for all
$w,w^{\prime}\in\text{dom}(f_{i})$

$\left\|\nabla f_{i}(w)-\nabla f_{i}\left(w^{\prime}\right)\right\|\leq
L\left\|w-w^{\prime}\right\|.$

Due to $\nabla f_{i}(w)-\nabla
f_{i}\left(w^{\prime}\right)=H\left(w-w^{\prime}\right)$ where
$H=\nabla^{2}f_{i}\left(\hat{w}\right)$ is the Hessian at midpoint $\hat{w}$,
the largest eigenvalue of $H$ is no more than $L$. In view of it, we obtain a
boundary for the square distance between the gradients, i.e.,

$\displaystyle\left\|\nabla f_{i}(w)-\nabla
f_{i}\left(w^{\prime}\right)\right\|^{2}$
$\displaystyle=\left(w-w^{\prime}\right)^{\top}H^{2}\left(w-w^{\prime}\right)$
(2) $\displaystyle\leq
L\left(w-w^{\prime}\right)^{\top}H\left(w-w^{\prime}\right)$
$\displaystyle=L\left(\nabla f_{i}(w)-\nabla
f_{i}\left(w^{\prime}\right)\right)^{\top}\left(w-w^{\prime}\right).$

The individual $L$-smoothness implies the overall $P(w)=\frac{1}{n}\sum_{0\leq
i\leq|\Omega|}f_{i}(w)$ is $L$-smooth as well. Equivalently, we derive the
following bound

$P(w)\leq P(w^{\prime})+\nabla
P(w^{\prime})^{T}(w-w^{\prime})+\frac{L}{2}\|w-w^{\prime}\|^{2}.$ (3)

For the $L$-smoothness of the overall $P(\cdot)$, we have another equivalent
claim, i.e.,

$P(w)\geq P(w^{\prime})+\nabla P(w)^{T}(w-w^{\prime})+\frac{1}{2L}\|\nabla
P(w)-\nabla P(w^{\prime})\|^{2}.$ (4)

Assumption 2a (Strong convexity I). $P(w)$ is $\mu$-strongly convex, i.e.,
there exists $\mu>0$ such that, for all $w,w^{\prime}\in\text{dom}(P)$

$(\nabla P(w)-\nabla
P(w^{\prime}))^{T}(w-w^{\prime})\geq\mu\|w-w^{\prime}\|^{2}.$ (5)

Assumption 2b (Strong convexity II). Each component $f_{i}$ is $\mu$-strongly
convex, i.e., there exists $\mu>0$ for each $f_{i}$ such that for all
$w,w^{\prime}\in\text{dom}(f_{i})$

$(\nabla f_{i}(w)-\nabla
f_{i}(w^{\prime}))^{T}(w-w^{\prime})\geq\mu\|w-w^{\prime}\|^{2}.$ (6)

Assumption 2b can imply Assumption 2a, but not vice versa. Assumption 2b is a
stricter premise which requires the strong convexity on each $f_{i}$.

For the $\mu$-strongly convex $P(\cdot)$, we have another equivalent claim as

$P(w)\geq P(w^{\prime})+\nabla
P(w^{\prime})^{T}(w-w^{\prime})+\frac{\mu}{2}\|w-w^{\prime}\|^{2}.$

Define the global optimum $w_{*}=\operatorname{argmin}_{w}P(w)$, it further
indicates that

$2\mu\left[P(w)-P\left(w_{*}\right)\right]\leq\|\nabla P(w)\|^{2},\forall
w\in\mathbb{R}^{d}.$ (7)

To see this, we have by strong convexity

$P\left(w_{*}\right)\geq P(w)+\nabla
P(w)^{T}\left(w_{*}-w\right)+\frac{\mu}{2}\left\|w-w_{*}\right\|^{2}.$

Then, through some basic derivation, we have

$\displaystyle 2\mu\left[P(w_{*})-P\left(w\right)\right]+\|\nabla
P(w)\|^{2}\geq 2\mu\nabla
P(w)^{T}\left(w_{*}-w\right)+\mu^{2}\left\|w-w_{*}\right\|^{2}+\|\nabla
P(w)\|^{2},$ $\displaystyle 2\mu\left[P(w_{*})-P\left(w\right)\right]+\|\nabla
P(w)\|^{2}\geq\left\|\nabla P(w)+\mu\left(w_{*}-w\right)\right\|^{2}\geq 0,$
$\displaystyle 2\mu\left[P(w)-P\left(w_{*}\right)\right]\leq\|\nabla
P(w)\|^{2}.$

##  3 Motivation

###  3.1 Barzilai-Borwein Method and the Random Versions

Barzilai-Borwein method, originally developed in the pioneer work of Barzilai
and Borwein [1], shows great preeminence in solving nonlinear optimization
problems and has been widely improved up to now.

We automatically hope $\theta_{k}I$ approximates $H_{k}$ in the $k$-th epoch,
where $\theta_{k}$ denotes the $k$-th step size, $I$ the identity matrix and
$H_{k}$ the inverse of Hessian. To minimize the residual of the secant
equations, i.e., $\left\|\left(1/\theta\right)s_{k}-y_{k}\right\|_{2}^{2}$ and
$\left\|\theta y_{k}-s_{k}\right\|_{2}^{2}$, we have the following step size
solutions respectively

$\theta_{k}^{\mathrm{BB}1}=\frac{s_{k}^{T}s_{k}}{s_{k}^{T}y_{k}},\quad\theta_{k}^{\mathrm{BB}2}=\frac{s_{k}^{T}y_{k}}{y_{k}^{T}y_{k}},$

where $s_{k}=w_{k+1}-w_{k}$, $y_{k}=\nabla f\left(w_{k+1}\right)-\nabla
f\left(w_{k}\right)$.

Among the well-known Barzilai-Borwein methods, either BB1 or BB2 is expected
to be computed at the start of each epoch and employed uniformly during the
entire consecutive stochastic stages (see [tan][yang3][ma][yangy] for a brief
reference). Yang et al. [yang][yang2] advocated to calculate the Barzilai-
Borwein step size timely by stochastic curvature (batch curvature, similar to
Castera et al. [castera]) within each stochastic stage. In mS2GD algorithmic
setting, they proposed a random version of BB1, we instinctively deduce the
random BB2 by analogy

$\tilde{\eta}_{k}^{\mathrm{RBB}1}=\frac{1}{|S_{1}|}\cdot\frac{\left\|w_{k}-w_{k-1}\right\|^{2}}{\left(\left(w_{k}-w_{k-1}\right)^{T}\left(\nabla
P_{S_{1}}\left(w_{k}\right)-\nabla
P_{S_{1}}\left(w_{k-1}\right)\right)\right)},$
$\tilde{\eta}_{k}^{\mathrm{RBB}2}=\frac{1}{|S_{2}|}\cdot\frac{\left(\left(w_{k}-w_{k-1}\right)^{T}\left(\nabla
P_{S_{2}}\left(w_{k}\right)-\nabla
P_{S_{2}}\left(w_{k-1}\right)\right)\right)}{\left\|\nabla
P_{S_{2}}\left(w_{k}\right)-\nabla
P_{S_{2}}\left(w_{k-1}\right)\right\|^{2}},$

where $\nabla P_{S_{1}}\left(w_{k}\right)=\frac{1}{|S_{1}|}\sum_{i\in
S_{1}}\nabla f_{i}\left(w_{k}\right),\nabla
P_{S_{2}}\left(w_{k}\right)=\frac{1}{|S_{2}|}\sum_{j\in S_{2}}\nabla
f_{j}\left(w_{k}\right)$. The $S_{1},S_{2}\subset\\{1,\ldots,n\\}$ are
randomly selected subsets with size $|S_{1}|$ and $|S_{2}|$. As in MB-SARAH
algorithmic setting, another parameter $\gamma$ is multiplied to adjust the
current RBB for a better Hessian approximation, i.e.,

$\eta_{k}^{\mathrm{RBB}1}=\frac{\gamma}{|S_{1}|}\cdot\frac{\left\|w_{k}-w_{k-1}\right\|^{2}}{\left(\left(w_{k}-w_{k-1}\right)^{T}\left(\nabla
P_{S_{1}}\left(w_{k}\right)-\nabla
P_{S_{1}}\left(w_{k-1}\right)\right)\right)},$
$\eta_{k}^{\mathrm{RBB}2}=\frac{\gamma}{|S_{2}|}\cdot\frac{\left(\left(w_{k}-w_{k-1}\right)^{T}\left(\nabla
P_{S_{2}}\left(w_{k}\right)-\nabla
P_{S_{2}}\left(w_{k-1}\right)\right)\right)}{\left\|\nabla
P_{S_{2}}\left(w_{k}\right)-\nabla
P_{S_{2}}\left(w_{k-1}\right)\right\|^{2}}.$

(a) a8a (b) covtype (c) phishing (d) a8a (e) covtype (f) phishing

Figure 1: (a)(b)(c): Trajectories of BB1 step sizes and RBB step sizes in MB-
SARAH. (d)(e)(f): Trajectories of BB1 step sizes and RBB step sizes in mS2GD.

Let us observe a series of trajectories in Fig. 1, where the $x$-axis denotes
the number of inner iterations within the $2$-nd epoch and the $y$-axis
represents the corresponding step sizes. All relevant parameters are set as
suggested in [yang][yang2]. For illustration, we discard the first epoch that
is insufficient to compute the BB1 step sizes (instead a specified constant is
applied during the first epoch). It can be observed that BB1 remains unvaried
in a single value throughout the entire epoch, while RBB1 and RBB2 update in a
timely manner, evolving with RBB1 taking the precedence or both intertwining
mutually.

In fact, the BB1-type is more preferable due to its aggressive finesse,
however in a wide range of cases, it still may not reach the extreme of
efficiency as well. Note that hedging is an innocuous way to mitigate risks in
the financial sector. Inspired by this, we ensure or deliberately expand the
effective magnitude of RBB1, while offsetting any over-utility from the
opposite direction via its ‘twin’ RBB2 to devise our RHBB step size rule. It
can enlarge the adaptive step sizes smoothly and controllably. Further, to
improve the adaptivity along iterative periods, we expect to include an
adaptor to adjust the step sizes. Specifically, RHBB is based on an affine
combination of RBB1, RBB2 via an adaptive parameter
$\alpha^{h(\sigma_{1}s+\sigma_{2}k)}$, where
$\alpha^{h(\sigma_{1}s+\sigma_{2}k)}>1$. The adaptor $h$ is an exponential
decay rate of the affine magnitude. In early epochs, the $h$ boosts the step
sizes at a low cost to accelerate training. As approaching the global optimum,
it then enforces them conservative to ensure the final convergence. In fact,
the adaptor $h$ is monotone decreasing with respect to the linear indicators
$\left(\sigma_{1}s+\sigma_{2}k\right)$, and iteratively satisfies

$|\alpha-h(\sigma_{1}s+\sigma_{2}k)|>\epsilon(s),\forall k\in m,\forall s.$

In practice, we make $\sigma_{1},\sigma_{2}\in\\{0,1\\}$ and
$\epsilon(s)>\frac{s}{m}$, the gap should be expanded in latter periods. It’s
of great distinction in the structural sense from the composite Barzilai-
Borwein method (CBB in [li]) and the composite adaptive Barzilai-Borwein
method (CABB in [li]), which utilize a convex combination with parameters
within range $(0,1)$. The CBB and CABB comprise two components, each of which
extracts partial resources from either BB1 or BB2. Essentially, all Barzilai-
Borwein methods enjoy the ‘calculation’ adaptivity, we attach another adaptor,
$h(\cdot)$, to enhance the adaptivity along iterations.

The employment of RBB2 introduces another type of stochastic curvature. Notice
that we capture the stochastic curvature from two probabilistic subsets
$S_{1}$ and $S_{2}$, suggesting that we use quite another subset to do
hedging. Due to the quasi-Newton property, RHBB indeed reduces the noise in
the second order level. For the sake of convergence, we pick the larger batch
size as the batch correction. Hence, we have the RHBB for MB-SARAH algorithmic
setting as (with total $\alpha^{h(\sigma_{1}s+\sigma_{2}k)}>1$):

$\displaystyle(\eta_{k}^{s})^{\mathrm{RHBB}}$
$\displaystyle=\frac{\gamma}{\max\\{|S_{1}|,|S_{2}|\\}}\cdot\left(\frac{\alpha^{h(\sigma_{1}s+\sigma_{2}k)}\cdot\left\|w_{k}^{s}-w_{k-1}^{s}\right\|^{2}}{\left(\left(w_{k}^{s}-w_{k-1}^{s}\right)^{T}\left(\nabla
P_{S_{1}}\left(w_{k}^{s}\right)-\nabla
P_{S_{1}}\left(w_{k-1}^{s}\right)\right)\right)}\right.$ (8)
$\displaystyle+\left.\frac{\left(1-\alpha^{h(\sigma_{1}s+\sigma_{2}k)}\right)\cdot\left(\left(w_{k}^{s}-w_{k-1}^{s}\right)^{T}\left(\nabla
P_{S_{2}}\left(w_{k}^{s}\right)-\nabla
P_{S_{2}}\left(w_{k-1}^{s}\right)\right)\right)}{\left\|\nabla
P_{S_{2}}\left(w_{k}^{s}\right)-\nabla
P_{S_{2}}\left(w_{k-1}^{s}\right)\right\|^{2}}\right).$

Next, we multiply an extra balance parameter, $\gamma_{2}$, to the RHBB step
sizes in mS2GD algorithmic setting, i.e.,

$\displaystyle(\tilde{\eta}_{k}^{s})^{\mathrm{RHBB}}$
$\displaystyle=\frac{\gamma_{2}}{\max\\{|S_{1}|,|S_{2}|\\}}\cdot\left(\frac{\alpha^{h(\sigma_{1}s+\sigma_{2}k)}\cdot\left\|w_{k}^{s}-w_{k-1}^{s}\right\|^{2}}{\left(\left(w_{k}^{s}-w_{k-1}^{s}\right)^{T}\left(\nabla
P_{S_{1}}\left(w_{k}^{s}\right)-\nabla
P_{S_{1}}\left(w_{k-1}^{s}\right)\right)\right)}\right.$ (9)
$\displaystyle+\left.\frac{\left(1-\alpha^{h(\sigma_{1}s+\sigma_{2}k)}\right)\cdot\left(\left(w_{k}^{s}-w_{k-1}^{s}\right)^{T}\left(\nabla
P_{S_{2}}\left(w_{k}^{s}\right)-\nabla
P_{S_{2}}\left(w_{k-1}^{s}\right)\right)\right)}{\left\|\nabla
P_{S_{2}}\left(w_{k}^{s}\right)-\nabla
P_{S_{2}}\left(w_{k-1}^{s}\right)\right\|^{2}}\right).$

In the latter section, we will elaborate the role of $\gamma_{2}$ in relaxing
the stochastic hedge effect. And related trade-off rules will as well be
studied.

We choose the batch scheme for possible opportunities of parallel processing.
Notice from Tan et al. [tan] that they applied a convex combination to
approximate the full gradient at the snapshots, hence the absolute operation
has been taken upon the denominators in the step size calculations. Here in
MB-SARAH and mS2GD algorithmic settings, we bear no worry on possible negative
step sizes. According to Castera et al.[castera], when the curvature condition
$(g_{k}^{s})^{T}\nabla^{2}P(w_{k}^{s})g_{k}^{s}$ does not keep positive
($g_{k}^{s}$ is an update direction), it’s advisable to set $\eta_{k}^{s}=c$
($c>0$). Different from AS in Liu et al. [liu3], we conduct and use real-time
estimations instead of the accumulated moving average.

###  3.2 Importance Sampling technique

In terms of theory, this technique leads to the improvement of leading
constants in the complexity estimates (Richtarik et al [richtarik], Needell et
al [needell]). The overhead associated with configuring distributions and
withdrawing samples is negligible, and hence the net effect [csiba] is
speedup.

Uniform sampling enables unbiased estimators but sacrifices potential
opportunities of variance reduction, algorithms nowadays have strived for the
opposite. Most sampling techniques have been applied to the gradient
estimates, which include but not limit to Prox-SVRG [xiao], Prox-SDCA [zhao]
and SARAH-I [liu]. Hence, we inventively exploit the sampling schemes in the
step size level, customizing probability distributions to filter stochastic
variations. In fact, we configure $Q\sim$ $\left\\{q_{1},q_{2},\ldots
q_{n}\right\\}$ according to Zhao et al [zhao].

The basic moment estimate of $\nabla P(\cdot)$, over the uniform distribution,
on the subset $S\subset\Omega$, is in the form of

$\nabla P_{S}\left(w_{k}\right)=\frac{1}{|S|}\sum_{0\leq i\leq|S|}\nabla
f_{i}\left(w_{k}\right).$ (10)

In light of (10), consider that from a general distribution we have

$\nabla P_{S}^{+}\left(w_{k}\right)=\frac{1}{|S|}\sum_{0\leq
i\leq|S|}\frac{\nabla f_{i}\left(w_{k}\right)}{n\cdot
q_{i}}=\frac{1}{|S|}\sum_{0\leq i\leq|S|}\nabla f^{+}_{i}\left(w_{k}\right).$
(11)

Uniformity delivers $q_{1}=...=q_{n}=\frac{1}{n}$, suggesting a special case
(10) of $(\ref{P+})$.

By substituting (10) with (11), we extend the RHBB (8) (9) to support general
distributions, resulting in the enhanced rule RHBB+, i.e.,

$\displaystyle(\eta_{k}^{s})^{\mathrm{RHBB+}}$
$\displaystyle=\frac{\gamma}{\max\\{|S_{1}|,|S_{2}|\\}}\cdot\left(\frac{\alpha^{h(\sigma_{1}s+\sigma_{2}k)}\cdot\left\|w_{k}^{s}-w_{k-1}^{s}\right\|^{2}}{\left(\left(w_{k}^{s}-w_{k-1}^{s}\right)^{T}\left(\nabla
P_{S_{1}}^{+}\left(w_{k}^{s}\right)-\nabla
P_{S_{1}}^{+}\left(w_{k-1}^{s}\right)\right)\right)}\right.$ (12)
$\displaystyle+\left.\frac{\left(1-\alpha^{h(\sigma_{1}s+\sigma_{2}k)}\right)\cdot\left(\left(w_{k}^{s}-w_{k-1}^{s}\right)^{T}\left(\nabla
P_{S_{2}}^{+}\left(w_{k}^{s}\right)-\nabla
P_{S_{2}}^{+}\left(w_{k-1}^{s}\right)\right)\right)}{\left\|\nabla
P_{S_{2}}^{+}\left(w_{k}^{s}\right)-\nabla
P_{S_{2}}^{+}\left(w_{k-1}^{s}\right)\right\|^{2}}\right),$
$\displaystyle(\tilde{\eta}_{k}^{s})^{\mathrm{RHBB+}}$
$\displaystyle=\frac{\gamma_{2}}{\max\\{|S_{1}|,|S_{2}|\\}}\cdot\left(\frac{\alpha^{h(\sigma_{1}s+\sigma_{2}k)}\cdot\left\|w_{k}^{s}-w_{k-1}^{s}\right\|^{2}}{\left(\left(w_{k}^{s}-w_{k-1}^{s}\right)^{T}\left(\nabla
P_{S_{1}}^{+}\left(w_{k}^{s}\right)-\nabla
P_{S_{1}}^{+}\left(w_{k-1}^{s}\right)\right)\right)}\right.$ (13)
$\displaystyle+\left.\frac{\left(1-\alpha^{h(\sigma_{1}s+\sigma_{2}k)}\right)\cdot\left(\left(w_{k}^{s}-w_{k-1}^{s}\right)^{T}\left(\nabla
P_{S_{2}}\left(w_{k}^{s}\right)-\nabla
P_{S_{2}}^{+}\left(w_{k-1}^{s}\right)\right)\right)}{\left\|\nabla
P_{S_{2}}^{+}\left(w_{k}^{s}\right)-\nabla
P_{S_{2}}^{+}\left(w_{k-1}^{s}\right)\right\|^{2}}\right).$

##  4 Algorithms

We shall first clarify the notations: $\\{v_{k}^{s}\\}$ denotes the estimate
sequence of $\nabla P(\cdot)$ in MB-SARAH-RHBB/RHBB+, each estimate is in a
recursive form with

$v_{k}^{s}=\nabla P_{S}\left(w_{k}^{s}\right)-\nabla
P_{S}\left(w_{k-1}^{s}\right)+v_{k-1}^{s}.$ (14)

In mS2GD-RHBB/RHBB+, we signify the estimate array as
$\\{\tilde{v}_{k}^{s}\\}$ and have

$\tilde{v}_{k}^{s}=\nabla P_{S}\left(w_{k}^{s}\right)-\nabla
P_{S}(\widetilde{w})+\nabla P(\widetilde{w}).$ (15)

Notably, $\\{w_{k}^{s}\\}$ represents the inner iterative sequence within the
$s$-th outer epoch, and we use $\\{\widetilde{w}_{s}\\}$ for the outer series
and $\widetilde{w}$ for the snapshots.

Input: $\widetilde{w}_{0}$, update frequency $m$, batch sizes $b,b_{1},b_{2}$,
constant sequence $\\{\eta_{0}^{s}\\}$, modification parameter $\gamma>0$,
hedge base $\alpha$ and monotone decreasing function $h$, probability
distribution $Q$.

Output: approximate solution $\widetilde{w}_{s}$.

1

2 Outer Loop: for  _$s=1,2,...,$_ do

3 $w_{0}^{s}=\widetilde{w}_{s-1}$

4 $v_{0}^{s}=\frac{1}{n}\sum_{i\in\Omega}\nabla
f_{i}\left(w_{0}^{s}\right)=\nabla P\left(w_{0}^{s}\right)$

5 $w_{1}^{s}=w_{0}^{s}-\eta_{0}^{s}v_{0}^{s}$

6 Inner Loop: for  _$k=1,2,...,m-1$_ do

7 Pick subset $S\subset\\{1,\ldots,n\\}$ of size $b$ uniformly at random

8 Update the stochastic recursive gradient $v_{k}^{s}$ by  $v_{k}^{s}=\nabla
P_{S}\left(w_{k}^{s}\right)-\nabla P_{S}\left(w_{k-1}^{s}\right)+v_{k-1}^{s}$

9 Compute the step size $\eta_{k}^{s}$ by (1) Option I or (2) Option II

10 (1) Option I : RHBB rule

11 Configure $Q$ as uniform probability distribution

12 Pick subset $S_{1}\subset\\{1,\ldots,n\\}$ of size $b_{1}$ randomly
according to $Q$

13 Pick subset $S_{2}\subset\\{1,\ldots,n\\}$ of size $b_{2}$ randomly
according to $Q$

14 Calculate $\eta_{k}^{s}$ according to (8)

15

16 (2) Option II : RHBB+ rule

17 Configure $Q$ to our needs

18 Pick subset $S_{1}\subset\\{1,\ldots,n\\}$ of size $b_{1}$ randomly
according to $Q$

19 Pick subset $S_{2}\subset\\{1,\ldots,n\\}$ of size $b_{2}$ randomly
according to $Q$

20 Compute $\nabla P_{S_{1}}^{+},\nabla P_{S_{2}}^{+}$ according to (11)

21 Calculate $\eta_{k}^{s}$ according to (12)

22

23 Update the iterate by  $w_{k+1}^{s}=w_{k}^{s}-\eta_{k}^{s}v_{k}^{s}$

24 end for

25 $\widetilde{w}_{s}=w_{m}^{s}$

26 end for

Algorithm 1 MB-SARAH-RHBB/RHBB+

Input: $\widetilde{w}_{0}$, update frequency $m$, batch sizes $b,b_{1},b_{2}$,
constant sequence $\\{\tilde{\eta}_{0}^{s}\\}$, balance parameter
$\gamma_{2}\geq 1$, hedge base $\alpha$ and monotone decreasing function $h$,
probability distribution $Q$.

Output: approximate solution $\widetilde{w}_{s}$.

1

2 Outer Loop: for  _$s=1,2,...,$_ do

3 $\widetilde{w}=\widetilde{w}_{s-1}$

4 $w_{0}^{s}=\widetilde{w}$

5 $\varphi=\frac{1}{n}\sum_{i\in\Omega}\nabla
f_{i}\left(\widetilde{w}\right)=\nabla P\left(\widetilde{w}\right)$

6 $\tilde{v}_{0}^{s}=\varphi$

7 $w_{1}^{s}=w_{0}^{s}-\tilde{\eta}_{0}^{s}\tilde{v}_{0}^{s}$

8 Inner Loop: for  _$k=1,2,...,m-1$_ do

9 Pick subset $S\subset\\{1,\ldots,n\\}$ of size $b$ uniformly at random

10 Update the semi-stochastic gradient by  $\tilde{v}_{k}^{s}=\nabla
P_{S}\left(w_{k}^{s}\right)-\nabla P_{S}(\widetilde{w})+\varphi$

11 Compute the step size $\tilde{\eta}_{k}^{s}$ by (1) Option I or (2) Option
II:

12 (1) Option I : RHBB rule

13 Configure $Q$ as uniform probability distribution

14 Pick subset $S_{1}\subset\\{1,\ldots,n\\}$ of size $b_{1}$ randomly
according to $Q$

15 Pick subset $S_{2}\subset\\{1,\ldots,n\\}$ of size $b_{2}$ randomly
according to $Q$

16 Calculate $\tilde{\eta}_{k}^{s}$ according to (9)

17

18 (2) Option II : RHBB+ rule

19 Configure $Q$ to our needs

20 Pick subset $S_{1}\subset\\{1,\ldots,n\\}$ of size $b_{1}$ randomly
according to $Q$

21 Pick subset $S_{2}\subset\\{1,\ldots,n\\}$ of size $b_{2}$ randomly
according to $Q$

22 Compute $\nabla P_{S_{1}}^{+},\nabla P_{S_{2}}^{+}$ according to (11)

23 Calculate $\tilde{\eta}_{k}^{s}$ according to (13)

24

25 Update the iterate:
$w_{k+1}^{s}=w_{k}^{s}-\tilde{\eta}_{k}^{s}\tilde{v_{k}}^{s}$

26 end for

27 $\widetilde{w}_{s}=w_{m}^{s}$

28 end for

Algorithm 2 mS2GD-RHBB/RHBB+

###### Remark 1

At beginning of each epoch $s$, the constant step sizes $\eta_{0}^{s}$ and
$\tilde{\eta}_{0}^{s}$ are used in the first deterministic step of full pass.
Our RHBB or RHBB+ are placed in following stochastic stages to match
stochastic recursive or semi-stochastic estimators and create smooth paths for
the convergence. Distribution $Q$ can be tailored to the needs of particular
data sets.

##  5 Convergence Analysis

Hereafter, we use following notations: let the batch correction
$\overline{b}=\max\\{|S_{1}|,|S_{2}|\\}=\max\\{b_{1},b_{2}\\}$. Under the
probability distribution $Q\sim\left\\{q_{1},q_{2},\ldots q_{n}\right\\}$, let

$L_{q}=\max_{i}\frac{L}{n\cdot q_{i}},\quad
L_{r}=\frac{L}{L_{q}},\quad\mu_{q}=\min_{i}\frac{\mu}{n\cdot
q_{i}},\quad\mu_{r}=\frac{\mu_{q}}{\mu}.$ (16)

Then, we have straightforward results: $L_{q}\geq L$, $\mu_{q}\leq\mu$ and
$L_{r}\leq 1$, $\mu_{r}\leq 1$. Besides, we obtain an approximate condition
number
$\kappa^{+}=\frac{L_{q}}{\mu_{q}}=\frac{\kappa}{L_{r}\mu_{r}}\geq\kappa$.

Since the adaptor $h(\cdot)$ is monotone decreasing, we need the following
uniform boundness assumption in the convergence analysis.

Assumption 3 (Uniform boundness). The iterative adaptor $h(\cdot)$ is
continuous over the bounded closed domain, monotone decreasing with respect to
the epoch count $s$ and the inner count $k$, i.e., there exist constants
$\hat{\alpha}$, $\tilde{\alpha}$ such that

$1<\tilde{\alpha}<\alpha^{h(\sigma_{1}s+\sigma_{2}k)}<\hat{\alpha}.$

###  5.1 MB-SARAH-RHBB and MB-SARAH-RHBB+

To begin with, we provide subsequent Lemma 1 to show the summative boundary of
$\mathbb{E}\left[\left\|\nabla P\left(w\right)\right\|^{2}\right]$ within the
$s$-th epoch (inner loop).

###### Lemma 1

Suppose that Assumption 1, 2a and 3 hold. The subsets $S$, $S_{1}$, $S_{2}$
are selected uniformly at random of size $b,b_{1},b_{2}$ respectively. Then,
for any $s\geq 1$ in MB-SARAH-RHBB, we have

$\displaystyle\sum_{k=0}^{m}\mathbb{E}\left[\left\|\nabla
P\left(w_{k}^{s}\right)\right\|^{2}\right]\leq\frac{2\mu\overline{b}L}{\hat{\alpha}\gamma
L+(1-\tilde{\alpha})\gamma\mu}\mathbb{E}\left[P\left(w_{0}^{s}\right)-P\left(w_{*}\right)\right]$
$\displaystyle+\sum_{k=0}^{m}\mathbb{E}\left[\left\|\nabla
P\left(w_{k}^{s}\right)-v_{k}^{s}\right\|^{2}\right]-\left(1-\frac{\hat{\alpha}\gamma
L^{2}+(1-\tilde{\alpha})\gamma
L\mu}{\mu\overline{b}L}\right)\sum_{k=0}^{m}\mathbb{E}\left[\left\|v_{k}^{s}\right\|^{2}\right].$

Furthermore, if Assumption 2b holds, and the subsets $S_{1}$ and $S_{2}$ are
sampled according to the probability distribution $Q$ of size $b_{1}$ and
$b_{2}$. For any $s\geq 1$ in MB-SARAH-RHBB+, we have

$\displaystyle\sum_{k=0}^{m}\mathbb{E}\left[\left\|\nabla
P\left(w_{k}^{s}\right)\right\|^{2}\right]\leq\frac{2\mu_{q}\overline{b}L_{q}}{\hat{\alpha}\gamma
L_{q}+(1-\tilde{\alpha})\gamma\mu_{q}}\mathbb{E}\left[P\left(w_{0}^{s}\right)-P\left(w_{*}\right)\right]$
$\displaystyle+\sum_{k=0}^{m}\mathbb{E}\left[\left\|\nabla
P\left(w_{k}^{s}\right)-v_{k}^{s}\right\|^{2}\right]-\left(1-\frac{\hat{\alpha}\gamma
LL_{q}+(1-\tilde{\alpha})\gamma
L\mu_{q}}{\mu_{q}\overline{b}L_{q}}\right)\sum_{k=0}^{m}\mathbb{E}\left[\left\|v_{k}^{s}\right\|^{2}\right].$

###### Proof

Deferred to the Appendix A.

Next, we prove that the deviation (expected distance) of the full gradient to
the recursive estimates is upper bounded within the $s$-th epoch (inner loop).

###### Lemma 2

Suppose that Assumption 1, 2a hold. The subsets $S$, $S_{1}$, $S_{2}$ are
selected uniformly at random of size $b,b_{1},b_{2}$, respectively. Within the
$s$-th epoch of MB-SARAH-RHBB, for any $1\leq k\leq m$, we have

$\mathbb{E}\left[\left\|\nabla
P\left(w_{k}^{s}\right)-v_{k}^{s}\right\|^{2}\right]\leq\frac{n-b}{b\left(n-1\right)}\left(\frac{\hat{\alpha}\gamma
L^{2}+\left(1-\tilde{\alpha}\right)\gamma\mu L}{\overline{b}\mu
L}\right)^{2}\sum_{j=1}^{k}\mathbb{E}\left[\left\|v_{j-1}^{s}\right\|^{2}\right].$
(17)

If Assumption 2b holds further and the subsets $S_{1}$, $S_{2}$ are sampled
according to the probability distribution $Q$. Within the $s$-th epoch of MB-
SARAH-RHBB+, for any $1\leq k\leq m$, we then have

$\mathbb{E}\left[\left\|\nabla
P\left(w_{k}^{s}\right)-v_{k}^{s}\right\|^{2}\right]\leq\frac{L_{r}^{2}\left(n-b\right)}{b\left(n-1\right)}\left(\frac{\hat{\alpha}\gamma
L_{q}+(1-\tilde{\alpha})\gamma\mu_{q}}{\overline{b}\mu_{q}}\right)^{2}\sum_{j=1}^{k}\mathbb{E}\left[\left\|v_{j-1}^{s}\right\|^{2}\right].$