Reconstructing nonlinearities with intermodulation spectroscopy

We describe a method of analysis which allows for reconstructing the nonlinear disturbance of a high Q harmonic oscillator. When the oscillator is driven with two or more frequencies, the nonlinearity causes intermodulation of the drives, resulting in a complicated spectral response. Analysis of this spectrum allows one to approximate the nonlinearity. The method, which is generally applicable to measurements based on resonant detection, increases the information content of the measurement without requiring large detection bandwidth, and optimally uses the enhanced sensitivity near resonance to extract information and minimize error due to detector noise.

A harmonic oscillator with high quality factor Q, i.e., large stored energy in relation to the energy lost in each cycle, is a very useful tool for precision measurement. When the oscillator interacts with an object of interest, the first order correction to the linear response is a shift of the resonant frequency ω 0 , which results in a large change of the amplitude and phase of the response near resonance. Resonant detection is employed in a wide variety of measurements, and particularly in nanoscience, where examples include the radio frequency single electron transistor [1], dispersive readout of superconducting qubits [2], measuring the deflection of a nanomechanical oscillator [3], and nanoscale imaging of surfaces with dynamic atomic force microscopy (AFM) [4]. In this letter we describe a method of resonant detection which exploits the nonlinear phenomenon of intermodulation, to generate a spectral response near resonance that enhances the information content of the measurement. We provide a general theoretical framework for analyzing the intermodulation spectrum and for reconstructing the nonlinearity.
The sensitivity of resonant detection can be understood from the large transfer gain, G(ω), of a high-Q oscillator. In the narrow frequency band B = ω 0 /Q centered at resonance, the oscillator gives large response to a small stimulus. When the linear oscillator interacts with the object of interest, the equation of motion becomes nonlinear and higher harmonics of the drive frequency will appear in the response. One could in principle analyze the amplitude and phase of these higher harmonics to extract the nonlinearity. The problem with this approach is: (i) it requires a very large detection bandwidth, B D = nω 0 to get n harmonics and (ii) the transfer gain for harmonics is far less than one, falling off as 1/ω 2 above resonance, so the response at higher harmonics is often buried in the detector noise. In the field of AFM, the latter problem was addressed by coupling more oscillators, utilizing additional torsional [5] and flexural [6,7] eigenmodes of the cantilever. In this context, different types of multi-tone excitation schemes have been put forward as a means of extracting qualitative information about the nonlinearity [8][9][10].
An alternative approach, sketched in Fig. 1, is to drive the oscillator with two pure tones ω 1 and ω 2 , so as to produce intermodulation products (IMPs) or "mixing products" which occur at the frequencies nω 1 + mω 2 (n and m integers), where the order of the IMP is defined by |n| + |m|. The drive frequencies can be chosen in a variety of ways, with the objective of creating many IMPs near resonance where the sensitivity is enhanced. These IMPs measured near resonance represent a partial spectral response of the system from which one would like to reconstruct the nonlinearity.
Our work is motivated by the desire to understand the information content of intermodulation spectra, as measured in superconducting microresonators at GHz frequencies [11,12], and in AFM cantilevers oscillating at several hundred kHz while interacting with a surface [13,14] (see Fig. 2). However, we wish to stress that the reconstruction algorithm described here can be applied to a wide variety of experiments which exploit resonant detection. To test the reconstruction algorithm, we use simulated data obtained by numerical integration of the equation of motion of the nonlinear oscillator, where, in . Measured intermodulation spectra taken from a nonlinear superconducting resonator [12] (left panel) and an AFM cantilever interacting with a surface [14] (right panel). In these experiments the two central peaks with largest response are the drive tones. The amplitude in dB is normalized to the maximum response.
contrast with experiment, we can work backward from the simulated intermodulation spectrum to the known expansion coefficients of the nonlinearity, thus allowing a more systematic study of errors.
The equation of motion of a single-mode harmonic oscillator, which is driven by two pure tones, subject to viscous damping, and perturbed by a nonlinearity, reads: Here, Ω 1,2 = ω 1,2 /ω 0 , dots denote derivatives with respect to the dimensionless time τ = ω 0 t, and the dimensionless coordinate ζ is the deviation of the oscillator from the equilibrium position, suitably normalized. The dimensionless nonlinear "force" F nl (ζ) is assumed conservative, with no explicit time-dependence. While the motion of a damped, driven nonlinear system of the type (1) can in general be either regular or chaotic, we restrict ourselves here to nonlinearities and drive strengths which are weak enough, such that for commensurate drive frequencies with greatest common divisor ∆Ω ≡ gcd(Ω 1 , Ω 2 ) the motion in the steady state is periodic in 2π/∆Ω. We may therefore expand the response in a discrete Fourier series, ζ(τ ) = k ζ k e ik ∆Ω τ , where the expansion coefficients ζ k are the complex numbers (amplitude and phase) which constitute the intermodulation spectrum. With Fourier expansions of the nonlinear force F nl (ζ(τ )) = k F nl,k e ik ∆Ω τ , and the drive force, where G k = 1/(1 − k 2 ∆Ω 2 + ik∆Ω/Q). Note that while (2) is an exact equation, it is not a solution of (1), since F nl,k depends on all the ζ k in a highly nonlinear way.
To proceed, we represent the nonlinear force as a polynomial in ζ(τ ), F nl (ζ(τ )) = ∞ j=1 g j ζ j (τ ), which allows us to express the Fourier components of the nonlinear force F nl,k as a linear combination of the polynomial coefficients, g j : where the matrix element H kj is the kth spectral component of the Fourier transform of ζ j (τ ), alternatively written as H k1 = ζ k and H k, For a known nonlinear force, we immediately find a first approximation ζ (1) k of the spectrum, by evaluating F nl,k in (2) at the free solution ζ In the next step we use ζ (1) k in the same way to obtain ζ (2) k , etc. Each iteration generates new spectral response, but only at the intermodulation frequencies. If this procedure converges to a solution ζ k , then our assumption of periodic orbits is certainly justified. However, here we are interested in the inverse problem of finding the nonlinear force, starting from ζ exp k , which is measured with limited bandwidth in the presence of noise. Inspection of Eqs. (3) and (2) reveals that one must invert the matrix H to find the expansion coefficients g j of the nonlinearity: The inverse problem is generally complicated by the fact that the measured spectrum ζ exp k has only a finite number N p of peaks that can be observed above the detector noise level, notably those close to the resonance frequency. Nevertheless, with a proper intermodulation drive scheme one can achieve enough intermodulation peaks with good signal-to-noise ratio to construct a matrix H, and calculate the coefficients g j . We use the following algorithm: 1. Choose the drive frequencies so that the frequency spacing of IMPs, ∆Ω, gives N p peaks in the measured intermodulation spectrum ζ exp k at frequencies k∆Ω for k ∈ {k 1 , . . . , k Np }, where N p is at least twice the number of expansion coefficients j max desired in approximating the nonlinearity. 2. Define the intermodulation response ζ k to be ζ k = ζ exp k for k ∈ {k 1 , . . . , k Np } and ζ k = 0 otherwise. 3. With the intermodulation response ζ k , calculate the matrix H, restricted to j ≤ j max and k ≤ k max · j max , where k max = max{k 1 , . . . , k Np }. 4. Calculate the j max expansion coefficients g j (restricted to be real) using Eq. (4), with the inverse replaced by a pseudo inverse, the sum restricted to k ∈ {k 1 , . . . , k Np }, and H restricted to k ∈ {k 1 , . . . , k Np } and to j ≤ j max . The pseudo inverse gives a least square fit of the j max coefficients g j to the N p equations ζ k = ζ (0) k + G k jmax j=1 H kj g j , thus reducing sensitivity to both measurement errors and systematic errors in the calculation of H.
We have tested this reconstruction algorithm using realistic parameters for two different experimental realizations. In each case, the response of the nonlinear oscillator was simulated by numerical integration of Eq.(1) using the Fortran solver DDASKR [15]. The output of this integrator was sampled appropriately and fast Fourier transformed to generate an intermodulation spectrum. Tests were made to ensure that there was no Fourier leakage of the spectral peaks, so that the background level in the spectrum corresponded to the error tolerance set in the integrator. Detector noise was simulated by adding random gaussian noise ∆ζ(τ ) to each time sample, before the Fourier transform.
As a first example we consider a superconducting coplanar wave-guide resonator with a Josephson tunnel junction or weak link, currently studied for low-noise amplification and the readout of quantum bits [11,[16][17][18][19]. The classical dynamics of the resonator can be mapped to our model equation by using ζ = I/I 0 , where I is the current in the center strip and I 0 the critical current. The quality factor Q can be engineered in a wide range, up to 10 6 . Typically a nonlinear kinetic inductance provides a nonlinear "force" which is constrained by g j = 0 for even j, because the inductance dF nl /dζ should not depend on the sign of the current ζ. We can thus employ a drive scheme which is only sensitive to IMPs of odd order, with drive frequencies closely spaced near resonance, ∆Ω ≡ Ω 2 − Ω 1 ≪ 1. This drive scheme was used for the experimental spectra shown in Fig. 2. For weak drive, spectra are well described by the lowest nonlinearity of Kerr-type, g 3 , but for strong drive, higher nonlinear coefficients (g 5 ,g 7 , ...) play a role.
In Tab. I we give the percent error in the coefficients g recon j which were reconstructed from data that was simulated using three coefficients g 3 , g 5 , g 7 . Results are shown for two different simulated noise levels (average of 1000 simulations) and two different quality factors. In the absence of noise, ∆ζ = 0, the small errors in the reconstruction are explained by the fact that we only used 20 closely spaced spectral peaks around the resonance, leading to a slight miscalculation of the matrix H. These systematic errors, which are typically small, can be reduced by using an improved reconstruction algorithm which estimates the unmeasured spectral peaks to achieve a selfconsistent solution of Eq.(2). In the presence of noise, reconstruction from the intermodulation spectrum gives excellent results for all coefficients. For comparison, we applied the reconstruction algorithm to a simulated spectrum of the same number of harmonics, generated with a single drive having strength F d = F d1 +F d2 , for the same parameters given in Tab. I. In the presence of noise, the spectrum of harmonics was unable to reconstruct the coefficients g 5 and g 7 . This nicely illustrates the advantage of resonant detection in an intermodulation scheme.
dynamic atomic force microscopy (AFM), the resonator is a cantilever oscillating about its equilibrium position with dimensionless amplitude ζ = z/z s , effective mass m, and stiffness k = m ω 2 0 . The nonlinearity is provided by the tip-surface force, F nl = F ts (ζ), which we wish to determine. Intermodulation spectroscopy can be performed while scanning over a surface, measuring the dominant IMPs at each image point [13,14]. The inversion algorithm described here can be used to rapidly determine the force-distance curve at each image pixel, while scanning at normal speeds for dynamic AFM, one of the major objectives of AFM development [5,20].
The nonlinear force typically encountered in AFM experiments has an attractive region close to a repulsive surface. We restrict our study to conservative forces, derivable from a potential function, F ts = −dU/dz, and desire to reconstruct a polynomial approximation of the dimensionless force F ts (ζ) ≡ F ts (z + z 0 )/(z s k), where z 0 is the probe height above the surface in the absence of drive signals, an experimentally tunable parameter. For a Morse potential, this force has the functional form F ts = 2α(e −2(ζ−σ)/λ − e −(ζ−σ)/λ ), which is described by the location of the potential minimum, σ, its width, λ, and an overall scale factor, α [21]. We used this nonlinear force to simulate data with a quality factor Q = 50 and two drives centered around resonance. Because this drive scheme generates near resonance only IMPs of odd order, intermodulation peaks around twice the resonance frequency were needed to gain information about the coefficients g j with even j. 20 peaks around resonance and 10 peaks around 2ω 0 were taken for the reconstruction. All other spectral information was set to zero.
The reconstruction algorithm was run using the amplitude and phase of these N p = 30 peaks, and j max = 15 coefficients g j were taken for the polynomial representation of the nonlinear force. The result of the reconstruction is shown in Fig. 3, where we see excellent correspondence between the original Morse force curve and the reconstructed polynomial.
To make the inversion more realistic, we added Gaussian noise corresponding to a maximum-spectral- amplitude-to-noise ratio ∼ 10 −4 (Fig. 3, inset), a dynamic range typical of the optical lever detectors in most AFMs. With this noise we see that the polynomial missrepresents the force at maximum distance from the surface, where small changes in the coefficients g j with large j have a large effect. We also note that the polynomial may in general miss-represent the force near the point of maximum curvature, as many coefficients are needed to reconstruct sharp kinks in the force-distance curve. To improve upon these results, we can carry the reconstruction one step further. Using prior knowledge of the interaction in a particular experiment, for example that F ts goes asymptotically to zero for large ζ and that a sharp kink exists at the contact point, one can make a weighted fit of the reconstructed polynomial to a particular force model containing far fewer than j max parameters.
Thus it is possible to extract a very nonlinear forcedistance curve in the fast scanning mode without using high frequency components of the cantilever response spectrum. Intermodulation spectroscopy allows one to collect information in the form of coefficients g j with large j, by measuring high order IMPs at much lower frequencies than the j th harmonic, where such information would occur with a single drive tone. The ability of two drives to down-convert this information to lower frequencies, where it can be put near resonance and acquired with good transfer gain, is the major strength of the intermodulation spectroscopy technique. The method will be accurate if the dominant part of the spectral information is contained in the peaks which are used to construct the matrix H, and if the traversed force curve is well approximated by j max coefficients. Note that both amplitude and phase data are necessary to calculate the matrix H.
We briefly note that we also successfully performed a simulated reconstruction of force-distance curves with two other drive schemes: One with a low-frequency drive and a drive at resonance, as sketched in Fig. 1, and another scheme with Ω 1 = 1/3 and Ω 2 = 2/3 + ∆Ω. These drive schemes have the advantage that they produce IMPs of both odd and even order near resonance, and are therefore better with the highest Q oscillators, and for nonlinearities containing both odd and even coefficients.
In summary, intermodulation spectroscopy is a measurement technique which employs a high Q oscillator and two drive frequencies to extract response caused by a nonlinear disturbance. The drive frequencies are chosen commensurate with a greatest common divisor ∆Ω, and placed so as to produce many IMPs near resonance, where large transfer gain allows for detection with good signal-to-noise ratio. A reconstruction algorithm was derived and tested against simulated data including noise. Excellent reconstruction of the nonlinearity was demonstrated with relatively few spectral peaks obtained at low frequencies. Together with the reconstruction algorithm, intermodulation spectroscopy should find wide-ranging use with the many types of sensors that are based on disturbance of a resonance.
We acknowlegde helpful discussions with Jack Harris, Arvind Raman, and Ricardo Garcia. This work was supported by the Swedish Research Council and The Foundation for Strategic Research.