Fast and Accurate Finite Transducer Analysis Method for Wireless Passive Impedance-Loaded SAW Sensors

Abstract

An accurate and fast simulation tool plays an important role in the design of wireless passive impedance-loaded surface acoustic wave (SAW) sensors which have received much attention recently. This paper presents a finite transducer analysis method for wireless passive impedance-loaded SAW sensors. The finite transducer analysis method uses a numerically combined finite element method-boundary element method (FEM/BEM) model to analyze non-periodic transducers. In non-periodic transducers, FEM/BEM was the most accurate analysis method until now, however this method consumes central processing unit (CPU) time. This paper presents a faster algorithm to calculate the bulk wave part of the equation coefficient which usually requires a long time. A complete non-periodic FEM/BEM model of the impedance sensors was constructed. Modifications were made to the final equations in the FEM/BEM model to adjust for the impedance variation of the sensors. Compared with the conventional method, the proposed method reduces the computation time efficiently while maintaining the same high degree of accuracy. Simulations and their comparisons with experimental results for test devices are shown to prove the effectiveness of the analysis method.

1. Introduction

Wireless passive surface acoustic wave sensors have been investigated since the 1990s because they have several advantages, such as, no need for power supply, good environmental durability, high sensitivity, etc. [1]. There are various types of surface acoustic wave (SAW) sensors, pressure sensors [2], gas sensors [3], torque sensors [4], etc. Compared with them, impedance-loaded surface acoustic wave sensors offer a much wider application range [5,6]. These sensors use the SAW chip as a transponder and the load sensing part can be any impedance-varying sensor that matches the scope of impedance variation. Impedance-loaded sensors with a wide application range have received much attention in recent years [7,8,9].

Conventional impedance-loaded sensors have two types of structures: The delay line structure and resonator structure [6,7,10,11]. The basic construction of delay line structure sensors is shown in Figure 1. The SAW transponder is composed of several interdigital transducers (IDTs), one of which is used as a transmitting/receiving IDT responsible for transmitting and receiving wireless signals and the others are used as reflectors. One of the reflectors is loaded with conventional sensors whose varying impedance can be of any type, resistance, capacitance, or inductance. The wireless interrogating electromagnetic wave signal is received by the antenna and transformed to a surface acoustic wave by the transmitting/receiving IDT. The surface acoustic wave propagates along the substrate surface and is reflected by the reflectors. At the same time, the reflection coefficient would be influenced by the externally loaded sensors impedance. The returned surface acoustic wave is transformed back into an electromagnetic (EM) wave by the transmitting/receiving IDT and sent by the antenna. By changing the impedance of outside conventional sensors (electrical boundary condition), the scattering parameters S₁₁ of the sensors also changed, which can be detected by the interrogators.

Impedance-loaded SAW sensors operate at high frequency and should be designed to maintain low loss in order to achieve wireless measurement. Therefore, the sensor design must be based on an accurate model. The numerical combined FEM/BEM model is now widely used to derive the influence of device parameters on the frequency response in the design of SAW devices [12,13]. It is more accurate than other conventional models like the equivalent circuit model (ECM) model, this is because it includes the second order effects which other models do not consider. The FEM/BEM model is comprised of a periodic model and a non-periodic model.

The periodic FEM/BEM model can be used to calculate the characteristics of periodic structures, from which more accurate and flexible coupling-of-modes (COM) parameters can be derived. In our previous work [14,15], the COM model of impedance-loaded SAW sensors was constructed based on the models of SAW IDT. The COM parameters, whose accuracy determines the accuracy of the COM model, were calculated using a periodic FEM/BEM model. Next, the SAW device admittance can be derived by the COM model. Moreover, the effects of the impedance matching circuit and the loaded conventional sensors are included in the model.

The non-periodic FEM/BEM model is derived from the basic laws of physics such as Newton’s second law, piezoelectric constitutive equations and Maxwell equations to form a set of system equations. The admittance of a SAW device can be directly derived from the solutions of the system equations with its geometrical characteristics, materials properties, and mechanical and electrical boundary conditions known. The non-periodic FEM/BEM model is more accurate than the COM model for radio frequency (RF) devices because it includes most of the second order effects, but it requires more computing time. In this model, calculating the coefficients of the system equations, which involves double numerical integrals and solving the equations, takes most of the central processing unit (CPU) time. The main calculation content can be divided into four parts, among which the bulk wave contribution takes a long computation time, for it requires the calculation of the general Green’s function at all frequencies. A huge amount of computing time will be needed if the adaptive quadrature method was adopted. Weibiao Wang et al. presented a method based on the sinusoidal weighted integral [16]. However, the case when the integral variable tends to zero was not mentioned. In our paper, this case was carefully considered and the whole numerical properties in the local area were treated accordingly. Our new algorithm leads to a half reduction in required computation time.

When the non-periodic FEM/BEM model is applied to impedance-loaded SAW sensors, the calculation is more time consuming. The loaded impedance is varying, which changes the electrical boundary condition of the IDT. This means that every time the loaded impedance is changed the model calculation needs to be performed. When the electrical boundary condition changes, the complete system equation should be solved once, which requires a long time for solving the large dimension linear equations. For the sensor performance calculation, the loaded impedance might be assigned values twenty times, which leads to the system equation requiring to be solved twenty times. Thus, the calculated cost will be considerably increased. In this paper, a finite transducer analysis method using a non-periodic FEM/BEM model was introduced to simulate impedance-loaded SAW sensors. A fast algorithm was introduced to calculate the coefficients of the system equations. It was shown that, after some modification to the non-periodic model, the sensor admittance can be obtained at any external loaded impedance by using the network theory after performing the calculation for only two states: Short circuit and open circuit. Section 2 provides a short introduction to the non-periodic FEM/BEM model. Section 3 describes in detail an accurate and fast algorithm to calculate the bulk wave part of the equation coefficient. Section 4 and Section 5 present the method that is used to calculate the sensor response at any outside sensor impedance. Section 6 illustrates the simulated results of a test device with different structures. Comparison of the simulated results to the measured results showed good accordance and Section 7 presents the conclusion that our finite transducer analysis of the impedance-loaded SAW sensors is valid.

2. Basics of the Non-Periodic FEM/BEM Model

As shown in Figure 2, the structure is composed of a semi-infinite piezoelectric substrate extending in the half space z < 0, whose surface is loaded with metal strip gratings. The grating electrodes are parallel to the y axis and are assumed to be long enough along the y axis such that any dependence on y can be ignored. Throughout this paper, the time harmonic factor exp (jωt) will be omitted.

It is supposed that there are N_e electrodes in the array. The first IDT has N₁ electrodes and the second IDT has N₂ electrodes. Obviously, N_e equals the sum of N₁ and N₂. The location of the kth electrode is X_k, the width of the kth electrode is 2a_k, and the applied potential on the kth electrode is V_k. Because stress and charge distributions are zero outside the electrode-substrate interface, the three components (mechanical displacement, u; stress vector, t_s; electrical potential, $\varphi$) and the surface charge density $\mathsf{\sigma }$ can be related using the general Green’s function G(x) [17]:

$$\left(\begin{array}{c}\mathrm{u}(x)\\ \varphi (x)\end{array}\right)={\displaystyle {\int }_{-\infty }^{+\infty }\mathrm{G}(x-x\prime )}(\begin{array}{c}{\mathrm{t}}_s(x\prime )\\ \mathsf{\sigma }(x\prime )\end{array})dx\prime ,$$

(1)

The distributions of stress and free charge are expanded using the first-kind Chebyshev polynomial T_n(x) for the kth electrode:

$$\left(\begin{array}{c}{\mathrm{t}}_s(x)\\ \mathsf{\sigma }(x)\end{array}\right)={\displaystyle \sum _{n=0}^{N_k^{ch}}\left(\begin{array}{c}{B}_{kn}^{\mathrm{t}}\\ B_{kn}^{\mathsf{\sigma }}\end{array}\right)\frac{T_n(x_k)}{\sqrt{1-x_k^2}}}, x_k=\frac{x-X_k}{a_k}, \begin{array}{l}k=1,\dots ,N_e\\ n=0,\dots ,N_k^{ch}\end{array},$$

(2)

where $N_k^{ch}$ is the number of Chebyshev polynomials used to expand stress and charge on the kth electrode. By incorporating Equation (2) in Equation (1) and timing the Chebyshev polynomials $\frac{T_i(x_j)}{\sqrt{1-x_j^2}}$ at both sides of Equation (1), we obtain

$${\left(\begin{array}{c}{c}_j^u\\ c_j^{\varphi }\end{array}\right)}_i={\displaystyle \sum _{k=1}^{N_e}{\displaystyle \sum _{n=0}^{N_k^{ch}}{Y}_{kn}^{ij}\left(\begin{array}{c}{B}_{kn}^{\mathrm{t}}\\ B_{kn}^{\mathsf{\sigma }}\end{array}\right)}}\begin{array}{c}j,k=1,\dots ,N_e\\ i,n=0,\dots ,N_k^{ch}\end{array},$$

(3)

where

$$Y_{kn}^{ij}=a_k{\displaystyle \underset{-1}{\overset{1}{\int }}{\displaystyle \underset{-1}{\overset{1}{\int }}G(a_j{x}_j-a_k{x}_k+{\Delta }_{jk})\frac{T_i(x_j)}{\sqrt{1-x_j^2}}\frac{T_n(x_k)}{\sqrt{1-{x_k}^2}}}d{x}_{j}d{x}_k}, {\Delta }_{jk}=X_j-X_k,$$

(4)

$${\left(\begin{array}{c}{c}_j^u\\ c_j^{\varphi }\end{array}\right)}_i={\displaystyle {\int }_{-1}^1\left(\begin{array}{c}\mathrm{u}(a_j{x}_j+X_j)\\ \varphi (a_j{x}_j+X_j)\end{array}\right)\frac{T_i(x_j)}{\sqrt{1-{x_j}^2}}d}{x}_j,$$

(5)

The key point in the numerical calculation of the non-periodic FEM/BEM model is calculating the equation coefficient $Y_{kn}^{ij}$ in Equation (3). Because G(x), which is Green’s function in the space domain, is numerically calculated, its double integral should also be numerically calculated. It is necessary to decompose the function into a sum of various contributions (electrostatic contribution, $G^0(x)$; surface wave contribution, $G^s(x)$; asymptotic contribution, $G^{\infty }(x)$; bulk wave contribution, $G^b(x)$) due to its extremely fast variations [16].

$$G(x)=G^0(x)+G^s(x)+G^{\infty }(x)+G^b(x),$$

(6)

Incorporating Equation (6) in Equation (4) leads to

$$Y_{kn}^{ij}={\left(Y_{kn}^{ij}\right)}^0+{\left(Y_{kn}^{ij}\right)}^s+{\left(Y_{kn}^{ij}\right)}^{\infty }+{\left(Y_{kn}^{ij}\right)}^b,$$

(7)

where ${\left(Y_{kn}^{ij}\right)}^0$, ${\left(Y_{kn}^{ij}\right)}^s$, ${\left(Y_{kn}^{ij}\right)}^{\infty }$ and ${\left(Y_{kn}^{ij}\right)}^b$ are the equation coefficients related to electrostatic contribution, surface wave contribution, asymptotic contribution and bulk wave contribution, respectively. All the above parts can be easily integrated by numerical methods [12,13,16] except the bulk wave contribution. The bulk wave contribution relates to the computation of the general Green’s function at all frequencies, of which the integrated expression requires a long computing time. In the next section, an accurate and fast algorithm for computing the bulk wave contribution part of $Y_{kn}^{ij}$ is presented.

3. Accurate and Fast Algorithm for Computing the Bulk Wave Contribution Part of ${\mathit{Y}}_{\mathit{k}\mathit{n}}^{\mathit{i}\mathit{j}}$

The bulk wave contribution ${\left(Y_{kn}^{ij}\right)}^b$ can be expressed as [16]:

$${\left(Y_{kn}^{ij}\right)}^b=j^{n-k}\pi a_j{\displaystyle {\int }_{-s_m}^{s_m}\frac{{\overline{H}}^b(s)}{\left|s\right|}\exp (-js\omega \Delta )J_k(s\omega a_i)J_n(s\omega a_j)}ds,$$

(8)

where $\Delta =X_i-X_j$ is the distance between the ith and the jth electrodes. s is the slowness, and $\omega$ is the angular frequency. In Equation (8), the calculation of $\frac{{\overline{H}}^b(s)}{\left|s\right|}$, which is the bulk wave part of the Green’s function, requires most of the computing time. The property that the Green’s function is only dependent on the slowness for homogeneous substrates makes it possible that the computation of $\frac{{\overline{H}}^b(s)}{\left|s\right|}$ can be done only once. Then the Green’s function can be obtained for all cases on the same substrate simply by interpolation using Chebyshev expansion.

The values of $\frac{{\overline{H}}^b(s)}{\left|s\right|}$ in the slowness domain vary within wide limits. For example, Figure 3 shows the variation of $\frac{{\overline{H}}_{44}^b(s)}{\left|s\right|}$ in the scope of $s\in (0,1\times {10}^{-3})$ for 42° LiTaO₃ on which the leaky wave exists. There are still some sharp changes and turns near the leaky wave pole point and branch points after removing the surface wave part of the Green’s function. The function must interpolate in several intervals in the slowness domain to ensure accuracy and lower interpolation expansion numbers. The larger the amplitude of function variation, the greater the interval number. As shown in Figure 4, $\frac{{\overline{H}}_{44}^b(s)}{\left|s\right|}$ is monotonic and asymptotic when $s>1\times 10^{-3}$. So, fewer sections can be divided for interpolation in that area. At each interval, $\frac{{\overline{H}}^b(s)}{\left|s\right|}$ is monotonic and smooth. The interval number is 2N + 1, denoted as −s_m = s − N < s − N + 1 <, …, <s₀ = 0 < s₁, …, <s_N = s_m.

The Chebyshev expansion of $\frac{{\overline{H}}^b(s)}{\left|s\right|}$ is performed in the Qth interval [s_Q, s_Q+1],

$$\frac{{\overline{H}}^b(s)}{\left|s\right|}={\displaystyle \sum _{p=0}^{N_Q}{\overline{c}}_p^Q}{(\frac{s-s_{oQ}}{{\Delta }_Q})}^p,$$

(9)

${\overline{c}}_p^Q$ are the expansion coefficients of Chebyshev polynomials and they only need to be calculated once. The expansion numbers $N_Q$ can be 8–10 for large intervals and 3–5 for small intervals.

Substituting Equation (9) into Equation (8), we can get

$${\left(Y_{kn}^{ij}\right)}^b={\displaystyle \sum _{p=0}^{N_Q}{\overline{c}}_p^Q{\displaystyle {\int }_{s_Q}^{s_{Q+1}}{(\frac{s-s_{oQ}}{{\Delta }_Q})}^p{J}_k(s\omega a_i)J_n(s\omega a_j)\exp (-js\omega \Delta )}ds},$$

(10)

Equation (10) can be decomposed into typical sine or cosine weighted integration forms. Then, the adaptive quadrature method could be applied to approximate the integration.

In the general case when the electrode number is greater than 100, $s\omega \Delta$ could be larger than 200 for long transducer RF SAW devices. The integral of Equation (10) is an oscillatory integral, where the computation cost is proportional to the value of $s\omega \Delta$. When the scale of $s\omega \Delta$ is larger than 1 × 10⁴, the accuracy of the adaptive quadrature is poor.

A robust and accurate algorithm is needed to calculate the appropriate cases for all values of $s\omega \Delta$. The best way to solve this problem is by expanding the integrand to an analytic series. Thus, another algorithm is presented according to the properties of the Bessel function. The products of the Bessel function are expanded to a series, which makes the integral analytic in nature and thereby integrable.

When $\omega as\le 8$, the limited series expression can be used as follows:

$$J_k(as)J_n(bs)=\frac{1}{\mathsf{\Gamma }(n+1)}{(\frac{as}{2})}^k{(\frac{bs}{2})}^n\times {\displaystyle \sum _{q=0}^{\infty }\frac{{{(-1)}^q}_2{F}_1(-q,-k-q;n+1;\raisebox{1ex}{$b^2$}\!\left/ \!\raisebox{-1ex}{$a^2$}\right.)}{q!\mathsf{\Gamma }(k+q+1)}}{(\frac{as}{2})}^{2q},$$

(11)

where $b\le a$, $\mathsf{\Gamma }(n)$ is the gamma function and ${}_2{F}_1(\alpha ,\beta ,\gamma ,x)$ is the hypergeometric function. When $\omega as>8$, the asymptotic expansion formula of the Bessel function can be used as follows:

$$J_k^{Asy}(s)\approx \sqrt{\frac{2}{\pi x}}[\cos (s-\frac{k\pi }{2}-\frac{\pi }{4}){\displaystyle \sum _{p=0}^{\infty }\frac{{(-1)}^p(k,2p)}{{\left(2s\right)}^{2p}}-\sin (s-\frac{k\pi }{2}-\frac{\pi }{4})}{\displaystyle \sum _{p=0}^{\infty }\frac{{(-1)}^p(k,2p+1)}{{\left(2s\right)}^{2p+1}}}],$$

(12)

According to the numerical analysis, it is accurate enough for p = 5 in Equation (12) when $k,n\le 6$.

Substituting Equations (11) and (12) into Equation (3), we can get the analytical and integrable forms of Equation (10) as follows:

$${\left(Y_{kn}^{ij}\right)}^b={\displaystyle \sum _{p=0}^{N_Q}{\overline{c}}_p^Q}{\displaystyle \sum _{l=0}^{N_{kn}}{J}_{kn}(l){\displaystyle {\int }_{s_Q}^{s_{Q+1}}{(\frac{s-s_{oQ}}{{\Delta }_Q})}^p}}{(\omega a_{i}s)}^k{(\omega a_{j}s)}^{n+2Q}\exp (-js\omega \Delta )ds,$$

(13)

or

$${\left(Y_{kn}^{ij}\right)}^b={\displaystyle \sum _{p=0}^{N_Q}{\overline{c}}_p^Q}{\displaystyle {\int }_{s_Q}^{s_{Q+1}}{(\frac{s-s_{oQ}}{{\Delta }_Q})}^p{J}_k^{Asy}(\omega a_{i}s)J_n^{Asy}(\omega a_{j}s)\exp (-js\omega \Delta )ds},$$

(14)

In Equation (13), the general form of the integral is ${\displaystyle {\int }_{s_Q}^{s_{Q+1}}{(\omega a_{i}s)}^t\exp (-js\omega \Delta )ds}$. In Equation (14), the general form of the integral is ${\displaystyle {\int }_{s_Q}^{s_{Q+1}}{(\omega a_{i}s)}^{-t}\exp (-js\omega \Delta )cos[(a_j\pm a_i)\omega s]ds}$ or ${\displaystyle {\int }_{s_Q}^{s_{Q+1}}{(\omega a_{i}s)}^{-t}\exp (-js\omega \Delta )\sin [(a_j\pm a_i)\omega s]ds}$. These can be calculated using integration by parts for the large value $\Delta$ case.

Thus, the accurate calculation algorithm of Equation (8) can be obtained using the series expansion of $J_k(s\omega a_i)J_n(s\omega a_j)$ and polynomial interpolation of $\frac{{\overline{H}}^b(s)}{\left|s\right|}$. Figure 5 shows the calculated results of $Y_{kn}^{ij}$ with different ${\Delta }_{ij}$ for different methods. The value of k as well as n is 0, as the first series represents the main item. The frequency is close to the center frequency of the transducer. The black line represents the results calculated using our approximation algorithm, and the red and blue lines are the results calculated from adaptive quadrature. The calculation time for our algorithm is one tenth that of the previous algorithm while maintaining the same degree of precision.

4. Non-Periodic FEM/BEM Model of the Impedance-Loaded SAW Sensors

After calculating $Y_{kn}^{ij}$, the FEM method was applied to derive the relation between $c_j^u$ and $B_{kn}^{\mathrm{t}}$ [12,13]:

$${\left(c_j^u\right)}_i=Y_e{B}_{kn}^{\mathrm{t}}.$$

(15)

Equation (3) can be rewritten in the matrix form as

$${\left(\begin{array}{c}{c}_j^u\\ c_j^{\varphi }\end{array}\right)}_i=\mathrm{Y}\left(\begin{array}{c}{B}_{kn}^{\mathrm{t}}\\ B_{kn}^{\mathsf{\sigma }}\end{array}\right)=\left(\begin{array}{cc}{Y}_u& Y_{u\varphi }\\ Y_{\varphi u}& Y_{\varphi }\end{array}\right)\left(\begin{array}{c}{B}_{kn}^{\mathrm{t}}\\ B_{kn}^{\mathsf{\sigma }}\end{array}\right),$$

(16)

where Y is the matrix form of $Y_{kn}^{ij}$, and $Y_u$, $Y_{u\varphi }$, $Y_{\varphi u}$, and $Y_{\varphi }$ are matrices related to the mechanical-mechanical, mechanical-electrical, electrical-mechanical and electrical-electrical couplings, respectively.

Substituting Equation (15) into Equation (16), we get

$${\left(\begin{array}{c}0\\ c_j^{\varphi }\end{array}\right)}_i=\left(\begin{array}{cc}{Y}_u-Y_e& Y_{u\varphi }\\ Y_{\varphi u}& Y_{\varphi }\end{array}\right)\left(\begin{array}{c}{B}_{kn}^{\mathrm{t}}\\ B_{kn}^{\mathsf{\sigma }}\end{array}\right),$$

(17)

Equation (17) is the system equation for the non-periodic FEM/BEM model. As shown in Figure 6, the electrical boundary conditions of the impedance-loaded SAW sensors will be applied to Equation (17).

The potential of the ground bus bar of the IDTs is noted as V₁ and the potential of the output and input IDTs are noted as V₂ and (1 + V₁), respectively. The omitted time harmonic coefficient is $\exp (j\omega t)$. The condition dictates that the total electric charge must be zero. Substituting the electrical conditions thus gives

$${\left(\begin{array}{c}0\\ c_j^{\varphi }\\ 0\\ 0\end{array}\right)}_i=\left(\begin{array}{c}{Y}_u-Y_e\\ Y_{\varphi u}\\ 0\\ 0\end{array}\begin{array}{c}{Y}_{u\varphi }\\ Y_{\varphi }\\ C_2\\ C_4\end{array}\begin{array}{c}0\\ C_1\\ C_3\\ C_5\end{array}\right)\left(\begin{array}{c}{B}_{kn}^{\mathrm{t}}\\ B_{kn}^{\mathsf{\sigma }}\\ {\mathrm{V}}_1\\ {\mathrm{V}}_2\end{array}\right),$$

(18)

where C₁ is a $2\times {\displaystyle \sum _{k=1}^{N_e}{N}_k^{ch}}$ dimension submatrix, C₂ and C₄ are ${\displaystyle \sum _{k=1}^{N_e}{N}_k^{ch}}$ dimensional row vectors, and C₃ and C₅ are two-dimensional row vectors.

For the elements related to the N₁ input IDT in $c_j^{\varphi }$, the expressions are $c_j^{\varphi }(N_1)=\pi$, $C_1(N_1,1)=1$, and $C_1(N_1,2)=0$. For the elements related to the N₂ IDT in $c_j^{\varphi }$, the expressions are $c_j^{\varphi }(N_2)=0$, $C_1(N_2,1)=\pi$, and $C_1(N_2,2)=-\pi$.

C₂ is determined by the charge conservation condition, C₃ = 0, and C₄ and C₅ are determined by ${\displaystyle \sum _{k=N_1+1}^{N_e}{a}_{k}j\omega \pi B_{k0}^{\sigma }}=-\frac{V_2-V_1}{Z_L}$, where Z_L is the impedance of the external loaded sensor.

If Z_L is known, Equation (18) can be solved and Equation (19) is obtained as the final solution.

$$\{\begin{array}{c}{U}_1=1\\ U_2=V_2-V_1\\ I_1={\displaystyle \sum _{k=1}^{N_1}{a}_{k}j\omega \pi B_{k,0}^{\mathsf{\sigma }}}\\ I_2=-{\displaystyle \sum _{k=N_1+1}^{N_e}{a}_{k}j\omega \pi B_{k,0}^{\mathsf{\sigma }}}\end{array},$$

(19)

The outside sensor impedance variation will have a great deal of influence on the electric charge and electric potential distribution. As seen from Equation (18), every time the outside impedance varies, part of the coefficient of the system equations (C₄ and C₅) will change. Thus, the whole system equations must be solved once, which requires a long time for solving the large dimension linear equations. A new strategy was presented in next section.

5. Modification of the Non-Periodic Model

In the above method, the SAW sensor was calculated as a unitary device (shown in Figure 6a). In our method, the SAW sensor was treated in a two-port network and the SAW transponder and the loaded impedance sensor part were calculated as a form of the ABCD matrix separately [18] (Figure 6b). This route can isolate the time-consuming large sparse equation-solving task with varying boundary conditions. The ABCD matrix [18] of the SAW chip can be derived using the non-periodic FEM/BEM model under two electrical boundary conditions: Short circuit and open circuit.

When we set the output IDTs as short circuit and open circuit, we can calculate the solutions of Equation (18) for the two situations. Then, we can obtain the expression of the ABCD matrix at all frequency points according to Equation (20). These values in the ABCD matrix are set because it represents the device without considering the external boundary conditions.

$$\left[\begin{array}{c}{U}_1\\ I_1\end{array}\right]=\left[\begin{array}{cc}A& B\\ C& D\end{array}\right]\left[\begin{array}{c}{U}_2\\ I_2\end{array}\right],$$

(20)

When the end of the output IDT is loaded with an impedance Z_L, we can get $I_2=U_2/Z_L$. Substituting $I_2=U_2/Z_L$ into Equation (20), we get

$$\{\begin{array}{c}{I}_1=(C+D/Z_L)U_2\\ U_1=(B/Z_L+A)U_2\end{array},$$

(21)

The input impedance of the U₁ port is $Y=\frac{I_1}{U_1}=\frac{C+D/Z_L}{B/Z_L+A}$. Then, the scattering parameters of this special two-port device can be obtained using the conventional formulae.

Using this method, Equation (18) can be solved twice which reduces computation time. This method can also be used in the situation where an impedance matching circuit has an influence on the electrical boundary conditions of the IDT, either at the front or the back ends of the device.

6. Experimental Verification

A test device was manufactured to verify the simulated results. The device parameters are as follows: The substrate of the sample is 128° YX-LiNbO₃. The electrode is aluminium. The periods of the transmitting/receiving IDT and the programmable reflector are both 12.3 µm. The transmitting/receiving IDT has a single finger and the number of finger pairs is 20. The programmable reflector is of the split finger type and the number of finger pairs here is also 20. The metallization ratio of both IDTs is 0.5 and the film thickness is 0.2 µm. The aperture of the IDTs is 421.7 µm. The gap between the IDTs is 3631.8 µm.

Using the methods mentioned above, the device performance can be simulated easily. The S₁₂ characteristic of the sensor was thus calculated and compared with the measured results in Figure 7.

It can be seen that the profiles of admittance in the simulation show good similarity with the experimental results. However, the amplitudes at some points are different.

The parameters derived from the above figure are shown in

Table 1. Results analysis between the simulation and measurement.

Results Type	Fc (MHz)	Insertion Loss (dB)	3 dB Bandwidth (MHz)	Band Ripple (dB)
Simulation	317.12	8.78	14.2	4.2
Measurement	316.83	9.16	15.1	1.98
Difference	0.09%	−4.1%	−5.9%	52.9%

. The center frequency and bandwidth of the simulation agree well with the measurements. The band ripple of the simulation is larger than the measured values. This is because the test device was measured on a test fixture, where the bond wire effect was not considered in the simulation.

The simulated results, including the loaded matching circuit and capacitive sensor are shown in Figure 8. The return loss is the amplitude of the response pulse in the time domain. The return loss variation curve shows a resonant behavior as presented before [15]. This property was verified by several groups [19,20,21]. The calculated curve corresponds with the trend of the measured points. In the range of 7–12 pF, there is a small ripple in the measured points which can also be seen clearly in the magnified view of the simulated curve. The simulated results are approximately 5 dB lower when the capacitance is 1 pF and 7.5 dB lower for a 2 pF value. This may be because when the capacitance value is small the bond wire effect can easily affect the measured results.

In the last paragraph the bond wire effect was not considered and the profile of the simulated and measured S₁₂ curve corresponds perfectly. On the other hand, the simulation result has larger errors for the smaller capacitance value. Both results came to a conclusion that the bond wire effect played an important role in the simulation results, especially for the small impedance value. The next step for our research will be taking this factor into our consideration in the simulation process.

Moreover, the result shows a resonance behavior of the return loss as a function of capacitance of the sensor. The sensor operation could be accepted only on a portion of the capacitance, i.e., in the region of 1–4 pF or 4–7 pF. For the sensing region below 4 pF, the capacitive sensor only can be integrated on the same chip with the SAW sensor to minimize the parasitic capacitance of connecting wires. In addition, it is hard to find a commercial sensor with appropriate capacitance variation. Thus, the need for sensor optimization is urgent in order to shift the monotonous range by trying different sensor configurations. Fast simulation methods will play very important roles in the optimization process.

7. Conclusions

The finite transducer analysis method based on the FEM/BEM model, which is one of the most accurate methods for simulating SAW filters, was successfully applied to simulate the impedance-loaded SAW sensors. An accurate and fast algorithm for calculating the bulk wave part of the equation coefficient was presented to reduce the required computing time. Modification to the non-periodic FEM/BEM model was made to adapt to the impedance-loaded sensors. Thus, we only need to calculate the short circuit and open circuit situations of the SAW sensors respectively in order to obtain the sensor response at each value of loaded impedance. The simulated results were compared with the measured results of the test device, which show good accordance.