Recent years have seen a significant expansion in the field of Micro-Electro-Mechanical System (MEMS). This growth can be attributed to the system’s capacity to integrate mechanical and electrical components at a microscopic scale, thereby enabling miniaturized devices with diverse applications across multiple disciplines [1–3].

In the field of mathematics, the development of MEMS has led to significant advancements in advanced modelling and simulation techniques [4–6]. These techniques assist researchers in the creation of complex mathematical models to describe the behavior of MEMS devices, taking into account factors such as micro-scale fluid dynamics, mechanical vibrations and electrical interactions. The development of such models is of paramount importance for the optimization of MEMS device design and performance.

In the domain of nanotechnology, Micro-Electro-Mechanical System (MEMS) play a pivotal role in the fabrication of nanoscale structures. The integration of nanowires and nanotubes into MEMS sensors has been demonstrated to enhance their sensitivity and selectivity, thereby resulting in highly efficient chemical and biological sensing devices [7, 8].

The field of material science has undergone significant transformations in the wake of MEMS. The demand for specific material properties has driven the development of new materials for MEMS applications. For instance, enhanced piezoelectric materials have been shown to facilitate more efficient energy harvesting and actuation in MEMS devices [9, 10].

MEMS oscillators are of particular importance within the domain of MEMS components [11]. Their utilization is pervasive in energy harvesters [12] and accelerometers [13]. In biosensors, they facilitate precise detection and analysis of biological substances by converting biomolecule-sensor surface interactions into electrical signals. In the context of wearable sensors, they enable uninterrupted monitoring of physiological parameters, including body movements and heart rate.

However, in the research of MEMS oscillators, several aspects require consideration. The periodic property is pivotal in determining the oscillator’s stable operation and frequency characteristics. It is imperative to comprehend the phenomenon of pull-in instability, as its onset can culminate in device failure when the oscillator approaches this state [14, 15]. Additionally, energy conversion efficiency affects power consumption, and temperature stability is vital due to the often varying operating temperatures of MEMS oscillators.

The harmonic balance method is a common technique for analyzing nonlinear systems [16, 17]. This method involves the representation of the response of a nonlinear system as a sum of harmonic functions, thereby providing a relatively simple means by which to analyze complex nonlinear behavior. In the context of MEMS oscillators, this method can occasionally yield analytical solutions, which are valuable for understanding the relationship between system parameters and oscillator performance.

However, it is important to note that the harmonic balance method has limitations when applied to MEMS oscillators. MEMS oscillators exhibit singular terms, which complicates the enhancement of solution accuracy in the vicinity of singularities. Furthermore, the method’s outputs are found to be highly sensitive to harmonic components, necessitating the inclusion of additional components to achieve higher levels of accuracy, which in turn increases the complexity of the calculations. Furthermore, as the applied voltage in MEMS oscillators changes, the system may shift from periodic motion to pull-in instability [18–20], further complicating the analysis.

It is evident that alternative methods, including the homotopy perturbation method [21, 22], the variational iteration method [23, 24], the variational principle [25, 26] and He’s frequency formulation [27, 28], have demonstrated efficacy in addressing singularities. These methods offer different ways to approximate solutions of nonlinear systems and may overcome the limitations of the harmonic balance method in MEMS oscillator analysis. The present study aims to explore the application and limitations of the harmonic balance method in MEMS oscillator analysis and to compare it with alternative methods, thereby providing insights for further research in this field.

The harmonic balance method is a well-established and powerful technique in engineering and physics, serving as a fundamental tool for analyzing nonlinear systems. The underlying principle of this method is based on the representation of the response of a nonlinear system as a sum of harmonic functions.

For a general nonlinear system described by an ordinary differential equation of the form: w ″ + f w , w ′ , w ″ = 0 (1) where f is a nonlinear function. Equation 1 has periodic solution when f/w > 0, and we assume the periodic solution has the following form: w = ∑ n = 0 N { a n ⁡ cos ⁡ ( n ω t ) + b n ⁡ sin ⁡ ( n ω t where ω is the frequency, a n and b n are coefficients to be determined late.

This method offers several advantages. It provides a relatively straightforward and intuitive way to analyze nonlinear systems. By decomposing the complex nonlinear behavior into individual harmonic components, engineers and researchers can more easily understand and study the fundamental characteristics of nonlinear oscillators. Specifically, it enables the analysis of frequency response and amplitude-related behavior independently. In numerous instances, the harmonic balance method has been shown to yield analytical solutions. These solutions establish explicit relationships between system parameters (e.g., capacitance, stiffness, and applied voltage in MEMS oscillators) and solution characteristics (e.g., amplitude and frequency). This enables a deeper understanding of how changes in these parameters affect the oscillator’s performance. To illustrate this point, consider the analytical expression for the amplitude of a MEMS oscillator in terms of capacitance and applied voltage. This expression can reveal the precise relationship between an increase in voltage and an increase in amplitude.

However, when applied to MEMS oscillators, the harmonic balance method also has some notable drawbacks. Firstly, MEMS oscillators contain singular terms. Conventional numerical techniques frequently employed in the harmonic balance method encounter difficulties in dealing with these singularities, thereby hindering the enhancement of the solution’s accuracy in the vicinity of these points. This can result in substantial errors in the approximation of the oscillator’s behavior. Secondly, the accuracy of the results obtained using the harmonic balance method is highly sensitive to the number of harmonic components considered. Achieving a high level of accuracy in the solution often necessitates the consideration of a substantial number of harmonic components. This not only increases the complexity of the calculation but also makes the computational process more time-consuming. Thirdly, it is important to note that MEMS oscillators exhibit a wide range of dynamical properties. As the voltage increases gradually, the system may transition from periodic motion to pull-in instability. The harmonic balance method may not be able to accurately capture this complex transition behavior, further limiting its effectiveness in analyzing MEMS oscillators.

The present paper focuses on a specific MEMS oscillator, the governing equation of which is given [4]: w ″ + w = k 1 − w , w 0 = 0 , w ′ 0 = 0 (2)

In this equation, w represents the dimensionless displacement and k is the voltage parameter. This equation is singular, and the presence of this singularity serves to render the problem considerably more complex [29].

In order to apply the harmonic balance method [16, 17], it is necessary to treat with the singular term. To this end, Equation 2 is rewritten in the form: w ″ + w − w ″ w − w 2 − k = 0 (3)

This transformation is a crucial initial step, as it aligns the equation with the principles of the harmonic balance method, thereby enabling the use of the approximation techniques of this method.

Subsequent to the establishment of the initial conditions of the oscillator, it is assumed that the solution may be expressed in the following form: w = A ⁡ sin 2 ⁡ ω t (4) where A and ω are unknown for further determination.

Substituting Equation 4 into Equation 3 yields the residual function R: R t = w ″ + w − w ″ w − w 2 − k = 2 A ω 2 ⁡ cos ⁡ 2 ω t + A ⁡ sin 2 ⁡ ω t − 2 A 2 ω 2 ⁡ cos ⁡ 2 ω t ⁡ sin 2 ⁡ ω t − A 2 ⁡ sin 4 ⁡ ω t − k = 2 A ω 2 − k + A 1 − 6 ω 2 sin 2 ⁡ ω t + A 2 4 ω 2 − 1 sin 4 ⁡ ω t which represents the net force imbalance at each phase ( ω t ). To satisfy the harmonic balance principle, we evaluate its residual function R at ω t = 0 and ω t = π , the points corresponding to the oscillator’s equilibrium position and maximum displacement, where the linear stiffness and nonlinear electrostatic forces are most dominant. Setting R = 0 at these locations, we obtain 2 A ω 2 − k = 0 (5) A 1 − 6 ω 2 + A 2 4 ω 2 − 1 = 0 (6)

Equations 5, 6 enforce force balance at critical phases, allowing us to solve for the unknown coefficients A and ω by projecting the nonlinearity onto the fundamental harmonic component: ω = 1 − A 6 − 4 A A = 1 + 2 k − 1 + 2 k 2 − 12 k 2 which provide explicit expressions for the coefficients in terms of the system parameters.

In this study, we will be comparing the approximate solutions obtained using the harmonic balance method with the exact solutions. As demonstrated in Figure 1, the figure reveals a substantial discrepancy between the approximate and exact solutions for different values of k. This discrepancy can be attributed primarily to the selection of the trial solution, as outlined in Equation 4, and the handling of the singularity during the solution process.

For the purpose of comparison, an alternative approach to treating the singular term is hereby presented. Utilizing the Taylor series, the following equation can be expressed as: w ″ + w − k 1 + w + w 2 + w 3 + w 4 + ⋯ = 0 (7)

Substituting Equation 4 into Equation 7 results in the following residual R t = 2 A ω 2 ⁡ cos 2 ω t + A ⁡ sin 2 ⁡ ω t − k 1 + A ⁡ sin 2 ⁡ ω t + A 2 ⁡ sin 4 ⁡ ω t + A 3 ⁡ sin 6 ⁡ ω t + A 4 ⁡ sin 8 ⁡ ω t + ⋯

By a similar operation as above, we have R t = 0 = 2 A ω 2 − k = 0 (8) R t = π / 2 ω = − 2 A ω 2 + A − k 1 + A + A 2 + A 3 + A 4 + ⋯ = 0 (9)

Equation 8 is equivalent to Equation 5, and Equation 9 can be expressed as follows: − 2 ω 2 + 1 − k A 1 − A = 0 (10)

By simultaneously solving Equations 8, 10, we obtain ω 2 = 1 − A 2 2 − A (11) k = A 1 − A 2 − A

Based on Equation 11, when A = 1, the pull - in instability emerges and the pull - in voltage is k = 0, a situation that deviates significantly from real - world conditions. Clearly, this method of dealing with the singular term is ineffectual. Given this shortcoming, we shift our focus to selecting an alternative approximate solution in the following form: w = B 1 − cos ⁡ ω t (12) where B and are unknown for further determination.

The alternative trial solution is proposed to address the limitations of the harmonic approximation of Equation 4, w = A ⁡ sin 2 ⁡ ω t = A 1 − cos ⁡ 2 ω t / 2 , though it is simple, and meets the initial conditions with a higher harmonic component ( cos ⁡ 2 ω t ), while Equation 12 has lower harmonics with also a simple harmonic component ( cos ⁡ ω t ). This simple modification is critical for improvement of its accuracy.

By a similar operation as above, we obtain the following residual: R t = B ω 2 ⁡ cos ⁡ ω t + B 1 − cos ⁡ ω t − B 2 ω 2 ⁡ cos ⁡ ω t 1 − cos ⁡ ω t − B 2 1 − cos ⁡ ω t 2 − k

We set R = 0 at two location points R 0 = 0 and R t = π / 3 ω = 0 , this results in B − B 2 − k = 0 ω 2 = 1 − 2 B + 1 2 B 1 + 1 2 B = 1 − 3 2 B 1 + 1 2 B

Figure 2 shows the comparison between the new approximate solution and the exact solutions for different values of k. It is evident that the selection of the trial solution and the treatment of singular terms are pivotal in the harmonic balance method. As discussed in the homotopy perturbation method [30], an initial guess that is well-founded invariably yields an ideal result. The singularity treatment ensures the reliability of the solution process.

Pull-in instability occurs when ω 2 = 0 , this results in B = 2/3, and the pull-in voltage k = 2/9, while the exact pull-in voltage is k = 0.20363, the relative error is 9.13%. This error may be acceptable for some engineering applications, but for more precise requirements, further improvement of the analysis method is needed. We can increase the accuracy by considering more harmonic components in the approximate solution or using more advanced numerical techniques to handle the singularities.

To improve the accuracy, the trial solution should include more harmonical components: w = ∑ n = 0 N B n ⁡ cos n ω t (13)

In view of the initial conditions, the coefficients in Equation 13 satisfy the following condition ∑ n = 0 N B n = 0

We can choose N location points to determine B n and ω , this can greatly improve the accuracy, but it also increases the complexity of the calculation.

The present work makes two distinct contributions to the analysis of MEMS oscillators with singularities. First, unlike traditional approaches that employ homotopy perturbation [31, 32] or variational iteration method [33–35] or numerical methods [36, 37] to address nonlinearity, we demonstrate that the strategic trial function selection in the harmonic balance method can mitigate singularity-induced errors with minimal computational overhead. This offers a design-friendly alternative for rapid parametric analysis, as shown by the 9.13% improvement in pull-in voltage prediction compared to the naive sinusoidal trial solution. The accuracy can be improved if more harmonic components are involved in the trial function. Second, we provide a suitable treatment of the singularity that quantifies the trade-off between trial function complexity and accuracy, a critical gap in prior HBM applications to MEMS. While methods like HPM may achieve higher precision, our approach prioritizes analytical transparency and computational efficiency, making it suitable for preliminary design stages or educational contexts. The insights herein thus expand the utility of HBM in MEMS research and establish a benchmark for future hybrid analytical-numerical methodologies.

Unlike prior studies that rely on complex transformations (e.g., homotopy or series expansions), we demonstrate that simply modifying the trial solution form can significantly reduce singularity-induced errors. This approach retains HBM’s analytical simplicity while improving accuracy, as shown in Figure 2.

In the in - depth study of MEMS oscillators, this paper has presented a detailed comparison between the harmonic balance method and the homotopy perturbation method. MEMS oscillators, with their complex structures and unique operating characteristics, pose significant challenges in the field of analysis.

The harmonic balance method, which is widely applied in engineering and physics, has limitations when used in the analysis of MEMS oscillators. Due to the complexity of MEMS oscillators, the harmonic balance method faces issues such as parameter errors and singularities, which significantly affect the accuracy of the results. Parameter errors stem from various sources. Manufacturing process variations can change the physical dimensions and material properties of oscillator components. For example, minor changes in the thickness of the vibrating structure or the dielectric constant of capacitive elements can notably impact the oscillator’s performance. Environmental factors, especially temperature fluctuations, also play a crucial role. As the temperature changes, the Young’s modulus of the structural materials varies, affecting the oscillator’s stiffness, resonant frequency, and other characteristics. These errors accumulate during the analysis process, making it difficult to accurately predict the behavior of MEMS oscillators.

The homotopy perturbation method shows great potential as an alternative. It can effectively handle the singularities in MEMS oscillators. By leveraging the concept of homotopy, it transforms complex problems into simpler ones, circumventing the difficulties posed by singularities and obtaining more accurate and reliable solutions. This helps in a deeper understanding of the behavior of MEMS oscillators, providing support for optimizing their design and performance. For instance, in the design of high - precision MEMS - based sensors, accurately predicting the oscillator behavior can enhance the sensor’s sensitivity and stability.

However, although the homotopy perturbation method has obvious advantages in dealing with singularities, the harmonic balance method still has its value. The harmonic balance method is intuitive when analyzing the basic characteristics of MEMS oscillators. In some cases, it can yield analytical solutions, facilitating a quick understanding of the relationship between system parameters and oscillator performance. Thus, it can serve as a reference in the initial stage of analysis or when a rough estimate is sufficient.

Looking ahead, this study motivates several avenues for advancing MEMS oscillator analysis. First, refining the harmonic balance method through higher-order trial functions or adaptive singularity modeling could further reduce errors in pull-in voltage predictions, particularly for systems operating near critical gaps. Second, developing hybrid HBM-numerical frameworks—such as using HBM solutions as preconditions for finite element simulations—holds promise for balancing computational efficiency and accuracy in complex geometries. Third, extending the methodology to multi-physics MEMS systems (e.g., thermal or fluid-coupled oscillators [38–40]) will address real-world operational challenges, such as temperature-induced drift. Finally, integrating machine learning [41], AI-powered problem solving [42] and ResNet Neural Network [43] and advanced numerical simulation methods [44] with analytical models could create self-calibrating frameworks that adapt to manufacturing tolerances and experimental data. These advancements will strengthen the harmonic balance method’s utility across precision engineering, nanotechnology, and bio-MEMS applications.