## 4 Parametric Amplification of the Thermal Fluctuations  
### Theoretical Model  
We consider a parametric modulation of the spring constant of a normal mode of the form \(k=k_{0}[1+\xi\sin(\omega_{p}t)]\), where \(k_{0}\) is the spring constant in absence of modulation, \(\xi\) the normalized modulation amplitude, which is proportional to the applied modulation voltage amplitude and \(\omega_{p}\) the pump modulation frequency. The classical dynamics of the amplitude around equilibrium, \(x(t)\), can be modeled by the following differential equation [60]  
\[\ddot{x}+\Gamma\dot{x}+\omega_{m}^{2}[1+\xi\sin(\omega_{p}t)]x=F_{\text{th}}/m, \tag{2}\]  
where \(\omega_{m}=\sqrt{k_{0}/m}\) is the bare mechanical resonance frequency of the normal mode considered, \(m\) the effective mass of the mode, \(\Gamma=\omega_{m}/Q\) its damping rate and \(F_{\text{th}}\) a thermal noise force arising from the coupling with the thermal environment. For a high quality factor oscillator (\(\omega_{m}\gg\Gamma\)) and a pump modulation frequency \(\omega_{p}=2(\omega_{m}+\delta)\) close to the second harmonic frequency, it is convenient to introduce the slowly varying envelope \(A(t)\) defined by  
\[x(t) =\left[A(t)e^{-i\omega_{m}t}+A^{*}(t)e^{i\omega_{m}t}\right]/2 \tag{3}\] \[=X_{1}(t)\cos(\omega_{m}t)+X_{2}(t)\sin(\omega_{m}t) \tag{4}\]  
where the quadratures \(X_{1}\) and \(X_{2}\) are respectively given by the real and imaginary parts of \(A\). Under the rotating wave approximation, Equation (2) yields  
\[\dot{A}=-\left(\frac{\Gamma}{2}+i\delta\right)A+\omega_{m}\frac{\xi}{4}A^{*}+ \tilde{F}, \tag{5}\]  
where \(\tilde{F}=F_{\text{th}}e^{i\omega_{m}t}/m\omega_{m}\). Fourier transforming Equation (5) using the convention \(f(\omega)=\int f(t)e^{-i\omega t}dt\) yields  
\[\left(\gamma+i\omega+i\delta\right)A(\omega)-\gamma eA^{*}(\omega)=\tilde{F} (\omega), \tag{6}\]  
where  
\[\gamma=\frac{\Gamma}{2}\qquad\text{ and }\qquad\epsilon=\frac{\omega_{m} \xi}{4\gamma}=\frac{Q\xi}{2}. \tag{7}\]  
Figure 4: (**a**) Measured (1,1) mode resonance frequencies of the A (red circles) and B (red dots) membranes in a double-membrane array as in Figure 1d, as a function of \(V_{dc}\). (**b**) Measured (2,2) mode resonance frequencies of the A (empty blue squares) and B (filled blue squares) membranes, as a function of \(V_{dc}\). From Equation (6) and its complex conjugate one readily obtains  
\[A(\omega) =\frac{(\gamma+i\omega-i\delta)\bar{F}(\omega)+\gamma\epsilon\bar{F }^{*}(\omega)}{(\gamma+i\omega+i\delta)(\gamma+i\omega-i\delta)-(\gamma\epsilon) ^{2}}, \tag{8}\] \[A^{*}(\omega) =\frac{(\gamma+i\omega+i\delta)\bar{F}^{*}(\omega)+\gamma\epsilon \bar{F}(\omega)}{(\gamma+i\omega+i\delta)(\gamma+i\omega-i\delta)-(\gamma \epsilon)^{2}}. \tag{9}\]  
The noise spectrum, \(S_{X_{\theta}}(\omega)=\langle X_{\theta}(\omega)X_{\theta}(-\omega)\rangle\), of an arbitrary quadrature \(X_{\theta}=X_{1}\cos(\theta)+X_{2}\sin(\theta)\) can be computed using the fact that  
\[\langle F_{\rm th}(\omega)F_{\rm th}(\omega^{\prime})\rangle=4k_{B}Tm\Gamma \delta(\omega+\omega^{\prime}), \tag{10}\]  
such that only \(\langle\bar{F}(\omega)\bar{F}^{*}(-\omega)\rangle\) and \(\langle\bar{F}^{*}(\omega)\bar{F}(-\omega)\rangle\) are nonzero and both equal to \(4k_{B}T\gamma/m\omega_{m}^{2}\), where \(k_{B}\) is the Boltzmann constant and \(T\) the thermal bath temperature.