Studying Soft Interfaces with Shear Waves: Principles and Applications of the Quartz Crystal Microbalance (QCM)

The response of the quartz crystal microbalance (QCM, also: QCM-D for “QCM with Dissipation monitoring”) to loading with a diverse set of samples is reviewed in a consistent frame. After a brief introduction to the advanced QCMs, the governing equation (the small-load approximation) is derived. Planar films and adsorbates are modeled based on the acoustic multilayer formalism. In liquid environments, viscoelastic spectroscopy and high-frequency rheology are possible, even on layers with a thickness in the monolayer range. For particulate samples, the contact stiffness can be derived. Because the stress at the contact is large, the force is not always proportional to the displacement. Nonlinear effects are observed, leading to a dependence of the resonance frequency and the resonance bandwidth on the amplitude of oscillation. Partial slip, in particular, can be studied in detail. Advanced topics include structured samples and the extension of the small-load approximation to its tensorial version.


Introduction
The use of the quartz crystal microbalance (QCM) in diverse areas keeps growing. Among its advantages is its simplicity. It is not difficult to mount the resonator plate in one way or another. Electrical interrogation of the resonance parameters is not difficult either. Simplicity entails versatility. The QCM can easily be combined with electrochemistry, optical reflectometry, scanning force microscopy, and other instruments of interface analysis.
Reviews on the QCM can be found in [1][2][3][4][5][6]. The QCM becomes slightly more complicated on second glance. That concerns intricacies in the operation as well as pitfalls in interpretation. This text is meant to give practitioners a quick start, still going to beyond the simple Sauerbrey picture. The concepts underlying the more advanced models are important when carrying the QCM to non-standard applications (such as the freezing of droplets or the impact of spheres).
The following list summarizes a few applications. (The list anticipates the later text in so far, as the topics and problems addressed have not been explained yet, but will be.) The numbers of citations are returns from a search in the Web of Science, where the dates were limited to January 2019April 2021. A total of 1067 entries contain the keyword "QCM". "QCM-D" returns 419 entries. "QCM-D" here denotes all instruments reporting frequency as well as bandwidth (or, equivalently, the dissipation factor) on a number of different overtones. In the majority of these publications, the QCM is one out of a few different instruments employed to study the respective samples. The QCM in these studies is a routine device.
Topics of particular relevance with regard to modeling are:  Numerous publications discuss the mass uptake of nanoporous and other rigid layers when exposed to a vapor of the analyte [7]. The porous layer takes the role of the receptor. The limit of detection of the QCM easily suffices for sensing building on this principle. (It does not easily suffice for similar sensors, building on adsorption to a planar surface.) These rigid structures swell and soften less than the polymer films, which took a similar role in the past [8]. While the emphasis in these works is on gravimetry, an analysis taking viscoelasticity into account (Equation (46)) will provide for more in-depth information. Also, it will yield a more accurate value for the mass uptake than the Sauerbrey equation.  The search term "EQCM" returns 137 citations. These are increasingly concerned with an analysis beyond gravimetry. The non-gravimetric effects in this context mostly originate from roughness (Equation (77)), from the viscoelasticity of the double layer (Equation (59)), and from the softness of an active polymer layer (if present, Equation (46)).  The keyword "QCM-D and brush" returns 37 entries. The brushes often undergo swelling/deswelling transitions or show electroresponsivity. Brushes should be modeled taking viscoelastic effects into account. The shear modulus varies between the bottom and the top, which necessitates the use of a viscoelastic profile (Equation (60)).  114 publications are returned for "QCM and particles". The interpretation of such QCM data is a topic of ongoing research (Sections 6.1 and 6.3). For instance, the amount of liquid mass contributing to the gravimetric signal (the "trapped mass") usually is not known, quantitatively. A few publications explicitly refer to the positive frequency shift induced by sufficiently large particles (Equation (70)).  69 publications mention bacteria, which often implies bacterial adsorption (reviewed in [9]). In these cases (and also for cell cultures and biofilms) the shear wave often does not reach to the top of the layer. The QCM then cannot measure the thickness.
If such a thick sample is homogeneous in viscoelastic terms, the QCM reports the shear modulus of this medium (Equation (31)). For reviews of applications in the life sciences, in general, see [10] and [11].  Interestingly, 247 hits are returned when asking for "QCM and protein". Protein adsorption is also routinely and successfully probed with optical techniques such as surface plasmon resonance (SPR) spectroscopy. The added information contained in the layer's viscoelasticity (Equation (52)) is a distinctive advantage of the QCM.  On the conceptual side, high-frequency rheology on polymers receives considerable attention. Some of these publications are returned when the keyword is "tribology" or "viscoelasticity". A recent review is contained in [12]. The equations applied for analysis mostly are similar to what is described in Sections 4.5.3. and 4.6.2. For thick films (microns), Equation (41) is a suitable fit function. Because the frequency shifts are large, temperature effects are irrelevant. For thin films (tens of nanometers), Equation (46) is more suitable than Equation (41). One can hope for data from more than 10 overtones available for analysis. However, the frequency shifts are smaller, which makes the analysis more susceptible to artifacts, for instance caused by changes of temperature.  Rather few (<10) publications mention large amplitudes and nonlinear behavior. While this is an interesting field in the authors' opinion, it has not been explored much, yet. Some recent reviews (such as [11]) cover acoustic sensors other than the QCM (SAW devices, FBARs, nanoresonators, …). These devices operate at higher frequency than the QCM, which implies improved sensitivity in gravimetry. Most of them are smaller than the QCM. Some models of the QCM can be transferred to the smaller devices, but vibration in a sufficiently clean thickness-shear mode is often in question.
The quartz crystal microbalance is about 60 years old. At that time, people working with quartz resonators knew that one can lower the frequency by scribbling onto the plate with a pencil. Günter Sauerbrey turned this practice into an equation and an instrument [13]. The frequency decreases, following the relation ω 0 ≈ (κ R /m R ) 1/2 with κ R some effective stiffness and m R some effective mass. By making the crystal heavier, one slows down its resonant vibration ( Figure 1).

Figure 1. (A)
The simplest possible harmonic resonator. A mass is linked to a rigid wall across a spring. The resonance frequency is ω 0 ≈ (κ R /m R ) 1/2 . If m R increases, the resonance frequency decreases in consequence. (B) A resonator abandoning the rigid wall. Two masses are linked across a spring. The resonance frequency is given as ω 0 ≈ (κ R /µ) 1/2 , where µ is the reduced mass (µ = m 1 m 2 /(m 1 + m 2 )). Again, increasing one of the two masses will lower the resonance frequency. The diagram in (B) contains discrete elements, similar to the equivalent lumped-element circuits discussed in Section 8.4. (C) Contrasting to the resonators in (A,B), the bell does not consist of discrete masses and springs. It is an elastic body with a certain shape, made from materials with a certain density and stiffness. Finding its resonance frequencies (plural) is a classical problem of acoustics. For any given resonance, one may construct an equivalent lumped-element model containing discrete elements (as in (A) or (B)), which reproduces this one resonance. In principle, one might tune the bell by gluing weights to its rim. The common practice rather is to remove metal in annular rings, usually from the inside [14]. That changes both the effective mass and the effective spring constant. The focus then is not usually on the absolute frequency of the fundamental mode, but rather on the ratios between the overtone frequencies and the fundamental frequency. These ratios govern the perception of the bell's sound.
Sauerbrey realized that the relation between mass and frequency shift is particularly simple if, firstly, the resonator is a plate oscillating in the thickness-shear mode and if, secondly, the sample consists of a thin film. Most resonators at this time indeed were thickness-shear resonators. Thin films were routinely coated onto these as electrodes. This insight allowed to develop film-thickness monitors for deposition processes of various kinds.
A side remark: Plates of α-quartz may resonate in the thickness-shear mode if the crystallographic x-axis is in the surface plane. There are certain angles between the crystallographic y-axis and the surface normal, at which the temperature-frequency coupling at room temperature almost vanishes. One of these temperature-compensated cuts is the widely used AT cut.
Sauerbrey's famous formula is: Δf/f ref is the fractional frequency shift. m f and m q are the mass per unit area of the film and the resonator, respectively. One might also talk about "mass" rather than "mass per unit area", but the latter term is more practical. For instance, the mass per unit area is easily converted to thickness, if the density is known. The resonator's mass per unit area, m q , may be replaced by Z q /(2nf 0 ), where Z q = 8.8 × 10 6 kg/(m 2 s) is the resonator's shearwave impedance, f 0 is the frequency of the fundamental (often 5 MHz), and n is the overtone order. These relations inserted into Equation (1) lead to: The reasoning behind the Sauerbrey equation is sketched in Figure 2. The interesting vibration modes are standing transverse waves with antinodes at the surfaces. Figure 2 shows the displacement patterns of the fundamental mode and the 3rd overtone as dashed lines. (An "overtone" here is an eigenmode, that is, a solution to the boundary value problem, not to be confused with second-harmonic generation or third-harmonic generation. The latter terminology is also common in acoustics.) The wavelength is 2d q /n with d q the thickness of the plate and n the number of nodal planes. Only odd overtones (n = 1, 3, 5, …) can be excited piezoelectrically, because the even overtones lead to a surface charge with the same sign on both sides. It is customary to label the overtones (the "overtone order") with indices equal to the number of nodal planes. The frequency is f res = nc q /(2d q ) with c q the speed of sound.

Figure 2. (A)
A film increases the wavelength of the standing wave, thereby decreasing the resonance frequency. (B) If the film is softer than the plate, the displacement pattern has a kink at the interface. In this case, the fractional frequency shift is proportional to the fractional increase in mass (rather than thickness).
The Sauerbrey mass, m f , is often quoted in units of µg/cm 2 . More intuitive would be a thickness (the "Sauerbrey thickness"), the calculation of which, however, requires knowledge of the density. With 5 MHz crystals and a density of 1 g/cm 3 , a 1 nm film shifts the frequency by Δf/n = −5.7 Hz. With this density, 1 µg/cm 2 corresponds to 10 nm.
Consider a film with the exact same acoustic properties as the resonator itself (top in Figure 2B). In acoustic terms, the film makes the plate thicker. The wavelength increases correspondingly and the frequency decreases. If the film is much thinner than the plate, one might expect the relation Δf f ref ≈ -d f d q ⁄ . This is not the Sauerbrey equation. The Sauerbrey equation makes a statement about mass, not about thickness. Above, the film was assumed to have the same acoustic properties as the plate (same density, ρ, same shear modulus, G). If this is not the case, the displacement pattern has a kink at the resonatorfilm interface (bottom in Figure 2B). Sauerbrey went through the mathematics and realized that the fractional thickness must be replaced by the fractional mass in order to let the relation be applicable to arbitrary materials.
Assume a resonator with a fundamental frequency of 5 MHz. Further, assume that the frequency can be determined with a precision of 0.5 Hz. This precision in frequency translates to a precision in mass of about 10 ng/cm 2 . Because the unit has the prefix "nano", people have advocated the term "quartz crystal nanobalance". In the end, the "quartz crystal microbalance" became the accepted term. With a density of about 1 g/cm 3 , the thickness resolution of this QCM is about 0.1 nm. The QCM has "submonolayer sensitivity". Monomolecular layers of typical bio-adsorbates are slightly thicker than 1 nm. The QCM has submonolayer sensitivity, but the sensitivity is not deep in the submonolayer range. In order to study the kinetics of adsorption in detail, one would wish for an even better limit of detection (LOD, Section 3.6.). The LOD of the QCM is good, but not strictly fantastic.
Two competitors are worth a mention. For gas sensing, the surface acoustic wave devices (SAW devices, [15]) have a better LOD. These are used in some electronic noses. For label-free biosensing in liquids, surface plasmon resonance spectroscopy (SPR spectroscopy [16]) also has an LOD better than the QCM. It is a question of both white noise and drift. Among the reasons to use the QCM is simplicity. Other reasons are connected to the depth of information. The QCM gives access to physical parameters beyond the mass per unit area. This is sometimes emphasized by calling the QCM an "instrument of surface analysis", rather than a "sensor".
Gravimetry in air or vacuum was the QCM's main use until the early 1980s [17]. At that time, Nomura in [18,19] and also Bruckenstein and Shay in [20] combined a QCM with an electrochemical cell and measured the mass transfer during electrodeposition. Attempts into that direction were made earlier but were less successful [21]. The first experiments with this "electrochemical QCM" (EQCM) were analyzed with the Sauerbrey equation. As was shown later, the Sauerbrey equation also applies in liquids as long as the layer is rigid (Equation (52), [22]). Of course, the liquid itself also has an influence. This influence is described by what today is called the Gordon-Kanazawa equation [23]. The Gordon-Kanazawa equation is a rediscovery. In slightly implicit form, it is, for instance, mentioned in [24]. Quite generally, a considerable body of experience on acoustic resonators in liquids was gathered in the 1930s to 1950s, using torsional resonators. Part of this knowledge is collected in Mason's book from 1948 [25]. The concepts, which underlie our current understanding of the liquid-phase QCM, mostly date from that period.
Later, there were two more additions to the techniques. Firstly, the resonance bandwidth was analyzed in addition to the resonance frequency [26] and the shifts of frequency and bandwidth were compared between overtones. Secondly, the oscillator circuits were largely replaced by passive interrogation (impedance analysis [27,28] and ring-down [29][30][31]). In recent years, the cost of impedance analyzers has come down [32]. These changes lead to the "advanced QCM". Another term is "QCM-D" for "QCM with Dissipation monitoring". "QCM-D" here is the name of a technique, not of one particular instrument. In this text, QCM is synonymous to QCM-D.

Forced Vibrations, Complex Resonance Frequencies
The following section motivates the complex frequency shift, Δf = Δf + iΔΓ [33]. The variable Γ denotes the half bandwidth at half height ("bandwidth" for short). The tilde denotes a complex parameter.
Start from the equation of motion of the forced resonator: m R is the mass. ξ R is the friction coefficient, also called "drag coefficient". In interfacial sliding, the "friction coefficient" is a ratio of two forces (tangential to normal). In liquid friction, it is a ratio of force to velocity. Renaming the force-velocity ratio as "drag coefficient" avoids this ambiguity. κ R in Equation (3) is the spring constant. We bring all terms containing x(t) to the left-hand side. The source term (the external force, F ext ) shall be of the form exp(iωt). The hat (̂) denotes a complex amplitude. Instead of exp(iωt), one might have also written exp(−iωt). That is a matter of convention, addressed in Box 1. For stationary oscillations of the form x(t) = x exp(iωt), the time derivative turns into a multiplication with iω: We divide by exp(iωt), divide by m R , rename ξ R /m R as 2γ, and rename κ R /m R as ω 0 2 : γ is the damping coefficient and ω is the natural frequency. Both have units of inverse seconds. The amplitude of displacement depends on ω as: -ω 2 m R x exp(iωt) + iωξ R x exp(iωt) + κ R x exp(iωt) = F ext exp(iωt) -

Box 1. Sign conventions.
Because the resonances of the QCM are extraordinarily sharp, the frequency of excitation, ω, is close to the natural frequency, ω 0 .
A side remark: For sharp resonances, the frequency of maximum displacement is the natural frequency. The natural frequency is called the "resonance frequency", here. For broad resonances, there is a slight difference between the natural frequency and the resonance frequency. The latter then is also called the "ringing frequency", equal to ω0(12γ 2 /ω0 2 ) 1/2 . One can always compute the ringing frequency from the natural frequency and the bandwidth. The difference is not of practical importance for the QCM.
If ω » ω 0 , the denominator can be simplified following ω 0 2 − ω 2 » (ω 0 +ω) When describing oscillations with complex numbers, one exploits Euler's relations, which imply that cos(ωt) = 1/2(exp(iωt) + exp(iωt)). In principle, all calculations containing the cosine should be carried out on the sum of exp(iωt) and exp(iωt). However, the two calculations with +iωt and with -iωt run in similar ways. One therefore carries out the calculation just once and eventually computes the (real) outcome of the calculation as Re(y ) = 1/2 (y + y *) where y is the outcome of the calculation for exp(iωt) and the asterisk denotes complex conjugation.
If entropy is supposed to always increase, the imaginary parts of certain complex response functions must have certain signs. The sign depends on whether the calculation is carried out with exp(iωt) or with exp(iωt). If exp(iωt) is chosen, the signs are: resonance frequency ω is real A wave propagating towards +z is written as exp(i(ωt -k )) = exp(iωt) exp(ikz) exp(kz).
The prefactor is often multiplied with an i and then hidden behind some normalization constant. Proceeding this way and separating the real and the imaginary part leads to: The first and the second term are shown as a black and a red line in Figure 3.
The complex resonance frequency plays out its strength when it comes to shifts thereof, called Δf in the following (Δf̃ = Δf + iΔΓ). The complex shift was proposed by Eggers and Funk [33]. Just about all equations predicting frequency and bandwidth can be formulated in terms of Δf . These equations cover Δf and ΔΓ at the same time.
The half bandwidth, Γ, is related to the energy dissipated per unit time, Ė, as: E is the energy contained in the oscillation.
In the authors' opinion, Γ is the best parameter for quantification of dissipative processes at the QCM surface. Γ puts frequency and bandwidth on equal grounds. For instance, the noise on Δf and ΔΓ is similar. Other parameters are in use. Some researchers use the full bandwidth, w = 2Γ, others use the Q-factor Q = f res /(2Γ), and still others use the inverse Q-factor Q −1 = 2Γ/f res and give it a new name and a new letter, namely "dissipation factor", D. Sometimes the "dissipation factor" is called "dissipation", for short. Δf may also be expressed in terms of the dissipation factor. The conversion is simplest for the overtone-normalized frequency shift: If ΔD is expressed in units of 10 −6 and if f 0 is 5 MHz, the conversion from ΔD [10 −6 ] to ΔΓ/n [Hz] amounts to a multiplication with 2.5.

Figure 3.
A typical output from impedance analysis. Panel (A) shows the conductance G el (black) and the susceptance, B el (red). Together, they form the complex electrical admittance, el = G el + iB el , which is equal to Z el 1 with Z el the impedance. The real part of the admittance forms the wellknown, symmetric resonance curve (assuming perfect calibration). This is different for the real part of Z el because of the parallel electrical capacitance, C 0 . G el (f) peaks at the series resonance frequency, f res . Panel (B) shows the polar diagram. Of interest in sensing are the shifts in frequency and bandwidth, f and ΔΓ (C).

Techniques of Read-Out
The methods of interrogation all rely on piezoelectricity and the plate's electrical impedance, Z el (ω), in the respective frequency range. Z el (ω) (and, also, the admittance, Y el (ω) = 1/Z el (ω)) form a resonance curve. The interrogation methods are sometimes grouped into "active" and "passive". In the active schemes, the resonator is part of an oscillator circuit. The amplifier contained in this circuit takes a certain influence on the Δf oscillation frequency. The other schemes are passive. Passive, however, does not mean that the apparatus would not take an influence on the resonance parameters, at all. (Even for the "grandfather clock" [34], the way of driving slightly affects the frequency.) For piezoelectric resonators, this influence is mediated by piezoelectric stiffening (Section 7.). The stiffness of a piezoelectric plate (and hence its resonance frequency) depends on whether the two electrodes are open, short-circuited, or connected across some electrical impedance, Z ext . The latter situation is realized in all electrical instrumentation controlling the resonator. Again, f res depends on Z ext . Calculating the impedance, which the crystal "sees", is nontrivial.
The different modes of interrogation differ from each other in cost, speed, and susceptibly to artifacts. They do not, actually, differ much from each other in precision (Section 3.6). References [35][36][37] cover the interface electronics in more detail.

Oscillator Circuits
An oscillator circuit is an amplifier with a resonator in the feedback loop. Because the resonator's impedance is small on the resonance frequency, the circuit spontaneously oscillates at this frequency. Oscillator circuits are the method of choice for clocks [38]. QCMs based on oscillators can be cheap, even after the frequency counter is included in the total cost. For advanced sensing, oscillators are problematic because the frequency of oscillation is not strictly equal to the frequency at the peak of the conductance curve. The latter frequency (the acoustic resonance frequency, also "series resonance frequency") is the frequency of relevance for interpretation. The parallel electrical capacitance (C 0 , Section 8.4.) takes an influence on the oscillation frequency. There are more intricacies in the details. These would not be a problem, if the small difference between the oscillation frequency and the series resonance frequency was constant, but this difference depends on damping and on details of the electronics. Oscillator circuits are available, which output bandwidth in addition to frequency (bandwidth being often converted to the dissipation factor) [39,40]. Oscillator circuits usually run on one harmonic only, often the fundamental.

Impedance Analysis
Impedance analysis [27] avoids the complications inherent to oscillator circuits. An impedance analyzer (synonymous to "vector network analyzer", "VNA") sweeps the frequency of excitation across the resonance. The resonance parameters are obtained from a fit of a resonance curve to the admittance trace. A suitable fit function is the phase-shifted Lorentzian, which is: The phase shift in Equation (13), φ, accounts for an asymmetry of the resonance curve. Imperfect calibration causes such an asymmetry. The asymmetry can be small, but it rarely vanishes. G max is an amplitude. The parameter G max does not contribute much to sensing. The product G max Γ is proportional to the effective area of the plate (Equation (113)). G max Γ sometimes varies slightly during experiment. How these variations depend on the sample's properties, is poorly understood Impedance analysis is among the passive techniques. "Passive", however, does not imply that the impedance analyzer would not affect the resonance frequency, at all. The analyzer's output resistance, its input resistance, and the length of the cables all take an influence on frequency and bandwidth because of piezoelectric stiffening. A second caveat: The resonance frequency as determined from the admittance trace depends on the (13) sweep rate. Impedance analysis is not quite as reliable as one would wish. Still: impedance analysis is rather transparent. The problems are noticed and their consequences can be quantified with moderate effort.
For measurements in liquids, the through ("thru") configuration is advantageous because it leads to a small current into the impedance analyzer. The small current is measured against zero background and may by amplified. The background is nonzero in the "shunt" configuration, which is also common and works well for experiments in air. In the shunt configuration (depicted in Figure 4C), a large impedance of the device under test lets the voltage from the output go straight to the input of the VNA. If the resonator's impedance is much larger than 50 , it causes small changes to this input against a large background. Because the background is amplified as well, amplification lets the detector run into overload. A resonator immersed in a liquid has a large impedance on resonance and should be wired in thru configuration. If grounding the front electrode is an issue, a transformer as shown in Figure 4B can be employed. Grounding the front electrode is advisable because the electrical properties of the sample may otherwise affect the resonance via piezoelectric stiffening.

Ring-Down
Resonant phenomena can always be probed in either the frequency domain or the time domain. As long as the dynamical equations are linear, the two modes of interrogation yield equivalent information ( Figure 5). One may either sweep the frequency of an AC excitation across the resonance (as in impedance analysis) or abruptly shut off the driving signal and watch the decay as a current trace on an oscilloscope (as in ring-down). The latter principle is implemented in the instrument marketed by Biolin Scientific (Västra Frölunda, Sweden). Figure 5. The time trace in ring-down and the resonance curve as a function of frequency are related to each other by a Fourier transform. The resonance parameters (mostly f res and Γ or D) can be obtained from both sets of data. The precision is similar.

Multi-Frequency Lock-in Amplification
The multi-frequency lock-in amplifier (MLA) stands between ring-down and impedance analysis. The MLA applies a comb of frequencies to the resonator. The resonance curve can be reconstructed from the signals returned to the MLA at these frequencies. The raw data are time-domain data, which are Fourier-transformed on the instrument's main board. In the time domain, the excitation amounts to a series of pulses. The current as displayed on an oscilloscope (lower right in Figure 6) visualizes the ring-down. The lefthand side and the right-hand side in Figure 6 describe the same process in the frequency domain and the time domain, respectively. The current response shows the resonance curve (bottom left). Transformed to the time domain (right), the excitation amounts to a sequence of pulses with a spacing in time of 1/(Δfcomb). The response is similar to a ring-down process (lower right). Downloaded from www.intermodulation-products.com/applications/lock-in on 13 February 2019.
The time between two pulses, Δtcomb, sets the time resolution. Δtcomb is equal to the inverse frequency spacing between two members of the comb, Δfcomb. The frequency spacing, in turn, must be smaller than the bandwidth of the resonance. Otherwise, the comb will miss the resonance. Using 32 frequencies, which are evenly distributed over a resonance with a width of about 3 kHz, one achieves a time resolution of 10 ms. However, one may also let the comb consist of only five frequencies and space those 500 Hz apart. The time resolution then improves to 2 ms.
Should a comb covering one particular resonance contain only, for example, 5 frequencies, the other 27 frequencies can be invested in the other overtones. The MLA can interrogate multiple overtones at the same time. In principle, one might worry about crosstalk between the different overtones excited in parallel, but this does not appear to be a problem, in practice [41].

Fast Measurements, Modulation Experiments
Data acquisition rate is critical for the study of transient phenomena. It is particularly important in analytical electrochemistry, which often exploits transients [42]. Most current advanced QCMs are not particularly powerful in terms of time resolution. Typical data acquisition rates are between 1 and 10 data points per second (for experiments in liquids) [43]. While the MLA in the comb mode is faster than most other instruments, one would still wish for more.
The data acquisition rate can be further improved, if the analysis is based on the electrical admittance at one, fixed frequency [44][45][46][47][48]. Following [47], we call this mode the "fixed-frequency-drive" (FFD) mode. As sketched in Figure 7, there is a one-to-one correspondence between the electrical admittance at this one frequency and the complex resonance frequency, f res + iΓ. This mode of data acquisition suffers from electrical artifacts, though. The conversion from G el + iB el to f res + iΓ assumes that the other parameters of the fit function in Equation (13) (given as G max , φ, G off , and B off ) are constant, which is not always the case. Figure 7. Impedance analysis requires a determination of the entire admittance trace (full and dashed lines), which takes time. Fast measurements can build on the conductance, G el , and the susceptance, B el , at one single frequency close to the center of the resonance, f c . Shifts in G el and B el are converted to shifts in f res and Γ (blue arrows). As the top panel shows, a shift in resonance frequency changes B el , but leaves G el unchanged. Conversely, an increase in bandwidth lowers G el , but leaves B el unchanged (bottom).
Even in the fixed-frequency-drive mode, it is difficult to achieve data acquisition rates beyond 2πΓ because the resonator remembers previous resonance conditions on a time scale of (2πΓ) . The memory is related to ringing up and ringing down after the resonator is turned on or off. The details are complicated. They can be studied with a varicap diode wired in series with the resonator. Switching the capacitance of the diode changes the stiffness of the quartz plate (piezoelectric stiffening, Section 7.), thereby rapidly switching the resonator's natural frequency. The resonance frequency as determined from impedance analysis follows with delay. Deconvolution of experimental data with a memory kernel should be possible but has not been done so far. Without deconvolution, the time per data point cannot be less than about 100 µs (for experiments in water, where Γ is a few kHz [47,48]).
Fast measurements may be combined with accumulation and averaging. This requires an experimental setting, where the sample responds to a periodic stimulus of some kind. Among others, a suitable parameter for modulation is the electrical DC potential of the front electrode when this electrode at the same time is the working electrode of an electrochemical setup. The instrument then operates as an electrochemical QCM (an EQCM, Section 9.1.). Figure 8 shows an example. The sample is an aqueous electrolyte. When the voltage of the front electrode is switched, Δf and ΔΓ respond, but they do so with a delay. The delay is linked to the kinetics of double layer recharging [48].
Modulation and accumulation avoid a critical problem of the QCM, which is drift. Typical QCMs drift by about 1 Hz/h when the crystal was mounted carefully and when all static stresses have relaxed. Otherwise, the drift can be much larger. The drift is mostly caused by migration of crystal defects, possibly also by insufficient control of temperature. The drift can hardly be prevented, but it can be circumvented by choosing the target of research suitably. The study of fast, repetitive processes does not suffer from drift because the average (taken over the period of the repetitive process) can be subtracted from the time traces. The average will drift, but the difference from the average can be accumulated over extended periods of time. The data shown in Figure 8 have been accumulated overnight. As a side remark: Oscillators also allow for fast data acquisition, as demonstrated in [49].

Noise and Drift
According to the conventions in sensing, the limit of detection, LOD, is three times the rms noise. Noise is one of the reasons, why surface plasmon resonance spectroscopy (SPR spectroscopy) is more widespread in label-free biosensing than the QCM. A second reason is drift. Noise has in depth been studied for clocks [50], but not to the same extent for the liquid-phase QCM.
To the best of the authors' knowledge, the different techniques and instruments driving the liquid-phase QCM reach a frequency noise in a similar range. A convenient way to calculate a drift-corrected noise from any data set builds on the Hadamard variance, which is: The Hadamard variance is zero for a straight, sloped line. After a linear fit is subtracted from a sloped line with added white noise, the root-mean-square noise (rms noise) of this data set is equal to the square root of the Hadamard variance of the original data set. Basing the definition of the noise on the Hadamard variance avoids the linear fit. The drift-corrected rms noise is δf In the authors' laboratory, measurements in water lead to an rms noise on Δf/n of about 30 mHz/Hz 1/2 . The noise depends on the time interval of data acquisition, hence the Hz 1/2 in the denominator. One can always lower the noise by averaging over longer times (assuming white noise). 30 mHz is the noise, if the instrument outputs one data point per second. Similar noise is seen in most figures in the published literature (which is a rough estimate, evidently). Reference [51] reports similar noise for resonators in liquids driven with an oscillator circuit (as opposed to impedance analysis or ring-down). Using a density of 1 g/cm 3 , an LOD in frequency of 90 mHz (which is 3 × the rms-noise) corresponds to an LOD in adsorbate thickness of ~0.05 nm.
In time and frequency control, the Allan variance is employed more commonly than the Hadamard variance. When frequencies are determined in time intervals of τ = 1 s, the Allan variance is: y is the fractional frequency shift. Good quartz clocks achieve σy(τ = 1 s) ≈ 10 -11 [52]. For the QCM in water, the noise is larger by a factor of about 1000. With a noise of 30 mHz on a 5 MHz signal, σy(1 s) is 6 × 10 9 . The Q-factor, on the other hand, decreases by only a factor of about 30. The noise is not proportional to Q −1 , as one might expect.
One can understand that the frequency noise increases stronger than Q −1 on a qualitative level. The discussion can build on the fixed-frequency-drive mode (Section 3.5.). In a liquid environment, the large damping increases Γ and it also decreases Gmax (Equation (13), Figure 9). More generally, the noise has a white component, which scales as the ratio of the thermal energy, kBT, to the power going into the device. When immersing a resonator into a liquid, the power into the device decreases because the resistance (R ) increases. This amounts to a first factor of about Q −1 . A second factor of Q −1 enters, when a noise in Y el is translated to a noise in frequency. With Q decreasing by about a factor of 30, the noise (following this rough argument) increases by about a factor of 30 2 . Figure 9. When a resonance becomes broader, this lowers the precision in the determination of the frequency shift twice. Γ increases and G max decreases (both in proportion to the inverse Q-factor). The noise on δB is assumed to be independent of Q in this argument.
There may be other sources of noise. For instance, the liquid-phase QCM is susceptible to acoustic vibrations. Slamming the door leaves a trace in Δf(t). The coupling is mediated by bending of the plate as described in [53]. More generally: There probably is room for improvements on the frequency noise of the liquid-phase QCM.

The Acoustic Multilayer Formalism and Its Consequences
The following sections describe the quantitative analysis of QCM data acquired on planar samples in detail. As long as the samples are homogeneous in the surface plane, the acoustic multilayer formalism achieves the modeling.

Qualitative Data Inspection
Before starting a fit, some qualitative considerations are worthwhile:


Is -Δf ≫ ΔΓ and is -Δf/n ≈ const.? Did the experiment occur in air? If so, the response is probably dominated by inertia in the sense of the Sauerbrey equation ("inertial loading"). With a density of 1 g/cm 3 and 5 MHz crystals, a layer thickness of 1 nm leads to -Δf/n = 5.7 Hz. Did the experiment occur in liquid? If so, the response is probably dominated by the formation of a thin layer. However, -Δf/n may be smaller than 5.7 Hz per nanometer in case the film is soft (Equation (52)).  Is -Δf » ΔΓ, is -Δf/n 1/2 ≈ const., and was the resonator immersed in a liquid? If so, the response is probably dominated by changes in viscosity (Equation (29), "viscous loading"). With 5 MHz crystals, -Δf/n 1/2 = 716 Hz corresponds to a viscosity of 1 mPa s (slightly more than the viscosity of water).  Is Δf > 0 and is Δf·n ≈ const.? If so, the response may be dominated by point contacts ("elastic loading", Section 5.1.).  Do Δf and ΔΓ show unexpected patterns? If plots of ΔΓ versus Δf show circles or spirals, the data may originate from a coupled resonance (Equation (79), Section 6.).

The Small-Load Approximation in 1D (Parallel-Plate Model)
In the following, we go beyond the equation ω 0 ≈ (κ R /m R ) and formulate a continuum model. The lumped-element description from Figure 1B is abandoned (no discrete springs, no discrete masses). We treat the resonator as a vibrating body, similar to the bell shown in Figure 1C. Piezoelectric stiffening is ignored, for now. Piezoelectricity at this level simply is a convenience, which allows to probe acoustic resonances by electrical means.
In the continuum picture, a resonance amounts to a displacement pattern, which occurs time-harmonically and which is easily excited to a large amplitude. The deformation pattern u(x,y,z,t) is: The amplitude u (x,y,z) is the mode of vibration. Because the displacement always occurs along x, u may be viewed as a scalar (rather than a vector). Further simplifying the problem, we let all gradients in the plane vanish. The "parallel plate" can be viewed as an "infinite parallel plate". It can also be a plate with finite area, A eff , but the edges must not affect the mode of vibration (which is unrealistic for AT-cut quartz because of its anisotropic elasticity). Within the parallel-plate model, the amplitude of displacement, u ≡ u x , is a function of z, only.
The resonant modes of vibration are solutions of a boundary value problem. The boundary condition here are surfaces, which are free of stress. Because the shear stress is proportional to the shear strain, the shear strain must vanish at the surface (at z = 0 and z = d q ): The origin of the z-axis in Equation (17) is at the back of the plate. For the parallel plate, the modes of vibration are standing waves: u(x,y,z) = u (x,y,z) exp (iωt) (16) σ z = 0, = G q du dz z = 0,d q = 0 (17) u n (z,t) = u S cos(k n z) exp(iωt) (18) u S is the displacement amplitude at the surface. The boundary condition fixes the wavenumber to discrete values, which are: where n is the overtone order. Expressed differently, the wavelength, λ, must be an integer fraction of twice the plate's thickness. Such discrete sets of solutions are characteristic of boundary value problems. Small deviations of the overtone frequencies from the integer multiples of the fundamental frequency are discussed in Section 8.1. Critical to the above argument was the fact that the resonator surface was stress-free. The surface must coincide with an antinode. When a sample exerts a periodic stress onto the surface, the resonance condition changes. Within linear acoustics, the stress, σ S , is proportional to the displacement, u S . In acoustics, is customary to not discuss the stress-displacement ratio, but rather the stress-velocity ratio, which is the impedance. Velocity, v S , and displacement, u S , are related as v S = iωu S . The stress-velocity ratio at the resonator surface is the load impedance, Z L . The load impedance is a key variable in the physics of the QCM. The displacement and the stress at z = d q are u S cos (k q d q ) and u S G q k q ( sin (k q d q )) sin (k q d q )), respectively. The stress-velocity ratio follows as: The first minus sign occurs because the stress is exerted by the sample onto the resonator surface (in the direction of −z). It follows that: The relations c q = (G q /ρ q ) 1/2 and Z q = (G q ρ q ) 1/2 were used. Z q , c q , and f ref are complex, in principle. However, the resonator's intrinsic losses are not of interest in sensing. Writing them as real parameters certainly affects the absolute value of the bandwidth, but not its shift induced by the sample. The use of an effective complex k q (rather than k q = /c q ) is justified in Appendix D.
The relation d q = c q /(2f 0 ) leads to: Equation (22) is an implicit equation in Δf , which can be solved numerically. It can also be turned into an explicit equation in Δf = Δf + iΔΓ by: -linearizing the tangent as tan(πΔf /f 0 ) ≈ πΔf /f 0 -evaluating the load impedance Z L (f) at the frequency of the unloaded crystal, rather than the resonance frequency in the presence of the load.
This explicit equation is: Equation (23) is the small-load approximation applied to the parallel plate. Given its importance, it is written down in slightly different form one more time: Δf f All terms have been normalized to overtone order, as is common in gravimetry. Also, the shift in the dissipation factor, ΔD, was used in step 3, replacing 2ΔΓ/f ref . This section deals with stratified layer systems. For those, the stress-velocity ratio follows from how the shear wave bounces back and forth inside the sample. Three simple cases are sketched in Figure 10. However, Equation (23) is more general. Should the sample be structured, laterally, the load impedance may be replaced by its area average: Area averaging is possible, because Equation (23) is linear in the load impedance. Within the parallel-plate model, ⟨…⟩ area is an unweighted area average. For more realistic resonators, the square of the local amplitude, |u S (r S )| 2 , must be included as a weight function (Section 8.1.): r S is a point on the resonator surface.
The calculation of the stress-velocity ratio is rather simple for thin rigid films and it is also simple for semi-infinite media. For layered systems, there are analytical equations (not all equally simple) [54][55][56]. Some of them are discussed below. For single contacts with small contact area, the stress can be replaced by the ratio of the restoring force, F, and the acoustically active area, A eff . Should the sample have a more complicated structure, the stress-velocity ratio needs to be calculated numerically, solving the equations of continuum viscoelasticity for the given geometry (Section 5.3.).

Inertial Loading
For a thin rigid film as shown in Figure 11, the stress at the resonator surface is governed by inertia. From Newton's third law (force = mass × acceleration) it follows that the stress exerted onto the surface is iωm f v S with m f the mass per unit area. The load impedance is Z L = iωm f . This leads to the Sauerbrey equation: Angle brackets denote area averaging as before.  Figure 11. The Sauerbrey equation applies to thin rigid films (A) and other rigid samples (B). In the latter case, the mass per unit area must be area-averaged.

Semi-Infinite Viscoelastic Media
For semi-infinite, homogeneous viscoelastic media, the load impedance is equal to the shear-wave impedance, Z bulk : The load impedance and the shear-wave impedance must not be confused. Z L is the area-averaged ratio of stress and velocity at the resonator surface. The shear-wave impedance, Z or Z bulk , is the stress-velocity ratio of a propagating shear wave. Z is a materials constant, given as Z = ρG 1/2 = ρ G /ρ 1/2 = ρc . The wave impedance governs the reflectivity at interfaces (Equation (37)).
The relations Z bulk = ρG 1/2 = ρc and G = iωη inserted into Equation (28) lead to: Equation (29) is the Gordon-Kanazawa relation [23,24]. If η is independent of frequency, Δf and ΔΓ scale as n 1/2 . The Gordon-Kanazawa relation can be inverted for viscosity as: In more compact notation, one may write The density was moved to the left-hand side in order to emphasize that the QCM measures the viscosity-density product (or, equivalently, the product ρG ). The density often is known and it often varies less than the shear modulus. For instance, adding a polymer to a solvent much increases the viscosity but leaves the density unchanged within a few percent. Still: ρ and η cannot be determined separately, using the Gordon-Kanazawa relation. (Such a separate determination can be achieved with porous coatings [57].) The reference state must be the resonator in air, if the target of the study is a small change in the viscosity (for instance caused by a change in pH). That is so, because ρG and ρη depend on the square of Δf . If a change in pH causes a slight change in viscosity, δη , this causes a change in frequency, δf , following δ ρη ≈ -1/(iω) πZ q /f 0 2 2Δf δf . The term in square brackets is the mixed term of the binomial. A term of the form δf was neglected. Δ and δ have different meanings. Δ denotes the difference from the dry state, δ denotes the small shift induced by the change in pH. For viscoelastic media, ΔΓ is larger than -Δf. For a purely elastic medium (η′ = 0, η′′ > 0, or G′ > 0, G′′ = 0) the frequency shift vanishes and ΔΓ is equal to (ρG′) / f 0 / πZ q . This result may appear as counterintuitive, given that the bandwidth is usually associated with dissipative processes. However, ΔΓ quantifies all forms of energy withdrawn from the resonator (Equation (11)). The energy may or may not be dissipated inside the medium. (Of course, it is dissipated eventually, somewhere.) For the semi-infinite elastic medium, the energy is radiated away towards z = +∞ and ΔΓ is nonzero for that reason.
The displacement pattern in a Newtonian liquid is shown in Figure 12. It is of the form: In the last step, a Newtonian liquid was assumed (η'' = 0, η' independent of frequency). The wave number is then given as k = (1  i)/δ, where δ is the depth of penetration: With ρ = 10 3 kg/m 3 , η = 10 −3 Pa s, and ω = 2πn × 5 MHz, the depth of penetration is δ = 252 nm/n 1/2 . These values inserted into the Gordon-Kanazawa relation predict −Δf = ΔΓ = 716 Hz/n 1/2 . Figure 12 clarifies what "semi-infinite" means for the QCM. The sheared layer seen by the QCM is around 200 nm thick (depending on overtone order and viscosity). The finite thickness of the sheared layer turns the liquid-phase QCM into a surfacespecific instrument. This is expressed diagrammatically in Figure 13. If exposed to a fluid, the QCM does not see the bulk outside the sheared layer. This is strictly correct to the extent that the QCM indeed vibrates in a pure thickness-shear mode. There are small flexural admixtures to the mode of vibration (Section 8.1.). Because of these, Δf and ΔΓ are slightly sensitive to objects outside the sheared layer, whenever these scatter compressional waves.
Can the quartz crystal microbalance be turned into a quartz crystal viscometer? Firstly, there are other simple ways to measure viscosity. Further working against the QCM are artifacts, which are caused by compressional waves ( Figure 13) and by the adsorption of debris to the resonator surface, acting as a Sauerbrey load. The Sauerbrey load and the Gordon-Kanazawa load can be separately quantified with the advanced QCMs, but only with these. Problems with mass deposition have a characteristic signature in QCM-based viscometry, which is an apparent negative η''. η'' is proportional to ΔΓ 2  Δf 2 (Equation (30)). If some adsorbate lowers the frequency following Sauerbrey, this may drive the apparent η'' into the negative range when data are analyzed with Equation (30).

Figure 13.
When shear waves dominate the resonator's response, the QCM response is surfacespecific. Flexural admixtures to the mode of vibration and the concomitant compressional waves may spoil surface specificity. Compressional waves may be reflected somewhere in the bulk and return to the crystal.
More conceptually, the QCM determines the viscosity at a frequency of a few MHz. For small-molecule liquids, the steady-shear viscosity and oscillatory-shear viscosity at a few MHz are similar. The more interesting fluids, however, often contain soft matter with varying degrees of complexity, which entails relaxation and viscoelastic dispersion. The high-frequency viscosity then may be different from what the engineer cares about.
Torsional resonators [58][59][60][61] mitigate these problems by virtue of their lower frequency. They are less sensitive to the deposition of mass than the thickness-shear resonators and their frequency (tens of kHz) is closer to the frequencies and time scales of practical relevance. Torsional resonators are commercially available as viscosity sensors [62]. References [63][64][65] report on the use of kHz resonators for an array of other purposes.
One may envisage a role for the QCM in viscometry, when it comes to small sample volumes. The problem is of much technical relevance and was addressed with other miniaturized sensors, as well [66]. One may deposit small droplets onto the resonator surface. The shifts in frequency and bandwidth are correspondingly small, but they are still well above the noise. One may determine the contact area, Ac, with a camera and attempt to derive the viscosity from the relation [67,68]: A eff is the acoustically active area of the plate and K A is a function of the droplet area, which takes the amplitude distribution, v S (r S ), into account (Section 8.1.). Typically, this function would be determined by calibration, using liquids with known viscosity. A similar analysis can be applied to a combination of a QCM with a JKR apparatus [68]. The JKR apparatus pushes a lens of a soft material against a substrate and determines the contact radius as a function of the normal force. The JKR apparatus targets the contact-mechanics of soft materials. The substrate may be a QCM, in which case the complex frequency shift reports the material's high-frequency shear modulus.
From Equation (34), one would expect the prefactors to be the same for frequency and bandwidth. If that was so, the ratio ΔΓ/(−Δf) would be independent of contact area and related to the material's loss tangent, tan(δL) = G′′/G′. Experiment shows, however, that the ratio slightly depends on contact area. The problem has to do with the fact that the degree of energy trapping changes when a sample contacts the resonator in the center, only. Energy trapping increases the resonance frequency (see the discussion in in Section 8.1.). Its effect on bandwidth is different from its effect on frequency [69]. These problems must be kept in mind when analyzing Δf with Equation (34). QCM-based viscometry amounts to high-frequency rheology on bulk samples. When applied to engine oils, the QCM's high frequency is a disadvantage. For other complex fluids, the high-frequency viscoelasticity actually is of interest because it depends on the fluid's internal organization. Pharmaceuticals for parenteral administration are often formulated as concentrated protein solutions, which display viscoelastic relaxation in the MHz range [70]. Protein-protein interactions (PPIs) may turn these solutions into weak gels. High-frequency rheology is among the techniques probing such interactions [70][71][72]. Figure 14 shows an example, taken form a study on pharmaceutical formulations. For the data in Figure 14A, the shear modulus as derived with Equation (31) displays viscoelasticity with a characteristic dependence on frequency. The data could be fitted with the Maxwell model (a lumped-element model consisting of a spring in series with a dashpot). The relaxation time was in the range of a few tens of nanoseconds. The data shown at the bottom did not show this kind of viscoelasticity. This formulation looked like a Newtonian liquid to the QCM.
Note that time-temperature-superposition (TTS) is not needed to interpret these experiments. For certain types of polymers, one may study the "high-frequency" viscoelasticity with conventional, low-frequencies rheometers by cooling the sample, such that the relaxations of interest slow down and then are accessible to the instrument. Complex liquids often are not thermorheologically simple in this sense. For these, acoustic instrumentation operating at high frequencies (such as the QCM, but not limited to the QCM [73]) is needed.

Films in Air
If the sample contains interfaces with some impedance contrast, the reflected wave contributes to the periodic stress at the resonator surface in proportion to -Z f v ← , where v ← is the amplitude of the reflected wave, evaluated at the resonator surface ( Figure 15). Z f is the film's wave impedance. There is a minus sign because the reflected wave travels towards -z. Dividing by the total velocity, v → + v ← , the load impedance is found to be: The ratio r̃S = v ← /v → is the complex reflectivity evaluated at the resonator surface. The QCM may be viewed as an acoustic reflectometer (and may compete with other reflectometers, for instance described in [74,75]). In particular, there is a close correspondence

G' G''
between the physics of the QCM and optical reflectometry (Section 9.2.). The reflectivity of the sample can be inferred from Equation (35), solved for r̃S: Equation (23) was used when expressing α as a function of Δf + iΔΓ. Figure 15. The stress at the surface of a resonator coated with a film contains a contribution from the reflected wave.
The calculation of r̃S for a film in air is sketched in Figure 15. r̃S is given as exp −2ik f d f × r̃ where the exponential covers the propagation through the film (twice, hence the factor of 2) and r̃ is the reflection amplitude at the film-air interface. The reflectivity of a wave at an interface between two media with different wave impedances, Z 1 and Z 2 , is: The proof exploits that the velocity and the stress are continuous at the interface. Equivalently, the reflectivity can be calculated from the conservation of energy and momentum. It is worthwhile to remind oneself of two related situations: -When an optical wave hits an interface at normal incidence, the reflectivity is n r,1  n r,2 / n r,1 +n r,2 . While one might think so, the refractive index, nr, is not strictly the same as the impedance of the optical wave, but it is related to this impedance.
-Upon a central elastic collision of two spheres, the velocity of the first sphere after collision is v ← = v → (m m )/(m +m ). The mass takes the role, which the impedance has for waves.
Because Z air = 0, the reflectivity at the film-air interface is unity. From Equations (35) and (37), the load impedance follows as: Euler's relation implies: This result inserted into Equation (22) yields: Equation (40) was first derived by Lu and Lewis [76]. The Lu-Lewis equation does not invoke the small-load approximation. It is an implicit equation in Δf + iΔΓ, which must be solved numerically. An analysis of frequency shifts based on the Lu-Lewis equation is implemented in some commercial film-thickness monitors. The algorithm is called "Z-match method" [77,78]. Film-thickness monitors often become heavily loaded when crystals are not replaced between deposition runs.
A side remark on film thickness monitors: The frequency shift in these instruments can be above 1 MHz. Resonators with plane-convex surfaces are used, which tolerate large loading but only work reliably on the fundamental. This is one of the cases, where the small-load approximation does not apply.
In order to apply the Z-match algorithm to these data (Δf determined on a single overtone), the wave impedance of the layer must be known. Some values for metals are tabulated in [79]. One might also use a numerical solution of the Lu-Lewis equation (Appendix C.2) as part of a fitting process, determining not only the thickness, but also the layer's viscoelastic parameters [80]. This analysis of course requires experimental values of frequency and bandwidth on a few overtones as input. Otherwise, the problem is underdetermined.
The Lu-Lewis equation does not make use of the small-load approximation. If the small-load approximation is employed (which amounts to inserting Equation (39) into Equation (23)), the following result is found: Figure 16 shows Δf and ΔΓ as predicted by Equation (41). The following sections address the four different regimes indicated with arrows in Figure 16.
Δf + iΔΓ f This is the Sauerbrey result. The relations Z f = ρ f ̃ and = ω/ ̃ were used. Δf/n is proportional to the film's mass per unit area. ΔΓ vanishes because the film does not undergo shear deformation to any appreciable extent under its own inertia. Again, the Sauerbrey result is more general than Equation (42) because area averaging may be applied.

Infinite Thickness
In the limit of infinite thickness, Z f turns into Z bulk . The tangent turns into i as long as k f '' > 0: The Gordon-Kanazawa relation is recovered.

Thin Viscoelastic Films
If the film is thin, but still thick enough to let viscoelasticity be noticeable, the tangent can be expanded to 3rd order as tan( . This regime is of much practical importance. The Taylor expansion leads to: The It is instructive to express Δf /n as a function of n 2 : Clearly, both Δf/n and ΔΓ/n depend on n 2 [81]. If (!) J f ' and J f '' themselves do not depend on frequency, Equation (45) describes a linear relation between Δf/n and ΔΓ/n, on the one hand, and n 2 , on the other. The slopes then are proportional to the elastic compliance, J f ', and the viscous compliance, J f '' ( Figure 17). However, J f ' and J f '' may depend on frequency, in which case the lines in Figure 17 have some curvature. Figure 17. Plots of Δf/n and ΔΓ/n versus n 2 as motivated by Equation (46). The film thickness as derived from the offset in panel (A) is about 1.6 µm. The slope d(Δf/n)/d(n 2 ) is almost constant, indicating that J f ' only weakly depends on frequency. There is such a dependence on frequency in (B). The data were taken on a spin-cast film of polyisobutylene. Adapted from [81]. offset ® m f When the film is not much softer than the crystal, a correction to Equation (45) is needed because the assumptions inherent to the small-load approximation produce a sizeable error. In the derivation of Equation (23), the term tan(πΔf /f 0 ) was linearized, while the term tan k f d f from Equation (41) was expanded to 3rd order. This is an inconsistency, which can be removed with a systematic perturbation calculation [82]. Box 2 addresses the issue in more detail.
Box 2. 3rd-order perturbation applied to the films on a parallel plate.
On the way to Equation (44)), a tangent contained in the load impedance was Taylor-expanded to 3 rd order, while a similar tangent contained in the Lu-Lewis equation (Equation (40))) was linearized. This is an inconsistency, which can be removed [82]. At the same time, one needs to deal with the fact that the load should be evaluated at the resonance frequency of the loaded crystal, rather than the reference frequency. Dealing with these complications, one should also take electrode effects into account, that is, treat the 2-layer system (electrode plus film).
For the 2-layer system, the Lu-Lewis equation turns into The indices e and f denote the electrode and the film, respectively. The equations become more compact if the following variables are used: Importantly, the coefficient to n 2 contains the thickness and the wave impedance of the electrode. If the electrode is neglected, the equation simplifies as The terms independent of n slightly modify the Sauerbrey equation. They are negligible in practice. The n-dependent term leads to Equation (46)).
In view of these complications, one might also go back to the Lu-Lewis equation, solve it numerically, and use this solution when fitting a model to experimental data. The Python code in Appendix C.2 solves the Lu-Lewis equation.
The above remarks mostly concern films in air. For films in liquids, Equation (49)) can be trusted. At least, it is not grossly invalidated in quantitative terms by the full numerical solution to the Lu-Lewis equation.
Piezoelectric stiffening is not covered by this formalism. Piezoelectric stiffening does affect the result, in principle, but the changes are small [80].
3rd-order perturbation leads to: The difference between Equations (45) and (46) is essential for stiff films (with small J ). Analyzing data from thin, glassy polymer films with Equation (45) can easily produce a negative apparent shear modulus. The analysis becomes even more complicated if the viscoelastic properties of the electrode are taken into account. For more details see [80].
Equation (46) is the basis of quantitative rheometry on thin films [83,84]. The QCM as an instrument is unique in this regard because it does not require a clamp on the other side of the film. Thin layers can be clamped from both sides [85], in principle, but these experiments (using the surface forces apparatus, SFA) are more demanding than the QCM. Films in air shear under their own inertia, hence the proportionality to m f 2 in the viscoelastic correction. A film thickness of a few tens of nanometers is needed to see viscoelasticity. Thinner films would have to be extremely soft to show such effects. If Δf/n does not agree between overtones for such films, this may go back to an overtone-dependent modal mass (Section 8.6.). An increase in bandwidth has been seen in experiments on monolayers of noble gases [86]. This increase in damping is not easily explained in the standard framework of molecular interactions. Superlubricity may be an explanation [87].

The Film Resonance
At d f » λ/4, Equation (41) hits the "film resonance" [8890]. The film resonance is an example of a coupled resonance (Section 6.) and is therefore labeled with subscript CR. In the thickness range of the film resonance, Δf increases with thickness and ΔΓ goes through a maximum. The dependence of Δf + iΔΓ on d f looks like a resonance (cf. Figure 3A) and we briefly convince ourselves that the algebra confirms that. If the imaginary part of k f is small, the real part of the tangent at the pole first goes to +∞ and later returns from ∞. Close to the pole, one may write tan(x) = 1/cot(x) = 1/cot(y + π/2). The variable x was substituted by y = x  π/2. Taylor expansion of the cotangent to 1st order in y leads to cot(y + π/2) = y. The tangent turns into 1/(x  π/2) and Equation (41) is approximated by The resonance frequency of the film, ω' CR , is governed by the condition k f ′d f = π/2. An experimental example of a film resonance is shown in Figure 18. For the study of soft films, the film resonance mostly is a problem. The data can rarely be fitted well by Equation (41).

Figure 18.
An experimental example of a film resonance. The sample consisted of a polyelectrolyte multilayer, the thickness of which was gradually increased by repeated dipping. When displayed in polar form, the data form a spiral. These authors display the load impedance rather than the complex frequency shift. The load impedance here has the unit "Rayl", in the honor of John William Strutt, 3rd Baron Rayleigh. 1 MKS-Rayl is equal to 1 kg/(m 2 s). Red and black data denote results obtained with an even number of layers and with an odd number of layers, respectively. The difference occurs because anionic and cationic layers alternate. Adapted from [91].
The film resonance is among the examples, where the small-load approximation is not quite good enough. The problem occurs if the film is not lossy, that is, if k f ′′ ≪ k f ′. In the range of the film resonance, the Lu-Lewis equation then has two solutions, corresponding to a "symmetric" and an "antisymmetric" mode ( Figure 19). For k f ′d f < π/2, the mode with negative Δf has the larger amplitude. When k f ′d f » π/2, the mode with positive Δf grows in amplitude and eventually takes over. This picture emerges in the frame of parallel-plate model (not using the small-load approximation). More specifically, it emerges, when the calculation of the electrical admittance is done with the Mason equivalent circuit, briefly mentioned at the end of Section 8.4. For the details, we refer the reader to the literature [54,92]. Figure 19. If the small-load approximation is abandoned, a film resonance produces two separate peaks in the conductance trace, G el (ω), corresponding to two separate modes of vibration. The shear gradients inside the film have opposite sign for the two modes, hence the labels "symmetric" and "antisymmetric". Far away from the coupling condition, ω 0 and ω CR are not affected by coupling (to the right and to the left in panel (A)). If the two frequencies match and if the two modes indeed couple, anticrossing results. For a more quantitative treatment see Chapter 4.63 in [5]. Panel (B) shows an enlargement of panel (A) in the region of anti-crossing. If the bandwidth is large, the two modes are not actually resolved (sketched in orange in panel (B)). When the film becomes thicker, ω CR approaches the coupling condition from above (green arrow). Δf is less than zero, following Sauerbrey. Because the two modes are not resolved at the coupling condition, the peak in G el (ω) is broad. The center of the peak gradually moves up, because the antisymmetric mode becomes stronger. Eventually, the peak sharpens again and returns to the original frequency from above. A similar behavior is seen in Figure 16, based on Equation (41).
Modes growing and shrinking in magnitude are seen in experiment when swelling polymer films with d f » λ/4 in solvent vapor. However, an antisymmetric mode as shown in Figure 19 is not easily identified. Most often, one of the anharmonic sidebands grows and eventually becomes the largest peak in the conductance trace. A full understanding of the film resonance would require a realistic model of the resonator in 3D. That is difficult, in the first place. Also, the patterns of the anharmonic sidebands (Section 8.2.) show some variability between crystals and batches. Crystal imperfections play a role.
If the details of the film resonance are so difficult, why even bother? Firstly, the film resonance is an instructive example of a coupled resonance. The coupled resonance will concern us further in Section 6. Also, the film resonance is occasionally seen in experiment. Figure 20 shows an example. These authors were interested in the dissolution of polymer films and in the preservation of old paintings. The QCM worked well, basically, but the frequency at some point jumped. Knowing about the film resonance, one understands the jump.

General
Many adsorbates from the liquid phase do not have a sharp interface with the bulk. Still, the viscoelastic box profile (that is, the homogeneous film with thickness d f ) is a good starting point. Using the reflectivity at an interface from Equation (37), the frequency shift is [93][94][95]: Applying Euler's relation to the right-hand side (similar to Equation (39)) yields While not immediately evident, Equation (49) is equivalent to the Voigt-model from [56] and, also, to [54,91].
Equation (49) also leads to a film resonance, but the resonance condition is different from k f ′d f = π/2 (which is the resonance condition for dry films). The film resonance can be seen while films grow from a liquid phase or while they dissolve into a liquid phase [96]. It is also observed when surface-attached gels [97] or polymer brushes [89] swell and deswell. In the latter case, however, the layer becomes softer as it swells. k f and Z f vary strongly.
For thick gels swelling in a solvent, Δf as a function of the swelling degree can also be qualitatively portrayed as a transition from Sauerbrey-type behavior to Gordon-Kanazawa-type behavior. The compact layer obeys the Sauerbrey relation, possibly with a small viscoelastic correction. -Δf/n becomes larger as the layer swells. At some point, the layer thickness is comparable to the depth of penetration of the shear wave. Beyond this point, the layer appears as a soft semi-infinite medium. Δf and ΔΓ then report the medium's complex shear modulus, regardless of its thickness. The shear wave no longer reaches to the top of the film. An example for this behavior is shown in [97].
Cell cultures behave like soft gel layers in this regard. The shear wave usually does not reach to the top. It probes the layer's shear stiffness at the bottom rather than the layer thickness. Cell cultures have been extensively studied with the QCM [98][99][100]. The interpretation is usually based on certain correlations between Δf and ΔΓ, on the one hand, and the conditions of the experiment, on the other. Quantitative modeling is difficult. Δf/n may certainly be converted to an apparent mass, but the emphasis here is on "apparent". Δf + iΔΓ f

Thin Adsorbates
We now turn to thin viscoelastic layers. When the tangent in Equation (49) is expanded to 1st order in d f , one arrives at: It is convenient to choose the resonator immersed in the liquid as the reference state, which results in: The relation k f Z f = ωρ f was used in line 2. Equation (51) can be rearranged as: If the film is much stiffer than the liquid (if |Z f | ≫ |Z bulk |), this relation reduces to the Sauerbrey equation. For rigid layers and even moderately rigid layers, the Sauerbrey contribution and the Gordon-Kanazawa contribution to the complex frequency shift simply are additive [101,102]. This analysis approach is applied to electrogravimetry [103]. Electrodeposition and electroetching can be analyzed with the Sauerbrey equation. (Roughness may take an effect, though [104].) The second term in the square brackets in Equation (52) is a viscoelastic correction. This term differs characteristically from the viscoelastic correction in Equation (46) (experiments in air). The difference goes back to the fact that a film immersed in a liquid feels a stress from the other side. It is partially clamped by the liquid. In air, films are sheared by their own inertia, only. For films in air, viscoelastic effects are seen after expanding the tangent to 3rd order as k f d f + k f d f 3 /3. In liquids, viscoelastic effects enter the picture in 1st-order Taylor expansion, already. Even molecularly thin films are sheared by the adjacent liquid and can be studied with regard to their softness.
For soft films in liquids, the apparent mass as derived with the Sauerbrey equation is smaller than the film's mass [105]. Voinova et al. call this the "missing-mass effect" [106]. Viscoelastic effects can be recognized by plotting Δf/n and ΔΓ/n versus n ( Figure 21). A positive slope indicates a finite J f ''. A nonzero ΔΓ/n indicates finite J f '. Why does the slope have positive sign? (It has negative sign in air, Figure 17.) The missing-mass effect is caused by the film being clamped from the other side. The stress exerted by the liquid is proportional to ωη, meaning, increases with overtone order.  (54)). The slope in this log-log plot is the power law exponent β' (Section 4.7.). The sample is a block copolymer adsorbed to the gold surface, where the soluble part forms a brush-like structure. Data kindly shared by Anna M.C. Maan, University of Groningen, FSE-Zernike Institute for Advanced Materials.
Ideally, one would wish to derive J f ' and J f '' (or equivalently, the moduli, G f ' and G f '') on each overtone. This would amount to a rheological spectrum, albeit in a limited frequency range. Unfortunately, the problem is underdetermined as long as the film thickness is not known. The elastic compliance alone, however, can be determined on each overtone, at least approximately. In the thin-film limit, m f can be eliminated from Equation (52) by taking the ratio of ΔΓ and -Δf: It was assumed that ρ bulk » ρ film . For such thin layers, the ratio ΔΓ/(Δf) is independent of layer thickness. It is a materials parameter. Equation (53) further simplifies if the denominator is about unity. That is often an acceptable approximation because polymer films (even when swollen) are much stiffer than the ambient liquid. (They are stiffer at MHz frequencies than at low frequencies.) The denominator in Equation (53) can be , the denominator is unity, leading to [107]: J f ' can be rather robustly inferred from QCM data. This will be important for the discussion of viscoelastic dispersion at the end of Section 4.7.
In linear rheology, viscoelasticity is commonly expressed in terms of the shear modulus G = G' + iG''. For the QCM, it is more convenient to instead use the shear compliance, J = 1/G = J'  iJ'', because the compliance occurs in the numerator in Equations (51) and (46). The trivial case (Sauerbrey-like behavior) corresponds to zero compliance. The conversion between G and J is nontrivial because they are complex: Writing J' » 1/G' often is grossly inaccurate. The conversion is easier for the absolute values and the loss tangent: The loss tangent is the same for the modulus and the compliance. Another relation worth remembering is η = G /(iω) with η = η − iη = G′′/ω − iG′/ω the viscosity.
Because the film thickness enters Equation (51) as a linear term, shifts in frequency and bandwidth resulting from multiple films are additive (assuming ≪ 1 for all layers): Equation (57) holds in an integral sense: ρ(z) » ρ bulk was assumed in lines 3 and 4. Equation (58) may also be expressed in terms of density and viscosity: The term in square brackets is a contrast function. The integral can be viewed as a "shear-wave acoustic moment" of a profile of the polymer volume fraction, ϕ(z), as displayed in Figure 22.
Separating real and imaginary parts in Equation (59) leads to: One can apply the Sauerbrey equation to adsorbates and derive an apparent mass, m app (n), which is: In principle, the parameter m app is a complex function of n. Application of the Sauerbrey equation is most meaningful, if the imaginary part of m app is small and if m app weakly depends on n.
Δf  (59). The z-range, in which ρ(z) and η(z) are significantly different from the corresponding bulk values, must be much below the depth of penetration of the shear wave, δ, in order to let the integral formulation be applicable. The displacement profile is sketched at the bottom. For most polymers, ρ(z) is similar to ρbulk. The dependence of the viscosity, η (z), on the polymer volume fraction, ϕ(z), is strong and not easily guessed based on simple models.
For sufficiently rigid layers, the apparent acoustic thickness, m app /ρ, is close to the geometric thickness. It is slightly smaller because of the missing mass effect (Equation (52)). There is, however, another situation which also lets m app be real and independent of n. If there is a near-surface layer of a purely Newtonian liquid with increased viscosity (possibly also with an increased density), the contrast function in Equation (59) is again real and independent of n. In this case, the apparent acoustic thickness may be much smaller than the range with increased viscosity because the contrast function is less than unity. This situation is encountered in electrochemistry. The diffuse double layer is viscoelastic, in principle, but sometimes modeling it as a layer with increased Newtonian viscosity (no elasticity) matches the experimental data well [48]. This frequency shift should not be misinterpreted as the consequence of adsorption and desorption.
A side remark: If the near-surface viscosity is lower than the bulk viscosity, Equation (59) predicts a negative apparent Sauerbrey thickness. This situation is further discussed in Section 4.8.
Similar to the apparent mass, an apparent elastic shear compliance can be derived following Equation (54)  It is instructive to remind oneself of the differences between Equation (52) (thin film in a liquid) and Equation (46) (thin film in air):  For the thin film in air, the Sauerbrey mass is larger than the true mass. It is smaller for the thin film in a liquid (because of the missing-mass effect).  The viscoelastic correction scales as n 2 in air, while it scales as n in liquids (constant compliance assumed).
 In both cases, J f ' and J f '' are the coefficients to the viscoelastic correction. In air, J f ' enters the correction for Δf/n, while J f '' enters the correction for ΔΓ/n. The roles of J f ' and J f '' are reversed in liquids.  In air, the viscoelastic correction scales as the square of the film's mass because the film shears under its own inertia. Viscoelastic effects are only seen for films with a thickness of at least a few tens of nanometers. The film in a liquid is clamped from the other side. Viscoelastic effects are seen even for layers with a thickness corresponding to a few molecules.
In air, the 3rd-order perturbation analysis (Box 2) makes significant correction to the small-load approximation. In liquids, this difference exists, in principle, but it is negligible in practice.

Thick Layers
We now turn to adsorbates, which are comparable in thickness to the depth of penetration, δ. One might apply the multilayer formalism and depict the profiles, G (z) and ρ(z), as a sequence of many thin films with the wave being reflected at the numerous interfaces. It is easier to solve the wave equation directly. If the functions G'(z), G''(z), and ρ(z) are given, one may calculate the displacement profile, u (z), with one of the software packages, which numerically solve ordinary differential equations (Scipy being among them). One infers the complex frequency shift from u (z = 0) and du /dz(z = 0) as: The function u (z) is obtained by solving the wave equation, which is Note that the shear modulus appears inside the outer derivative because G (z) itself is a function of z. Appendix C.1 contains Python code solving this problem. The bottom of Figure 23 shows a velocity profile, v (z), obtained with this code. The top shows the shear modulus, G'(z) + iG''(z), which entered this calculation. Δf + iΓ f

Viscoelastic Dispersion and High-Frequency Rheology
Given that the QCM can determine shear stiffness, one wonders whether it can also do viscoelastic spectroscopy. The term "viscoelastic spectroscopy" here is synonymous to "high-frequency rheology". To cut a longer argument short, the QCM cannot explicitly do viscoelastic spectroscopy on thin films, even if one makes peace with the power laws from Equation (65) or Equation (66). It can, in principle, but the error bars are large.
The viscoelastic constants depend on frequency if the sample undergoes relaxations with rates comparable to the frequencies of excitation. (G and J then are complex for the same reason.) The QCM can yield spectroscopic information, but only to a limited extent. Firstly, the frequency range is only about one decade wide. Also, it is impossible to determine J' and J'' on every overtone separately, because this problem is underdetermined. As long as the thickness is not a priori known, 2n + 1 parameters would have to be derived from Δf and ΔΓ on n overtones.
At this point, one can exploit the narrow frequency range and the fact that rheological spectra usually are smooth. Rheological spectra are displayed on a logarithmic scale (Figure 24). Over a single decade, the frequency dependence of G f ' and G f '' can be approximated by power laws with power law exponents γ' and γ'': The subscript "cen" refers to a frequency in the center of the range accessible to the QCM. Typical is f cen = 30 MHz. J f '(ω) and J f ''(ω) are approximated as: The power law exponents in Equation (65) and Equation (66) differ. Worse, power laws in G f ' and G f '' do not turn into power laws in J f ' and J f '' after transformation from G to J with Equation (55). A set of power laws in G f ' and G f '' is not strictly equivalent to the corresponding set in J f ' and J f ''. Following from the Kramers-Kronig relations, β', β'', γ', and γ'' must be in certain ranges. If viscoelasticity is expressed in terms of compliance (J f ' and J f ''as in Equation (66)), one has 2 < β' < 0 and 1 < β'' < 1. If moduli are used (G f ' and G f ''), one has 0 < γ' < 2 and 1 < γ'' < 1. The software packages supplied by Biolin use the variables "µ" and "η". µ is equal to G f ' and η is equal to G f ''/ω. The power law exponent for η is between 2 and 0.
In the context of the QCM, power laws applied to the compliance, J , are closer to experiment because one of the corresponding power law exponents (only one, β') can be determined from experiment with good accuracy. It is often difficult to obtain a robust fit with all five fit parameters free (thickness, , tan(δ ), β', β''). Robust fits are obtained, though, if β'' is fixed. β' alone often is fitted with fair accuracy, because the film's mass, m f , enters the imaginary part of Equation (52) as a prefactor, only. (It enters the real part of Equation (52) as a prefactor and as an additive term.) Rather than fixing β'', one might also fix the difference between β' and β''. For polymers, the Rouse model and the Zimm model both predict β' » β'' in the high-frequency regime, see the right-hand side in Figure 24.
A side remark: For experiments in air, the situation is reversed. β'' rather than β' can be determined with good confidence from the curvature in Figure 17B. One may read a curvature from Figure 17A, as well, but the accuracy suffers from the unknown offset (proportional to the mass). Determination of both β' and β'' would be attractive because thin films are not easily studied with conventional rheology. Demanding equipment is needed [85 and others]. Again: the unknown layer thickness is the problem. If the layer thickness can be determined independently, this will help. The power law exponents give access to spectroscopic information. One of them (β') can be derived from the fits. With a model at hand (Zimm, Rouse, reptation, …) the value of β' can be interpreted.

Slip
"Slip" here denotes slip of a simple liquid at a solid wall. Slip in that sense is the exception rather than the rule. Often, the "no-slip condition" is a suitable boundary condition for liquid flows. In the following, slip does not denote sliding between solid surfaces and, also, does not denote wall slip of complex liquids, induced by shear thinning under large stress [108]. Slip in simple liquids is associated with a near-surface layer of reduced viscosity. The above statements are phrased in terms of the continuum picture. Arguably, a molecular description would be more appropriate. That would make no difference for the experiment.
Given constant shear stress, the shear gradient in the low-viscosity region is larger than the shear gradient in the bulk (Figure 25). When extrapolating the linear portion of the velocity profile to v(z) = 0, the intercept is negative. The negative intercept is the "slip length", b sl [109]. The reverse situation with a near-surface viscosity higher than the bulk viscosity causes a positive intercept, which may be viewed as the hydrodynamic thickness of the respective adsorbate. More quantitatively, the slip length is defined as: Slip of this kind might appear as exceptional because the density of a liquid close to a wall tends to be larger than the density of the bulk due to the attractive forces exerted by the wall. The question has caused much debate, but the evidence for slip in special situations has now solidified [110]. These situations include water at hydrophobic surfaces [111] and various flows inside hollow carbon nanotubes [112]. The slip length typically is of the order of a few diameters of the respective molecules, meaning, a few nanometers. Experimentally determining the slip length therefore is a challenge. At this point, the QCM plays out an advantage, which is the small depth of penetration of the shear wave. To the QCM, slip looks like an apparent negative Sauerbrey thickness. A Sauerbrey thickness (positive or negative) is easily determined with an accuracy of 1 nm or better. This being said: The slip length as determined with the QCM is different from b sl as defined in Equation (67). Following Equation (59), the negative Sauerbrey thickness (termed b sl,ac here, "ac" for acoustic) is: For the QCM, the density enters. The density of a slipping layer may be lower than the density of the bulk if slip is caused by nanobubbles or nanopancakes [113]. There are two more complications with nanobubbles:  Nanobubbles constitute a sample with lateral structure, while Equation (68) assumes lateral homogeneity.  This discussion ignores the surface energy of air-water interfaces (between the nanobubbles and the bulk liquid). Surface tension does play a role on the nanoscale. Surface tension turns nanobubbles into stiff objects [114,115]. (For macroscopic droplets or bubbles, the surface energy does not affect the resonance frequency because the associated oscillatory capillary pressure is small compared to the viscous stress.)

Point Contacts with Large Objects Clamped in Space by Inertia
By touching the resonator with a sharp tip, one increases the resonance frequency. This was first recognized by Dybwad in 1985 [116]. On an elementary level, the increase can be explained with the relation ω 0 » (κ R /m R ) 1/2 . When operated in the gravimetric mode, κ is about constant and the added mass lets the frequency decrease. When the resonator is touched with a tip (or with a sphere), the restoring force exerted by the contact lets the effective stiffness increase more strongly than the effective mass, thereby increasing the frequency.
The small-load approximation makes this understanding more quantitative. Within the parallel-plate model, the load impedance of a point contact is: n P is the number of contacts and A eff is the plate's effective area. The subscript elas stands for elastic loading (to be distinguished from inertial loading). The transverse stress is replaced by the transverse force, F , multiplied by the number density of the particles, n P /A eff . The contact is modeled as a Hookean spring (F = κ P u S = κ P v S /(iω)) with u S the displacement). The spring constant, κ P , can be complex, in which case the contact increases the bandwidth. The term κ P /iω is the mechanical equivalent of a capacitor's electrical impedance, which is 1/(iωC). Within Mindlin-theory, the contact stiffness of a sphere-plate contact is κ P = 2G *r C with G * an effective modulus, similar in magnitude to the shear modulus, and r C the contact radius [117].
The contact stiffness results from a small-scale deformation in the immediate vicinity of the contact. The deformation involves both the external object and the substrate (close to the contact). The relation κ P = 2G *r C requires a small contact area, that is, a "point contact". Contacts are small if the displacement pattern has spherical symmetry and decays as about 1/r 2 with r the distance from the contact. This requires r C to be much smaller than the local radius of curvature of the external object (often called R) and it also requires r C to be much smaller than the wavelength of sound, λ. In the opposite limit of r C ≫ λ ("sheet contact" in Figure 26) the resonator launches a plane wave into the external object and Δf should be analyzed with Equation (34). where Δf is positive and scales as 1/n, which is indeed observed [118]. The contact stiffness as derived with Equation (70) was compared to the expectations resulting from JKR theory in [119].
Equation (70) neglects effects of added mass. It also assumes that the object supporting the contact from the back is fixed in space. This assumption is reasonable for contacts with sufficiently large spheres (R » 200 µm). These are clamped by inertia.

Large Amplitudes, Partial Slip
Many piezoelectric devices can act as both sensors and actuators. The piezo effect works in both directions. The QCM, on the contrary, rarely is an actuator. Concerning static actuation, the piezoelectric coefficient of quartz is small compared to the competing materials such as lead zirconate titanate, PZT. One might still hope that the vibration would take some effect. A typical application would be sonolubrication [120][121][122]. Sliding of powders induced by a QCM vibrating at large amplitude has been studied experimentally in [123]. Generally speaking, kHz vibrations are better suited to sonolubrication than MHz waves [120].
Amplitude effects have been studied with the QCM on a few occasions. An incomplete list follows:  Unbinding of virus particles at high amplitudes was studied in [124].
 Adsorption was prevented at high amplitudes in [125] and other publications by the same group.  Cell adhesion as a function of amplitude was studied in [126]. Cell adhesion was delayed by high amplitudes, but cells, which had already adhered, did not detach when shaken vigorously.  High amplitudes can induce steady streaming, as shown in [127]. More generally, the Reynolds number at high amplitude can be large enough to let the nonlinear term in the Navier-Stokes equation (of the form ρ(v • ∇)v) be significant. This term can cause an oscillatory Bernoulli pressure. There may be a net attractive force onto colloidal particles, mediated by a high-frequency version of the Magnus force [128].
The following text is concerned with a quantitative discussion of nonlinearities in high-frequency contact mechanics as evidenced in QCM experiments. For linear resonators (also: "simple harmonic resonators", SHOs), the resonance frequency is independent of amplitude. A "linear resonator" is characterized by the restoring force being proportional to displacement. An example for a slightly nonlinear resonator is the pendulum. The restoring force is m R g sin(θ), where g = 9.81 m/s 2 is the earth's gravitational acceleration and θ is the angle of swing. If θ ≪ 1, one has sin(θ) » θ and the pendulum operates as a linear resonator. The "grandfather clock" (named that way in Wikipedia) employs a long pendulum in order to ensure small angles. The amplitude is constant. Should the amplitude become large, this will decrease the clock's frequency. The effect can also be observed, watching children on a swing.
Speaking of clocks, the current quartz clocks are slightly nonlinear, as well. The stress inside a deformed quartz plate is not strictly proportional to the strain. The deviations are small, but they are noticed when driving quartz clocks with too much electrical power. The problem carries the name "drive level dependence", DLD [129]. The DLD constrains the maximum voltage when interrogating the crystal's resonance properties. The nonlinearity is such that the resonance frequency increases quadratically with the driving voltage. The DLD is of considerable concern in time and frequency control. Actually, the nonlinear elasticity of crystalline quartz is only of concern in this context because it can only be measured based on the resonance frequency. This is an example, where frequencybased metrology plays out its strength.
A nonlinear force displacement relation may also originate from the sample. This is not usually the case with films, with liquids, or with adsorbates from the liquid phase because the amplitude is too small. Typical amplitudes are a few nanometers at most (Section 8.5.). With a penetration depth of 200 nm, the shear angle is less than 1%. Shear gradients of this magnitude typically are in the linear regime (stress proportional to strain). This is good news in some ways (no need to worry) and is a disappointment in others. Nonlinear rheology is of much interest in polymer science, but high-frequency rheology on polymer films using the QCM [12] is bound to be linear rheology.
Nonlinear force displacement relations are commonplace in contact mechanics [118,130]. In contact mechanics, the local stress at the points of contact is large. Also, the contact area can vary in response to the load. Nonlinear high-frequency contact mechanics can be studied with the QCM. The experiments are rather simple. Frequency and bandwidth are determined as a function of amplitude. Figure 27 sketches one particular mechanism (partial slip), which softens contacts at large amplitudes. Partial slip lets the resonance frequency decrease. Figure 27. At the edge of a contact under a transverse load (red in (A)) one often finds a sharp peak in transverse stress (B). If the local stress exceeds a certain critical stress, sliding sets in. This is referred to as partial slip [131]. Panel (C) sketches the force-displacement relation of a sphereplate contact with partial slip as calculated with the Cattaneo-Mindlin model [132]. The differential contact stiffness is indicated as a dashed line. It decreases with increasing displacement.
When rough surfaces make contact, the stress distribution is strongly heterogeneous. The stress is large at the tips of the asperities, giving rise to nonlinear force-displacement relations and even plastic flow, also called "asperity creep" [133]. Less well known is the heterogeneous stress distribution at contacts between smooth spheres and plates (more generally, at the edges of extended contacts between smooth surfaces). When these contacts are sheared or loaded in some other way, a stress singularity develops at the edge. A similar singularity exists at crack tips in fracture mechanics. The stress-under certain conditions-scales as r −1/2 with r the distance from the crack tip. The open wedge outside the sphere-plate contact in Figure 27B can be viewed as a crack.
When the resonator exerts a periodic transverse force onto a sphere-plate contact, the stress at the edge may be so high that the contact starts sliding at the edge. Partial slip is useful when heavy objects hit piles of gravel. Such piles rarely fall over because the energy of the impact is dissipated efficiently in partial slip (and, also, in gross slip, which is still local) [134]. Partial slip can be detrimental in engineering. Contacts, which have seen prolonged vibrations, may suddenly fail because partial slip has turned into what is called "fretting wear" [135].
Partial slip was modelled in the 1950s by Mindlin [131], building on earlier work by Cattaneo [136]. The calculation leads to a friction loop as shown in Figure 28B (a lensshaped loop, rather than an ellipse). The target is to turn this function into a prediction for Δf and ΔΓ as a function of oscillation amplitude. The oscillation amplitude will be called u S in the following (no hat, because its complex nature is unessential). One may guess that Δf and ΔΓ will decrease and increase with u S , respectively, because the force-displacement relation in Figure 27C bends downward and because energy is dissipated in sliding. The Cattaneo-Mindlin model will confirm that. Further, it will predict these changes to be proportional to amplitude. Some background is needed. If the relation between stress and displacement is nonlinear, the area average (always inherent to the small-load approximation) must be complemented by a time average, following: The derivation of Equation (71) makes use of the two-timing approximation [137]. As in Equation (69), stress was replaced by the transverse force acting onto the contact, F(t), normalized to area. F(t) is assumed to be periodic with the frequency of excitation, but not necessarily time-harmonic. If F(t) is time-harmonic, Equation (71) reduces to Equation (70). In the calculation of F(t) from the force-displacement relation, it is assumed that the motion of the substrate is close to time-harmonic (of the form cos(ωt)). This type of displacement control (rather than force control or mixed control) is inherent to the smallload approximation.
Relations similar to Equation (71) are exploited in scanning force microscopy [138]. In that context, the function F(z) (with z the distance between the tip and the surface) can be explicitly reconstructed from the resonance properties of the cantilever as a function z. Such an explicit reconstruction of a force-displacement relation is not possible for partial slip. Partial slip is hysteretic; hysteresis violates the conditions needed for the explicit reconstruction. Again, the force-displacement relation cannot be explicitly derived from Δf (u S ), but certain models can be formulated and can be used to predict the functions Δf(u S ) and ΔΓ(u S ), which can be checked against experiment. Figure 28 sketches three such models. For the viscoelastic contact ( Figure 28A), Δf and ΔΓ are independent of amplitude. Figure  28B shows the prediction from the Cattaneo-Mindlin model for partial slip [131]. (There are other models of partial slip, one of them described in [139].) Figure 28C depicts the transition to gross slip [140]. When gross slip sets in, ΔΓ decreases at large amplitudes and Δf levels off to a small value.  [141], not by a quantitative model. ΔΓ decreases at large amplitude because the friction force in steady sliding weakly depends on velocity. Dividing by u S 2 causes ΔΓ to decrease with .
The following discussion is concerned with partial slip ( Figure 28B). In order to calculate Δf and ΔΓ from Equation (71), the function F(t) must be derived from the functions F(u N , u S ,ω) and F(u N , u S , ω) as predicted by the Cattaneo-Mindlin model. (For their algebraic form, see, for instance, [132].) The parameter u N is the displacement normalized to the peak displacement, u S . The subscripts  and  denote increasing and decreasing uN. Because these two forces are different, the force-displacement relation forms a loop.

Area inside loop
Peak force The term u N 1  u N 2 1/2 in line 1 can be viewed as statistical weight. Δf is proportional to a weighted average of |F + F|, where the weight function has a sharp peak at the turning point. Following this argument, Δf is roughly proportional to the force-displacement ratio at maximum displacement (blue dots in the top row in Figure 28). Following line 2, ΔΓ is proportional to the area inside the friction loop, divided by . (One divides by u S 2 because du N in Equation (72) is equal to d(u/u S ).) The latter result is exact, there are no approximations involved. For partial slip following Cattaneo and Mindlin, the calculation leads to [141]: κ P is the contact stiffness in the low-amplitude limit. FN is the normal force. The model's one free parameter is µ (the "friction coefficient"), which is the ratio of the critical tangential stress for sliding to the normal stress. µ should not be naively identified with the friction coefficient in macroscopic sliding. It turns out, though, that µ is of order unity, similar to the conventional friction coefficient [119]. Equation (73) was found to apply in some experiments [141,142]. Others rather show a quadratic dependence, which actually was predicted by another model of partial slip, put forward by Savkoor [139]. Still other experiments (mostly on granular media) show an increase of Δf with u S , which can be explained with shear stiffening.
There is a caveat: Equation (73) results from an integration over the friction loop. All models making use of Equation (71) involve such an integration. The shape of the friction loop cannot be inferred from the dependence of Δf and ΔΓ on amplitude. Also, the forcedisplacement relation leading to Equation (73) was derived assuming a quasistatic situation. The dynamics at MHz frequencies might be different. Even a response, which is linear on the sub-µs time scale and therefore leads to an elliptical friction loop ( Figure 28A), can produce an amplitude dependence of Δf and ΔΓ if the width and angle of the ellipse depend on amplitude.

Structured Samples, Numerical Calculations
Many samples of interest in soft matter physics have some in-plane structure ( Figure  29). This includes proteins [143], dendrimers [144], biological cells [145], and colloidal particles [146,147]. Adsorbed vesicles, which may or may not rupture and flatten out into supported lipid bilayers (SLBs), have been studied in much detail [148][149][150]. Predicting shifts in frequency and bandwidth induced by such samples from the structure and the viscoelastic parameters requires a numerical model. As long as the area average in Equation (25) can be applied, the shifts in frequency and bandwidth can be computed numerically. The input to the small-load approximation is the area-averaged amplitude of the periodic stress at the resonator surface. The stress field at the surface can be extracted from a solution of the equations of continuum viscoelasticity.
The medium obeys the Stokes equation expressed in the frequency domain, which is: where ν = η /ρ is the kinematic viscosity. The pressure, p , follows from the modulus of compression, K , and the divergence of the velocity field as: The density was assumed as constant in Equation (74), but may also depend on position. The complex dynamic viscosity, η (r), can take any value. Elastic objects have η'' ≫ η'. Rigid objects are represented as objects with large |η | or, equivalently, with small η . Such objects are hardly deformed by the shear stress. For rigid objects, the question of viscous or elastic response does not occur, because it does not matter whether a negligible deformation occurs in-phase or out-of-phase to the stress. This argument reiterates the previous statement that the nontrivial samples to the QCM are the soft samples.
The Stokes equation is a linearization of the Navier-Stokes equation. The Reynolds number is assumed to be so small that the nonlinear term (ρ(v · ∇)v ) in the Navier-Stokes equation can be neglected. Actually, effects of finite Reynolds number can be seen in QCM experiments at high amplitude. Because these effects are weak, they can be modeled, based on the algorithms described below. Call the solution to the linear problem v (0) . The nonlinear term then generates a Bernoulli pressure of the form ρ v (0) · ∇ v (0) . This pressure vanishes for pure shear flow because it contains the dot-product of the velocity and its gradients. It does not vanish for structured samples. It occurs at = 0 (steady streaming) and at 2 (2nd harmonic generation) because of the relation cos 2 (ωt) = 1/2 1 + cos(2ωt) . The Bernoulli pressure drives a weak 2nd-order flow, which can be computed from v (0) . Figure 30 sketches a simulation volume. It is a few tens of nanometers wide and contains a few adsorbed particles (if this is the problem of interest). Periodic boundary con- ditions ("b.c.") apply at the side walls. A Dirichlet boundary condition applies at the bottom (v (z = 0) = (v S ,0,0)). The boundary condition at the top should be an impedance boundary condition (also: "Robin boundary condition" [151]), meaning that the ratio of the velocity gradient to the velocity should be such that the stress-velocity ratio inside the boundary is equal to the wave impedance outside the boundary. There is ongoing research turning these concepts into practical and efficient algorithms. Computational techniques applicable to this problem include the finite element method (FEM) [152,153], the finite volume method (FVM) [154], and the Lattice-Boltzmann method (LBM) [155,156].
If the sample is much thinner than the wavelength of sound, it may look to the QCM like a film. The QCM then does not actually recognize the structure. One may infer an apparent mass, m app , and an apparent compliance (J app ', J app '') from plots of Δf/n and ΔΓ/n versus n (as in Figure 21). The problem to the experimentalist is to interpret these parameters. For instance, the apparent mass is different from the true mass because a certain amount of liquid, which is trapped in the space between the particles, takes part in the resonator's motion [157][158][159]. Likewise, the apparent compliance has a contribution from the flow of the liquid around the particles.
The numerical simulation calculates Δf and ΔΓ from a known structure. The task to the experimentalist usually is the reverse problem, namely inferring the structure from experimentally determined values {Δf/n, ΔΓ/n} is another. The reverse problem is underdetermined, in general. The formalism can only run backwards, if the structure is known to a significant extent with few free parameters, so that the experimental values of {Δf/n, ΔΓ/n} can be exploited to determine these parameters. Such a situation is (for example) encountered, when rigid spheres of known size adsorb to the QCM surface, following random sequential adsorption [144]. If coverage can be estimated independently, the only free parameter of this problem is the stiffness of the sphere-plate contacts. Arguably, contact stiffness involves a total of four free parameters, which are the real and the imaginary parts of the shear stiffness, κ Sh , and the bending stiffness, κ b (related to translation and rotation of the sphere). For spheres in air, the problem is amenable to an analytical treatment (Section 6.2.). In liquids, hydrodynamic comes into play. For this problem, a calibration curve relating contact stiffness to the sets of {Δf/n, ΔΓ/n} can be obtained from simulation. These calibration curves can be used to analyze experiments. A related problem is the fractional trapped mass as a function of coverage. The trapped mass may also be estimated from calibration curves obtained with simulations on structures, which are welldefined, on the one hand, but still reasonably close to experiment, on the other.
Interface to bulk, impedance b.c. side walls periodic b.c. substrate, v S = 1

Roughness
Roughness in QCM experiments is always a reason to worry. Gold surfaces prepared by physical vapor deposition (PVD) have an rms roughness of about 1 nm. Spin-cast polymer films usually are considered as smooth. One might model roughness effects numerically as in Section 5.3., but there is a wide range of possible geometries.
We limit the discussion to small-scale, shallow roughness as modelled analytically in [160,161]. The vertical scale of roughness is assumed to be smaller than the lateral scale, which is realistic for gold surfaces prepared by PVD. The lateral scale is assumed to be smaller than δ, which again often is realistic. The model has two free parameters, which are the vertical scale, h r , and the lateral scale, l r . One might also use the vertical scale, h r , and the aspect ratio, h r /l r . The aspect ratio is assumed to be smaller than unity ("shallow roughness").
There is an interesting experimental statement in the literature, which is that the bandwidth is less affected by small-scale roughness than the frequency [162]. This finding is corroborated by the analytical treatment of shallow roughness following [161]. These authors Fourier-decompose the height profiles into sinusoidal corrugation waves with wave vector q. They solve the hydrodynamics problem for the different Fourier components, separately, and calculate Δf and ΔΓ. The total frequency shift follows from integration over all wave vectors, where the weight function is the square of the respective amplitude. A Gaussian distribution is assumed, the center of which is much larger than the inverse penetration depth ("small-scale roughness"). The following equations are obtained: Equation (76) can be reorganized as: The term in square brackets modifies the Gordon-Kanazawa term and it does so in proportion to h r 2 . This term affects both frequency and bandwidth. The second term only affects the frequency shift. It can be thought of as a Sauerbrey-like term, covering trapped mass. In principle, this terms scale as h r 2 , as well. However, if the aspect ratio is constant, one factor of h r is absorbed in the aspect ratio (given as h r /l r ), which leaves the other h r as a linear term. Roughness is often created in a way, which leaves the aspect ratio constant. The amount of trapped mass then depends linearly on the cluster size. The authors of [162] did experiments of that kind. Because effects in bandwidth are proportional to h r 2 , they are hardly seen at small roughness.

The Sphere with Moderate Mass
One example of a coupled resonance (the film resonance) was already discussed in Section 4.5.4. Here, we start from another example, which is the sphere attached to the resonator. (Particle adsorption to a QCM was discussed as early as 1971 in [163].) Contrasting to Section 5.1., the sphere has moderate size. It is neither clamped in space by inertia (elastic loading), nor is it a nanoparticle in the sense, that it would be rigidly attached to the resonator and constitute a Sauerbrey load. Figure 31 shows a lumped-element circuit. The link between the particle and the resonator is depicted as a spring with stiffness κ CR . A dashpot with drag coefficient ξ CR was placed in parallel to the spring, accounting for dissipative processes. The subscript CR denotes the coupled resonance. Written that way, the circuit suggests that the spring constant and the drag coefficient were independent of frequency. They may well depend on frequency. The spring constant may also be written as a frequency-dependent complex parameter κ CR (ω) = κ CR '(ω) + iκ CR ''(ω) . κ CR (ω) then is a response function similar to G (ω). κ CR (ω) has different units, though (force/displacement rather than stress/strain).

Figure 31.
A lumped-element representation of a load giving rise to a coupled resonance. This mechanical circuit implicitly invokes the small-load approximation insofar, as the resonator has been drawn as a plate. It might also have been drawn as two large spheres, coupled to each other by another spring ( Figure 1B). Had it been drawn as in Figure 1B, the solution to the full set of dynamical equations would have predicted anticrossing (as for the film resonance, Section 4.5.4.).
Anticrossing is not further discussed in the following. Coupled resonances caused by particles will usually be too broad to let anticrossing be visible.
The elements in Figure 31 represent mechanical impedances, which are ratios of force to velocity (rather than stress to velocity). Mechanical impedances and the electromechanical analogy [164] are discussed in Box 3. In mechanics, the total impedance of two elements arranged in parallel is the sum of the two individual impedances, following the "mechanical Kirchhoff rules". Application of these rules to the circuit in Figure 31 leads to: where F is the transverse force exerted by the sphere onto the resonator surface. Inserting Equation (78) into Equation (25) leads to: where n p is the number of particles.
A side remark: One might be tempted to write ω CR 2  ω 2 as (ω CR + ω)(ω CR  ω) and (ω CR + ω) » 2ω CR as in the mathematics leading to Equation (7). This is problematic because the coupled resonance is not usually a narrow resonance. For the same reason, polar diagrams of the complex frequency shift as in Figure 32 usually show spirals (as opposed to circles, Figure 3B). Equation (79) contains the Sauerbrey case and the elastic-load case in the limits of ω ≪ ω CR ′ and ω ≫ ω CR ′, respectively ( Figure 33). When ω ≪ ω CR ′, the right-hand side in Equation (79) is about n P ωm CR /(A eff πZ q ), which is equivalent to the Sauerbrey equation.
There actually is a frequency range, in which the response is close to Sauerbrey-like, but in which the dynamics of the sphere is still seen as a small deviation from Sauerbrey behavior (compare to Figure 16 and Section 4.5.3.). If ω is less but not much less than ω CR ' Δf / deviates from Sauerbrey-type behavior. This becomes evident when expanding Equation (79) to 3rd order in ω as: The parameters m CR and κ CR will need interpretation because the motion of the particle combines translation and rotation. m CR and κ CR are a modal mass and a modal stiffness. Section 6.2. elaborates on that problem.
The deviations from Sauerbrey behavior are proportional to the softness of the contact (the inverse contact stiffness). In this regard, Equation (80) is the analog to Equation (45), where the latter equation describes a soft film, rather than a soft link to a particle.
When the frequency of the coupled resonance is in the range covered by the QCM, Δf crosses from negative to positive [165,166]. An example is shown in Figure 32, adapted from [167]. The resonators were in contact with silica spheres with a radius of 2.5 µm. The spheres were attracted to the surface by gravity and by van-der-Waals forces. Both the spheres and the resonator surface carried negative charge, resulting in an electrostatic repulsion, which competes with the attractive forces. Variation of ion strength as indicated in the legend tuned this repulsive force [168]. Adding salt has two consequences. Firstly, the frequency of zero crossing increases. This can be understood as the consequence of increased contact stiffness, which in turn is the consequence of reduced electrostatic repulsion. Also, the radius of the semi-circle (more precisely, the semi-spiral) decreases. A closer look at Equation (79) shows that this radius is proportional to the ratio κ CR '/ξ CR . Why electrostatic screening decreases this ratio is nontrivial.   Figure 33. The frequency shift induced by a coupled resonance, following Equation (79). At low and high ω, the Sauerbrey limit and the elastic-load limit are recovered. The overtone scaling in these limits is Δf/n » const and +Δf·n » const. In a range next to the Sauerbrey regime, the softness of the contact can be inferred from a non-trivial the overtone dependence of Δf/n + iΔΓ/n.

Box 3. The electromechanical analogy.
Differing from electricity, mechanical impedances are additive when two elements are arranged in parallel. Inverse impedances are additive, when the elements are arranged in series. These are the mechanical Kirchhoff rules. In electricity, the reverse rules apply. Current (the analog of velocity) is additive, when two resistors are placed in parallel. For dashpots arranged in parallel, the forces (the analogs of voltage) are additive. Ground U = 0 Open end F = 0 Open end I = 0 Wall V = 0

Influence of Rotation on the Frequency Shift
Incorporating rotation into the formalism is worth the effort. Similar work was reported by Tarnapolsky et al. in [169]. Differing from the text below, these authors were concerned with large spheres (bacteria) in a liquid environment. The text below addresses small spheres (ω CR ' ≳ ω) and avoids hydrodynamics by sticking to a dry environment.
The sphere has two degrees of freedom, which are translation (with velocity u P ) and rotation about the sphere center (with rate Ω ). Inertia is balanced against the force and the torque exerted by the contact, F C and M C , by the relations: The inertial terms are iωm P u P for translation and iωm P (2/5)R 2 Ω for rotation. The subscript P denotes the particle. m P (2/5)R 2 is the moment of inertia of a solid sphere rotating about its center. The transverse restoring force exerted by the contact is F C = κ Sh /(iω) u P  RΩ  u S . κ Sh is the shear stiffness of the contact. Following Mindlin theory, κ Sh is given as 2G *r C where G * is an effective modulus and r C is the contact radius. The term (u P -RΩ )u S is the difference between the velocity of the sphere and the velocity of the substrate at the point of contact (Figure 34). The torque has two components. The first component is F C R. The second component follows from the contact's bending stiffness. Following [169], we write the torque as M C =-κ b /(iω)R 2 Ω , where κ b is the bending stiffness. Defined this way, κ b has the same units as κ Sh . Defined this way, κ b is not a ratio of torque to angle (as elsewhere, the ratio of torque to angle here is κ b R 2 ). Following Dominik and Tielens [170], the bending stiffness is: With these relations, the force and the torque are: Inserting Equation (83) into (81) and reorganizing terms, the following equation system results: The terms in brackets in Equation (84) form a matrix, which must be inverted to obtain the vector ( P ,Ω ) from the source term (κ Sh (iω) ⁄ u S , κ Sh R (iω) ⁄ u S ). A resonance occurs, when (u ,Ω ) is nonzero even in the absence of a source term. The resonance frequencies are found by setting determinant of the coefficient matrix to zero. There are two resonance frequencies, because there are two dynamic variables. They are given as: Both modes combine translation and rotation, but one is predominantly translational, while the other is predominantly rotational. Note that the two resonance frequencies are not vastly different unless r C ≪ R. A search for two coupled resonances in experimental QCM data is difficult because the experiment must involve an ensemble of spheres and because there will be a distribution in the parameters κ Sh and r C . The heterogeneous line broadening caused by the distribution of parameters often lets the two resonances merge.
The above remarks concerned situations, in which the frequency of the resonator is close to one of the frequencies of the coupled resonances. In general, the particles will be driven off-resonance. Setting u S on the right-hand side in Equation (84) to unity fixes u and Ω (to off-resonance values). The force, F C , follows from Equation (83). Inserting F C into the small-load approximation yields Δf and ΔΓ. Plotting Δf and ΔΓ versus frequency, yields two resonances. This graph is not shown here because it looks very similar to Figure  36, which was obtained with a numerical simulation.
In the following, we are concerned with small spheres, which are still large enough to show some dynamics. These spheres are almost Sauerbrey loads (are almost rigidly attached to the plate). Similar to Equation (80), we expand Δf to 3rd order in ω. In Equation (80), the coefficient to ω 2 was named m CR /κ CR and it was left open, what these parameters mean. This gap can now be closed. The 3rd order expansion in ω of Δf calculated from F C leads to: Clearly, m CR is equal to m P (in this limit). With regard to contact stiffness, Equation (86) shows that κ CR 1 from Equation (80) follows κ CR 1 = κ Sh 1 + κ b 1 (as also noted in [169]). Figures 35 and 36 illustrate the situation at hand of an FEM simulation. The simulation modeled a cylinder (rather than a sphere) because the simulation occurred in two dimensions. Also, the environment was a liquid. In agreement with Equation (85), two coupled resonances are found. For the higher frequency, the axis of rotation is close to the center of the sphere. For the lower frequency, the sphere rotates about the point of contact. Such a rotation is equivalent to a superposition of a rotation about the sphere center and a translation.
Δf , the sphere rotates about the contact, in (C) it rotates about its center. Adapted from [171]. Figure 36 shows Δf and ΔΓ versus frequency as derived from the FEM calculation. Because this is a simulation, the frequency is not limited to the odd integers of the fundamental. One observes two coupled resonances, one stronger than the other. Figure 36 looks similar to an absorption spectrum form IR spectroscopy, which is no coincidence [171]. Molecules can vibrate in different ways, just like the sphere on the plate may combine translation and rotation in different ways. When probed with IR waves, some vibrational modes (like the carbonyl stretching mode) show prominently in the spectrum, while others are weaker. Some vibrations are not seen in IR light at all. These are "forbidden", meaning that the vibration is not linked to an oscillating dipole. The vibrations of homonuclear diatomic molecules are forbidden. They are only seen in Raman spectroscopy. They are not seen in IR spectroscopy because the vibration "does not couple" to the electromagnetic wave. Some vibrations are "weakly allowed". They appear to be not coupled to IR light, at first glance, but there is a loophole. The overtone vibrations in IR spectroscopy are forbidden in the frame of the harmonic resonator. The selection rule is Δv = ±1 (with v the vibration quantum number). Anharmonicity weakens the selection rule. The overtones are seen, albeit with a small integrated line strength.
The coupled resonance to the left in Figure 36 is weakly allowed in that sense. The vibration mostly is a rotation about the point of contact. The restoring force mostly amounts to a torque. Torques are not seen by the QCM. The QCM only sees transverse forces. If the axis of rotation would intersect the point of contact, exactly, this coupled resonance would be "forbidden". The vibration would not couple to the substrate's transverse motion. Looking at the mode a little closer ( Figure 35B), one finds that the axis of rotation is slightly displaced from the contact. For that reason, the vibration does exert a small transverse force and is seen by the QCM, although with a small line strength. (Coupling across the liquid phase may also play a role.)  Figure 35. There are two separate coupled resonances. These are the consequence of two dynamical variables, which are translation and rotation. Translation and rotation are always coupled, but there are two special linear combinations of translation and rotation, which both produce a resonance. Only the transverse force at the point of contact shifts Δf and ΔΓ, a torque does not take an influence. The relative importance of the force compared to the torque is larger for the high-frequency mode. The sphere rotating about the point of contact mostly exerts a torque. This low-frequency mode therefore appears as less pronounced when probed with the QCM. Adapted from [171].

Other Types of Coupled Resonances
This text elaborates in some depth on coupled resonances caused by films and coupled resonances caused by particles. More generally, a "coupled resonance" can denote any load, which lets Δf be described by Equation (79). Two other noteworthy examples are the following:


Standing compressional waves (Section 8.1.) give rise to coupled resonances, when the distance to the opposite cell wall is an integer multiple of half the wavelength. At these distances, the compressional wave is a standing wave and the damping is large. This phenomenon can be exploited to check for the magnitude of compressional wave effects. The experiment is simple. One lets the water level in an open cell slowly decrease by evaporation overnight. Figure 7 in [172] shows data of this kind. In this example, the compressional wave effects were much stronger on the fundamental than at 15 MHz. This is a general rule and one of the reasons, why data from the fundamental often are discarded from the analysis.  The vibration of interest may couple to other modes of vibration of the crystal, where the exact mechanism of coupling is unclear and where even the nature of the other mode is unclear. These so-called "activity dips", which sometimes occur when ramping temperature up or down, can be a problem in time and frequency control [173]. An activity dip lets the bandwidth increase at a certain temperature and lets the frequency go through a corresponding antisymmetric pattern. Activity dips are not discussed further here, even though they are occasionally seen in sensing, mostly during temperature sweeps.

Piezoelectric Stiffening
Piezoelectric stiffening is an advantage in time and frequency control. It is closer to a problem in sensing. Piezoelectric stiffening also is instructing from a fundamental point of view.
When a bare piezoelectric plate is sheared, an electrical polarization is created in addition to the strain. The electrical energy contained in the polarization contributes to the overall energy of the sheared plate, thereby increasing the plate's stiffness. When the crystal surfaces are plated with electrodes and when these electrodes are short-circuited, the piezoelectrically induced polarization is compensated by a corresponding charge in the electrodes. Shearing the crystal then becomes easier because no electrical work is done. The spring constant decreases and so does the resonance frequency because of ω 0 » (κ R /m R ) 1/2 . Situations between those two limits are created by connecting the electrodes across a capacitor or some other electrical impedance. The situation is reminiscent of the difference between the two heat capacities of a gas, Cp and CV. Cp and CV are determined at constant pressure and constant volume, respectively. Cp is larger than CV because the volume expansion, which occurs when pressure is kept constant, amounts to a mechanical work. The open-circuit stiffness of the piezoelectric plate and the stiffness with short-circuited electrodes correspond to boundary conditions of constant charge (total charge in the bulk and at the surface) and of constant electrical potential, respectively. If the voltage between the two electrodes is zero, the creation of charge at the surfaces is not linked to an electrical work.
One can include a voltage-tunable capacitor (a varicap diode) into the circuitry and tune the resonator's frequency this way. This device is the voltage-controlled crystal oscillator (VCXO [52] or VCO). An approximative relation connects the frequency shift to the external capacitance, C ext , as: C 0 is the resonator's electrical capacitance (the "parallel capacitance") and C 1 is the motional capacitance (Section 8.4.). The subscript PE stands for piezoelectric stiffening, ∞ stands for infinite C ext . An infinite external capacitance amounts to short-circuit electrodes (constant potential). Note that C ext = ∞ is not a practical reference state because the current is then entirely absorbed in C ext and does not cause a vibration of the crystal.
C 0 and C 1 are related as: k t is the electromechanical coupling coefficient, e 26 = 9.65 × 10 2 C/m 2 is the piezoelectric stress coefficient, ε q = 4.54 is the dielectric constant, and G q = 29 × 10 9 Pa is the shear modulus. For AT-cut quartz, k t 2 is about 0.8%. Piezoelectric coupling is not particularly strong for quartz. The tuning range of quartz resonators is smaller than the tuning range of resonators made from competing materials such as langasite [174]. Equation (88) distinguishes between an "ideal" coupling coefficient (derived from the material constants) and an "effective" coupling coefficient (derived C 0 and C 1 ). The two are different because the geometry (electrode shape, energy trapping, …) takes an influence. Inserting numbers into Equation (87), one finds the pulling range to be around 1 kHz. In a sensing context, all kinds of stray capacitances (more generally, all kinds of electrical boundary conditions) enter C ext . For sensing, piezoelectric stiffening mostly is an annoyance. In particular, the cables must not move while the experiment is running. If they do, this will affect the capacitance seen by the crystal. A second problem are electric fields permeating the sample from the surface of the resonator. One avoids that by making the front electrode larger than the back electrode and by grounding the front electrode well. This problem also affects electrode-less resonators [175].
For the sake of formal consistency, we formulate a modified version of Equation (88), which lets the sample's electrical impedance look like any other load entering the smallload approximation. This equation is: where Z ext is the sample's electrical impedance. The term in square brackets converts between an electrical and an acoustic impedance (Section 8.4.). The hypothetical reference state with short-circuited electrodes was replaced by a more realistic reference state with some external electrical impedance, Z ext,ref . For more details see Chapter 5 in [5].
In principle, piezoelectric stiffening might provide for a scheme to measure the sample's electrical impedance in addition to its shear-wave impedance. That has turned out to be difficult [176,177]. In practice, one will usually determine electrical impedances with electrical equipment (as in electrochemical impedance spectroscopy, EIS, or, more generally, electrical impedance spectroscopy).

Energy Trapping, Compressional Waves
The limitations of the laterally infinite parallel plate as a model for the QCM come in two forms. Firstly, the edges take an influence on the mode of vibration for AT-cut quartz. Secondly and more importantly, the practical resonators do not actually have parallel surfaces, at least in acoustic terms ( Figure 37). In order to mount the resonator between O-rings with little damping, the resonator's vibration amplitude at the edge must be as small as possible. That is achieved by making the resonator slightly thicker in the center than at the edge [178,179]. The mechanism is called energy trapping. The resonator may be viewed as an acoustic lens. The concave surfaces focus the acoustic energy to the center. The situation has been analyzed with analytical theory in considerable depth [180,181].
Resonators designed for use on the fundamental indeed have concave surfaces. The resonators used in sensing mostly achieve energy trapping with keyhole-shaped electrodes. Often the back electrode is smaller and thicker than the front electrode, so that the back electrode defines the amplitude distribution. If convex surfaces are employed, the amplitude distribution is similar to a Gaussian (circular or elliptical) [182]. With key-hole electrodes, there is small-scale variability (bottom in Figure 37). The displacement pattern can be visualized in a few different ways, none of them being simple [172,[183][184][185]. These images reveal rather complex patterns with a considerable amount of irregularity [186,187]. In particular, the patterns may deviate from what symmetry would dictate in the absence of crystal defects. Figure 37. Energy trapping allows to mount resonators at the edge with little damping. When the resonator is thicker in the center than at the edge, the shear vibration is confined to the center. Ideally, the amplitude distribution is smooth (top). Experiments often evidence small-scale variability (bottom, [172,186]). u 2 (r) is the weight function to be applied in area averaging following Equation (26).
Energy trapping has two important consequences, which are an increased resonance frequency and flexural contributions to the mode of vibration.
Why would energy trapping increase the frequency? The apparent stiffness of a vibrating body depends not only on the material's elastic modulus but also on the steepness u(r) u(r) O-rings of the gradients in displacement. For instance, the effective spring constant of the parallel plate is given as: High overtones have a large effective spring constant because the shear gradients are strong. The more deformation is contained in a mode, the higher is the energy in the elastic deformation.
Increased energy trapping increases the steepness of the in-plane gradients, thereby increasing frequency. This effect is noticed when comparing the resonance frequencies of the different overtones. These do not occur at the exact integer multiples of the fundamental because of energy trapping (and, also, because of piezoelectric stiffening, but the latter influence is smaller) [188]. The influence, which energy trapping takes on the resonance frequency, is problematic when a sample contacts the resonator at the center, only, because the sample then improves the energy trapping, thereby increasing frequency (Figure 38). This effect is superimposed onto the sample's load impedance (see the discussion below Equation (34)).
A second consequence of energy trapping are flexural contributions to the displacement pattern. Because the amplitude of shear varies between the center and the edge, the resonator bends ( Figure 39A-C). On the high overtones, bending is reduced for the reasons discussed in [172]. Because bending is particularly strong on the fundamental, data from the fundamental often are discarded from the analysis. Bending is reduced in a liquid environment because the liquid itself is compressed in the regions with in-plane gradients of the transverse displacement. It exerts a pressure onto the plate, reducing the flexural displacement (green in Figure 39B).  Figure 39. (A,B) When a parallel plate experiences a hypothetical pure thickness-shear deformation with in-plane gradients, volume is not conserved at those places, where the in-plane gradient is large. The pure thickness-shear mode therefore is not realized: The plate bends. In liquids, the bending is reduced because the liquid is compressed as well. (C) A map of the vertical displacement as determined with laser doppler velocimetry (LDV). Adapted from [172]. (D) Δf and ΔΓ determined in a temperature ramp. The opposite wall of the liquid cell was parallel to the resonator surface, which generated standing compressional waves with λ depending on temperature.
Flexural motion launches compressional waves and standing compressional waves can cause coupled resonances. Figure 39D shows an example (from [189], see also [190,191]). The crystal and the wall of the liquid cell opposite to the crystal form a cavity for compressional ultrasound. The cavity resonates when the distance between the two surfaces is equal to nλcomp/2 with n an integer. The wavelength and the dimensions of the cell vary with temperature, which causes the nightmare shown in Figure 39D. In order to not let that happen, one may design the cell such that the resonator surface is inclined relative to the opposite wall. That does not remove the compressional waves, it only avoids the coupled resonances. This recommendation was actually formulated as early as 1987 by Eggers and Funk in the same paper, in which they propose to analyze QCM data in terms of the shift of the complex resonance frequency [33]. Eggers and Funk placed a "spoiler" (an irregular piece of Teflon) in front of the crystal in order to deflect the compressional waves. When working with open cells, one may place a paper clip onto the airwater interface (held there by surface tension). Few researchers have reported on temperature sweeps employing liquid cells. It is always difficult to say why certain experiments were not done or not reported, but one may guess that compressional waves were a problem.
Effects of compressional-waves can never be eliminated in measurements of viscosity because one compares frequencies measured in air to frequencies measured in liquid. With regard to adsorption experiments in liquid, one can hope that the effects of compressional waves remain constant during experiment and therefore disappear from the frequency shift.
A historical note: Soft matter was studied with acoustic resonators in the 1930s to 1950s, already. The frequencies mostly were in the kHz range. Torsional resonators and reflection devices [75] were used, mostly. Thickness-shear resonators were available at that time as clocks. Mason and McSkimin did use these when they expanded their frequency range to beyond 1 MHz, but they used them in an indirect way. In 1949 they glued AT-cut crystals to the ends of cylindrical rods of fused quartz. The crystals were excited by a radio pulse, thereby launching a transverse acoustic wave. The wave travelled down the cylinder, was reflected at the other end (at the interface to the liquid under study, see the remarks around Equation (37)) and returned to the crystal, which transduced the reflected wave's shear stress to a voltage. The liquid's viscosity was inferred from the reflected amplitude. In the introduction, the authors discuss using the plates as such: "Consideration was given to the use of a thickness vibrating shear crystal of the AT or BT type, but it was found that the shear motion was too closely coupled to other modes of motion, such as flexure 20 40 1 kHz DG Df, DG

Temperature [°C]
Df D C modes, to give reliable results. Hence another method had to be used." In hindsight, these problems turned out to be less severe than the authors had believed.
With some diligence, crystals not employing energy trapping can be mounted such that they are only weakly damped by the O-rings. Such crystals have electrodes covering the entire area rather than key-hole electrodes. However, such resonators immersed in liquids still display compressional-wave effects. This is evidenced with open cells, the water level of which slowly decreases due to evaporation [172]. Coupled resonances as shown in Figure 39D are still seen. Evidently, the edges of the resonator alone give rise flexure modes because of the anisotropy of the elastic constants.
Experiments in the dry can occur even with no electrodes at all. An electrodeless resonator blank can be placed on a rough metal surface. The other electrode can be placed above the resonator, leaving an air gap. The rough surface supports the blank across small asperities, which hardly dampen the resonance [81]. The problem with this arrangement is that the resonance frequency depends on the width of the air gap. The upper electrode must be mounted rigidly. Once this problem is solved, the experiment is rather clean and well-defined.

Anharmonic Sidebands
Anharmonic sidebands are solutions to the acoustic boundary value problem, which contain nodal planes perpendicular to the surface. Figure 40 shows examples. The maps of the displacement amplitude were produced by Sauerbrey, based on a somewhat intricate optical method [192]. Anharmonic sidebands have been exploited for sensing [193], but the examples are scarce. The important modes for sensing do not have nodal planes perpendicular to the plate surface. They might be called n-00-modes, where n is the number of nodal planes parallel to the surface and the two other indices count the number of radial and elliptical nodal planes. Figure 40. Vibrational patterns of the 1-00-mode ("1") and a few anharmonic sidebands. The top shows |Z el | versus frequency. The mode shape was determined with a variant of optical reflectometry. Adapted from [183].
Anharmonic sidebands are a problem if they overlap with the n-00-modes. They then couple to these modes with detrimental consequences. By and large, the 13-00-mode, shown in Figure 41 is useless for sensing. A resonance curve can certainly be fitted to the admittance trace, but Δf and ΔΓ as derived from these fits are bound to vary erratically in experiment. Overlap with anharmonic sidebands is much less of a problem in air than in liquid because of the smaller bandwidth. In air, overtone orders up to 19 can be evaluated.

Towards 3D-Modelling: The Small-Load Approximation in Tensor Form
The following section borrows from quantum mechanics. In quantum mechanics, a small perturbation (for instance caused by a magnetic field or by neighboring molecules) is sometimes superimposed onto a stronger, unperturbed Hamiltonian. If the perturbation is small, the Schrödinger equation does not have to be solved again from scratch. One starts from the solution to the unperturbed Schrödinger equation and computes small corrections to the orbitals and energies. Importantly, the 1st-order shift in energy does not require knowledge of the 1st-order shifts of the corresponding orbital. It can be computed from the unperturbed orbital.
Let the total Hamiltonian consist of an unperturbed operator, H 0 , and a small perturbation, H 1 : Let the solution to the unperturbed Hamiltonian, ψ (0) , be unique (no degeneracy). Perturbation theory predicts the 1st-order correction of the energy eigenvalue, E (1) , as The superscript (0) denotes the solution to the unperturbed Hamiltonian. "" denotes complex conjugation. If the wave function is normalized, the denominator is unity (and is not actually needed). It was included for the sake of generality. (The 0th-order solution will be non-normalized further down.) Again, the 1st-order correction to the energy can be computed from the unperturbed wave function.
In the acoustics problem, the Hamilton operator is replaced by the ω 2 -operator, ω 2 . Standing elastic waves in the unloaded crystal follow the relation: The indices i, j, k, and l run over x, y, and z. As before, superscript (0) denotes the solution to the unperturbed problem. Equation (93) is an eigensystem, similar to the timeindependent Schrödinger equation. The wave function is replaced by the displacement field, u (a vector, rather than a scalar). σ ( ) /dr j is the force density. σ ij is the stress tensor, given as: 14 (93) where ϵ kl is the deformation tensor and cijkl is the stiffness tensor. Because Equation (94) holds generally, the superscript (0) was omitted. The stiffness tensor of the unperturbed problem must be real (no internal friction) because perturbation theory requires the unperturbed operator to be self-adjoint. An unperturbed operator including dissipation would violate self-adjointness. Self-adjointness is not required, however, for the perturbation operator. Internal dissipation can therefore be included in the formalism as a perturbation. This perturbation is separate from the perturbation by the sample (Equation (95)). Because the unperturbed operator is real, the eigenvalue to the unperturbed problem, (ω ) ( ) , and the 0th-order displacement field, ( ) , are real, as well. A similar situation was encountered below Equation (21), when Z q , c q , and f ref were expressed as real parameters. When a sample exerts a periodic stress onto the surface, this stress amounts to a perturbation. The eigensystem with surface traction included is: The i in iω is the imaginary unit (not an index). δ(r  S) is the Dirac δ-function, where S is a position on the resonator surface. nj are the components of the surface normal. The eigenvalue, −ω , is complex because Z L,ijk is complex. Z L,ijk is the load impedance in tensor form, defined by: As in Equation (20), there is a minus sign in front of the stress because the stress is exerted into the direction of -z.
The displacement field in Equation (95), ≠ ( ) , is the solution to the total ω 2 -operator. One might write ≈ ( ) + ( ) with ( ) the 1st-order correction, but the 1st-order correction to the displacement field does not concern us any further. We are only interested in the shift of the eigenvalue. The eigenvalue is approximated as the sum of a 0thorder and a 1st-order term: Following perturbation theory, the 1st-order term can be computed from the 0thorder displacement field ( ) as: The denominator is needed because ( ) is not normalized. (It is not even dimensionless.) In the numerator, the volume integral has turned into a surface integral because of the δ-function. Complex conjugation as in Equation (92) is not needed because ( ) is real. In line 3, the center-dot is the vector product (contraction over one index, "dot product"). ":" denotes contraction over two indices and "⊗" is the outer product. This concludes the perturbation calculation. We are left with the task to compute Δf from ω 2 (1) and (ω 2 ) (0) . This calculation proceeds as: Taylor expansion ((1 + ε) 1/2 » 1 + ε/2 for ε ≪ 1) was applied. Combining Equation (99) with Equation (98) and using ω » ω ( ) leads to: This is the small-load approximation in tensor form. It was first written down in slightly different form by Pechhold in [194].
We briefly convince ourselves that Equation (100) reduces to Equation (25) for the parallel plate. For the parallel plate, the displacement occurs along x and is of the form u S cos k q z . The integral in the denominator turns into u S 2 A eff d q /2. The only nonzero component of Z L,ijk is Z L,xzx , called Z L in Equation (25). The surface normal is along z. The displacement in the numerator is u . For the parallel plate, Equation (100) turns into: With ρ q d q = m q = Z q 2f 0 this relation is equivalent to Equation (25): Limitations of the formalism are:  Piezoelectric stiffening is not included. That can be done (in tensor form). It is simply a matter of not letting oneself be intimidated by large equation systems.  Some perturbations may actually be large perturbations. Among these are the compressional waves, because the plate's stiffness under bending (not shear) may be too small to let the normal pressure exerted by compressional waves be a small perturbation [127].  The above mathematics covers the 1st-order perturbation, only. 3rd-order perturbation is sometimes needed (Box 2).  Calculating the vibration pattern of the unloaded crystal with electrodes is a challenge. If such a calculation is not feasible, the mode of vibration can still be determined experimentally with laser Doppler vibrometry (LDV, Figure 39C).
Among the benefits linked to Equation (100) are:  Equation (100) is general. It also applies to other resonators (such as torsional resonators [60,61,65] or nanoresonators [195]).  It clarifies, why the statistical weight in area averaging in Equation (26) is the square of the local amplitude of oscillation.  It explains why the Sauerbrey relation is slightly incorrect on the low overtones, even for rigid films in dry environments. The problem is linked to the modal mass (Section 8.6.).  It allows to quantitatively estimate the effects of increased energy trapping discussed around Figure 38. In particular, it explains why the freezing of a liquid drop on the resonator surface (water, hexadecane, …) lets the resonance frequency increase [196].
The limitations acknowledged: Equation (100) is the starting point for full-fledged 3D-modelling.

The 4-Element Circuit and the Electromechanical Analogy
An equivalent circuit (also "lumped-element circuit") contains discrete elements, networked together in some way. This involves two separate approximations. Firstly, the real-world device may or may not consist of discrete elements. Printed circuit boards (PCBs) often contain discrete elements, linked by conductive tracks with small resistance. Even then, stray capacitances (often needed for a faithful representation with an equivalent circuit) have no corresponding elements on the PCB. In real-word acoustics and mechanics, discrete elements are rare. Two masses linked together with a spring are an example ( Figure 1B). The bell is a counter example ( Figure 1C). A second approximation concerns discrete elements in the real world, which are represented with two impedances in the equivalent circuit. A real-world inductor has some ohmic resistance and is therefore depicted as an inductor in series with an ohmic resistor. In mechanics, a contact often dissipates energy and is therefore depicted as a spring in parallel to a dashpot ( Figure 31). Figure 42A shows an electrical equivalent circuit of a quartz resonator (or some other piezoelectrically driven resonator). When drawn without the load (in red on the righthand side), this is the Butterworth-van-Dyke circuit (BvD circuit, also: "4-element circuit"). The upper branch (the motional branch) contains an inductance, a capacitance, and a resistance. This circuit exploits the electromechanical analogy, which maps the mass, the spring, and the dashpot onto the inductance, the capacitance, and the resistance. The motional branch can be modeled with electrical impedances because piezoelectricity acts as an impedance converter. While the discrete impedances certainly are an idealization, Figure 42A reproduces the resonator's overall impedance close to the resonances well. Note that the values of L , R , and C 1 differ between overtones. There is a more general circuit (the Mason circuit [92]), which covers the entire frequency range (on-resonance and off-resonance, all overtones, parallel-plate model). In the Mason circuit, the impedances of the elements are not just inductances, resistances, or capacitances. For instance, one of them is iA eff Z q tan k q d q /2 . The Mason circuit also covers piezoelectric stiffening. (The 4-element circuit from Figure 42A does not.) Models of the QCM loaded with planar films can entirely be based on the Mason circuit [54,197]. The Mason circuit is interesting in a few ways, but working one's way from the diagram to the frequency shift is laborious [198].
The capacitance at the bottom in Figure 42 is a genuinely electrical element. It is the "parallel capacitance", formed by the electrodes on both sides of the plate. C 0 is larger than C 1 by about a factor of 1000 (Equation (88)). When the (electrical) Kirchhoff rules are applied to Figure 42A, an electrical admittance as shown in Figure 3 results.
It is instructive to also draw the corresponding mechanical circuit. A mechanical circuit of that kind was already shown in the introduction (Figure 1). The mechanical circuit contains a spring, a mass, and a dashpot. These mechanical impedances (such as iωm R for the mass and κ R /(iω) for the spring) have dimensions of force to velocity. The acoustic impedances, on the other hand, are ratios of stress to velocity. Within the parallel-plate model, one converts between the two by multiplication with the area, A eff . (The 4-element circuit is not necessarily outside the parallel-plate model. It is general.) Figure 42B makes the mechanical nature of the elements in the upper branch explicit. Importantly, the elements must be arranged in parallel (see Box 3). Following the mechanical Kirchhoff rules, impedances are additive, when the elements are arranged in parallel, while inverse impedances are additive, when the elements are arranged in series.
A side remark: Not everyone draws mechanical elements as in Figure 42B. If drawn as in Figure 43 on the right-hand side, the electrical Kirchhoff rules apply. The difficulty here is to always remember that the elements in this circuit are arranged differently from the elements in the real world. That would not problem in the context of Figure 42B because the discrete elements do not exist in the real world. It would be a problem in Figure  31. Figure 43. Electrical engineers sometimes draw mechanical equivalent circuits, such that the electrical Kirchhoff rules apply [52]. This is not the convention followed in Figure 42.
The values of the circuit elements in Figure 42B are: Within the parallel-plate model, δ q is the plate's loss angle, G q ''/G q '. For practical resonators, the energy dissipated in the electrodes also contributes to ξ R . Figure 42B represents piezoelectric coupling as a transformer. A key parameter of a transformer is the turns ratio, ϕ = n 2 /n 1 . n 1 and n 2 are the numbers of turns. Voltage, U, and current, I, are transformed as: Following Equation (103), an impedance on one side of the transformer takes a different value when seen from the other side. The converted impedance is: The transformer in Figure 42B is an impedance converter in this sense. It converts between an electrical and a mechanical impedance (Z mech = ϕ 2 Z el , with ϕ = e 26 A eff /d q ).
The following equations derive the small-load approximation from Figure 42B. In the absence of a load, the resonance condition is Neglecting the dashpot, the resonance frequency is In the presence of the load, the resonance condition is Before solving this equation for ω, we turn it into a more familiar form, which is: Again neglecting the dashpot, the shifted resonance frequency is: In the second step, the meaningful solution (out of the two solutions of the quadratic equation in ω) was selected. The quadratic term under the root was neglected because the load impedance is small. With m q = ρ q d q and Z q = ρ q G q / , Equation (109) turns into: This is the small-load approximation.

Amplitude of Oscillation, Effective Area
The effective turns ratio, ϕ, of the transformer in Figure 42B leads to a relation between the current into the electrodes and the velocity at the resonator surface. On resonance, the velocity, v S , and the current, , are related as [199]: A factor of 1/2 enters because the current is proportional to the difference in velocity between the front and the back (which is 2v S ). The relation between the current into the electrodes, I , and the nominal external voltage, U , may be nontrivial. Further complicating the situation, vector network analyzers often control the power (in units of dBm), rather than voltage or current. With 5 MHz resonators on the fundamental, an area of A eff = 10 mm 2 , a voltage of 100 mV, and a motional resistance of R = 500 Ω (typical for experiments in liquids), an amplitude of oscillation of 0.1 nm results. The shear angle is u S /δ, which is below 10 3 . A more detailed calculation of the amplitude would have to account for energy trapping. The current through the motional branch may also be affected by the analyzer's output resistance and by the current through the parallel capacitance, C 0 . For that reason, "»" was written on the right-hand side in Equation (111).
If the resonator is under voltage control, one may also remember the relation: Δω d 26 = 3.110 12 m/V is the piezoelectric strain coefficient of AT-cut quartz. Equation (112) was confirmed experimentally in [200]. The transformer's effective turns ratio, ϕ = e 26 A eff /d q , also leads to an equation for the plate's effective area, which is: For the proof see Section 7.4 in [5]. The same approximations as in Equation (111) apply.

Modal Mass, Sauerbrey Equation for Plates with Energy Trapping
Flexural admixtures not only affect measurements in liquids (because of compressional waves), but also measurements of the mass of a film in the dry. Start from the classical formula for the resonance frequency of a harmonic oscillator, ω » (κ R /m R ) / . This formula holds for discrete objects. Applied to elastic bodies vibrating in some vibration mode, equivalent parameters must be defined. The definition of κ R is unessential here. Write the total mass as m R + A eff m f with m R the "modal mass" and m f the film's mass per unit area. Assuming A eff m f ≪ m R and applying the Taylor expansions (1 + ε) 1/2 » 1 + ε/2 and (1 + ε) -1 » 1 -ε leads to: At first glance, Equation (114) appears to differ by a factor of 1/2 from the Sauerbrey result, which is: m q is the mass per unit area of the resonator plate. The comparison shows that m R is only half the mass of the resonator plate, A eff m q . This is so because a plate in thickness-shear motion contains nodal planes. The kinetic energy is 1/2A eff m q 〈v 〉 t,V . Averaging occurs over time, t, and volume, V. Volume averaging produces a factor of 1/2 because 〈cos ( )〉 = 1/2, which implies m R = A eff mq/2.
Importantly, this factor of 1/2 turns into some other numerical factor when flexural modes are present because the flexural motion does not vanish at the nodal planes of the shear motion. The modal mass increases and -Δf/n decreases in consequence. This effect slightly modifies the Sauerbrey equation. Figure 44 shows an experimental example. The sample was a polymer film in air. Because this film is viscoelastic, Δf/n decreases in proportion to at large n (see Equation (46)). From the slope, one can infer ′. If interpreted with Equation (45), the positive slope would indicate negative shear compliance, but Equation (46) solves that problem.
Regardless of the slope: The low overtones must be excluded from this analysis, because the respective vibration modes have an increased modal mass.

Figure 44.
Overtone-normalized shifts of frequency and bandwidth obtained after depositing a polymer film with a thickness of about 500 nm onto the resonator plate. The area in green denotes values of -Δf/n which deviate from the prediction by Equation (46). The deviations go back to flexural motion. Adapted from [80].
The above argument started from ω 0 » (κ R /m R ) / . The argument can also start from the tensor form of the small-load approximation. The denominator in Equation (100) is proportional to the modal mass.

Combined Instruments
Because the QCM is simple, it can be easily combined with other instruments. For instance, it is quite trivial to push a QCM sensor against the prism of an ATR-IR spectrometer and obtain an IR spectrum from the sample on the resonator with monolayer sensitivity (ex-situ, evidently). One can rather easily do AFM, XPS, and Raman spectroscopy on a QCM, either in-situ or ex-situ. The examples are numerous. Two problems are worth a mention:


When combining a QCM with an AFM [201], the QCM does not usually respond to the contact with the AFM tip because the contact is too small. That can be understood from Equation (70) together from the Mindlin result for the transverse contact stiffness, which is κ P = 2G*r C with G* an effective modulus and r C the contact radius [117] (see Section 5). Inserting values (G* » 10 GPa, r C » 5 nm) leads to a frequency shift below 0.1 Hz (at 5 MHz). An AFM tip tapping onto the resonator amounts to a nanoscopic object perturbing the motion of a macroscopic object. It does so, in principle, but the effects usually disappear in the noise. In the reverse direction (the QCM perturbing the motion of the tip), there is a strong influence [202]. However, experiments of this kind can also be done with any other actuator. Nanoscale dynamicalmechanical studies based on an AFM tip in contact with a vibrating substrate are commonplace ( [203] and others).  In-situ combination with dielectric spectroscopy [204][205][206][207] or electrical cell-substrate impedance spectroscopy (ECIS [208]) is possible. A difficulty arises, when the sample requires an oxidic substrate, such as SiOx, because the commercially available SiOx coatings may be too thick. The electric field then does not reach to the sample. More technically, the coating's capacitance, CSiOx, is so small, that its impedance dominates the sample's overall electrical impedance. The properties of the sample are then masked by the term 1/(iωCSiOx). Thin dielectric coatings are needed. The text proceeds with two particularly well studied combinations, which are the electrochemical QCM and the combination with optical reflectometry.

The Electrochemical QCM (EQCM)
Among the combinations with other techniques, the combination with electrochemistry is most advanced. For early works see [1619,209]. More recent reviews are provided in [210][211][212].
Electrogravimetry can be rather complicated and we do not elaborate on the details. Matters are transparent as long as one sticks to electrodeposition [213] or electroetching [214] with layer thicknesses of many nanometers. The QCM then operates in the gravimetric regime. The charge passed through the electrode surface may be converted to a mass, using Faraday's law [103,215]. If this mass equals the mass as inferred from the Sauerbrey equation, the "current efficiency" is unity. If this current efficiency is less than unity, one suspects side reactions or soft layers. For soft layers, the bandwidth and the overtone dependence of Δf/n can be analyzed to derive the layer's viscoelastic constants (Equation (52), [104]). If the current efficiency is larger than unity, roughness may be a reason [216].
The difficulties admitted, there is additional information available for interpretation. Firstly, there is the electrical current as a function of time (or the voltage, in case the experiment occurs under current control). If the experiment occurs outside the gravimetric range, charge and frequency shift are not simply related by Faraday's law, but they should still be related in some other way. Also, kinetics comes to help. The fast QCM discussed in Section 3.5. resolves the kinetics with a time resolution down to 100 µs. The response times determined this way give clues to the nature of the processes involved. Finally, the frequency resolution can be improved by running the experiment repetitively and accumulating data. If a certain process occurs in response to a change in electrode potential on the time scale of, for instance, one second, accumulation overnight improves the resolution in frequency by a factor of about (40000) 1/2 (by N 1/2 with N the number repeats). The frequency resolution readily drops to below 1 mHz. A caveat: For the conventional QCM, dirt, scratches, or even small bubbles are not necessarily detrimental because the QCM averages over the entire active area (~10 mm 2 ). A minute response (of a few mHz) to a jump in electrode potential, however, may entirely result from such local heterogeneities.
The Paris group has exploited the combination of electrochemical impedance spectroscopy (EIS) with the QCM in numerous papers [103 and others]. They call the technique "AC electrogravimetry". The electrode potential is modulated sinusoidally. The QCM is driven by an oscillator circuit, the output of which is fed into a frequency-tovoltage converter. The respective voltage oscillates with the modulation frequency. It is displayed together with the electrical current, often in the same diagram, often showing similar features.
The molecular details of electrochemistry can be enormously complex. The EQCM (similar to EIS) only measures integral quantities (such as mass transfer, area-averaged softness, and area-averaged roughness) and their kinetics. It can still aid the structural investigations (examples in [221,222]) and provide for constraints to the models.

Combination with Optical Reflectometry
With regard to label-free biosensing, SPR spectroscopy [15] outperforms the QCM insofar, as SPR has a lower limit of detection and smaller baseline drifts. Irrespective of the competition, there are interesting conceptual similarities between the QCM and optical reflectometry. Also, there are benefits in the combination of the two [93,151,223].
Similar to the QCM, SPR spectroscopy exploits a shift of a resonance condition, where the frequency is replaced by the kx-vector of the surface plasmon. The analogy is illustrated in Figure 45.
If the adsorbed layer is thinner than the decay depth of the plasmon, δ opt , (which happens to be similar to the decay depth of the shear wave of the QCM), the shift in k x is proportional to the 1st invariant of the refractive index profile, J 1 , given as [224]: Subscripts sub and bulk denote the substrate and the bulk liquid, respectively. n r is the refractive index. n r 2 is equal to the dielectric constant at optical frequencies, ε. While ε and n r 2 are complex, in principle, they are treated as real parameters, here. In SPR spectroscopy, the substrate is a metal (often gold), which has a large, negative ε. The term n sub 2  n r 2 (z) may be approximated as n sub 2  n bulk 2 (equal to ε sub  ε bulk ) and be pulled out of the integral. Neglecting prefactors, the shift in k x induced by the adsorbate obeys [225]: Figure 45. Both the QCM and SPR spectroscopy exploit the shift of a resonance condition. The peak in the conductance curve is the analog of a dip in the reflectance curve. The resonance frequency is replaced by the coupling angle (that is, the angle of incidence with minimum reflectance). The coupling angle is related to the plasmon's kx-vector. The shift in the coupling angle is proportional to the adsorbed amount.
Both Δf + iΔΓ (in QCM experiments, Equation (58)) and Δkx (in SPR spectroscopy, Equation (117)) depend on an integral of a response function. Algebraically, the response functions are similar. Differing from SPR spectroscopy, the response function to be used for the QCM saturates to unity. This happens whenever the viscosity in the layer is much larger than the viscosity of the bulk. The contrast function in optics, on the contrary, is mostly proportional to concentration, because the refractive index increment, dn r /dc, is small. The SPR response is roughly proportional to the adsorbed amount, while the QCM probes the distance to the surface, at which the contrast function finally falls out of saturation (d ac in Figure 46). When the layer expands, the SPR response changes by a small amount because the increase in thickness is compensated by a decrease in concentration. The QCM mostly notices an increase in acoustic thickness. These remarks apply to all variants of optical reflectometry. Because SPR and QCM see the adsorbed amount and the thickness, respectively, the degree of swelling may be inferred from the combination of SPR and QCM (more generally, from the combination of optical reflectometry with shearwave acoustic reflectometry). SPR can be combined with the QCM in-situ, but that requires grating coupling on the surface of a resonator crystal, which is expensive [226,227]. Ellipsometry is an option [151 and others], but ellipsometry can be challenging in the details, especially when the beam passes through windows. Given that the electrodes often consist of gold, variants of optical reflectometry can be employed, which are simpler than full-fledged ellipsometry [228,229]. If it shall be SPR on the side of optical reflectometry, it can be Love wave sensors on the side of shear-wave acoustics [230]. Love wave sensors are a variant of the surface acoustic wave devices [231,232]. For reasons in the details, these are more easily combined with SPR than the QCM. The remarks above on saturation of the contrast function apply to Love waves and to shear waves as excited by the QCM in essentially same way [233].

Appendix D. Justification of a Complex k q in the Lu-Lewis Equation
The original form of the Lu-Lewis equation was based on a real-valued wave vector in the resonator, k q . In order to account for lossy samples, k q was turned in to an effective complex k q in Equation (20). The imaginary part is not linked to the resonator's intrinsic losses. If the sample is lossy, the displacement at z = dq is of the form: A phase contained in the complex resonance frequency ω res was lumped into k q . The phase is caused by the lossy sample. Again, k q is an effective wave number, introduced to keep the algebra simple. Because the apparent k q '' is small, the errors can be tolerated. (One is not forced to do that. The calculation also runs along with complex ω res .) u (d q ,t) = u S cos k q 'd q exp(iω res t) = u S cos k q 'd q 1 + iϕ exp(iω res 't) = u S cos k q 'd q + i cosh k q ''d q exp(iω res 't) = u S cos k q d q exp(iω res 't) (A22)