A Novel Method for Estimating and Balancing the Second Harmonic Error of Cylindrical Fused Silica Resonators

The cylindrical resonator gyroscope (CRG) is a type of Coriolis vibratory gyroscope which measures the angular velocity or angle through the precession of the elastic wave of the cylindrical resonator. The cylindrical fused silica resonator is an essential component of the CRG, the symmetry of which determines the bias drift and vibration stability of the gyroscope. The manufacturing errors breaking the symmetry of the resonator are usually described by Fourier series, and most studies are only focusing on analyzing and reducing the fourth harmonic error, the main error source of bias drift. The second harmonic error also is one of the obstacles for CRG towards high precision. Therefore, this paper provides a chemical method to evaluate and balance the second harmonic error of cylindrical fused silica resonators. The relation between the frequency split of the n = 1 mode and the second harmonic error of the resonator is obtained. Simulations are performed to analyze the effects of the first three harmonic errors on the frequency splits. The relation between the location of the low-frequency axis of n = 1 mode and the heavy axis of the second harmonic error is also analyzed by simulation. Chemical balancing experiments on two fused silica resonators demonstrate the feasibility of this balancing procedure, and show good consistency with theoretical and simulation analysis. The second harmonic error of the two resonators is reduced by 86.6% and 79.8%, respectively.


Introduction
The Coriolis vibratory gyroscope measures angular velocities or angles by the precession of the standing wave [1]. The cylindrical resonator gyroscope (CRG) and the hemispherical resonator gyroscope (HRG) are two types of Coriolis gyroscopes, and they are widely used in inertial navigation systems. HRGs have unparalleled precision guaranteed by a highly symmetrical resonator with extremely small dissipation. The HRG Crystal TM DUAL CORE made by Safran Electronics & Defense (Massy, France) has achieved a bias stability of 0.0001 • /h over 100 h [2]. Inspired by these breakthroughs, HRGs have become a topic of interest within the inertial community, with devoted researchers investigating traditional size HRGs [3][4][5][6] and micro-version HRGs [7][8][9]. However, the complexity of the manufacturing of HRGs has limited the mastering of this technology to a small number of manufacturers.
Compared to HRGs, the CRGs are simpler to manufacture. Similar to HRGs, the key component of a CRG is the cylindrical resonator, the quality factor and the residual imperfections of which generally determine the performance of the CRG. Many researchers have been devoted to improving the quality factor of cylindrical resonators through the selection of resonator material, optimizing the structure, improving the processing technology and applying special post-processing treatments [10][11][12][13]. Based on these improvements, the quality factor of cylindrical fused silica resonator has achieved 2.89 million and still

The Influence of the Unbalanced Mass on Quality Factor
A cylindrical resonator with a quality factor around 5 million is adequate for the performance up to 0.01 • /h. The real obstacle in the way of inertial grade cylindrical resonator gyroscopes comes from the unbalanced mass. As the first three harmonics of the unbalanced mass lead to motion of the center of mass, the real quality factor depends on the properties of the connection between the resonator and the base [31]. The quality factor Q u induced by the unbalanced mass is given by: where m u is the unbalanced mass of the resonator, M d is the modal mass of the resonator, ω h and Q h is the frequency and the quality factor of the harmonic mode respectively, ω 2 is the natural frequency of the n = 2 mode. According to Equation (1), the unbalanced mass and the characteristics of fixation between resonator and base both affect the quality factor Q u .
The real quality factor Q r for a assembled resonator is then: where Q n is the nominal quality factor of the n = 2 mode. For the resonator with nominal quality factor Q n over one million, the real quality factor Q n approximately equal to unbalanced quality factor Q u . The position of the standing wave hinges on the central angle of the second harmonic error under external vibrations, and the increasing rate of the vibration amplitude of the standing wave is proportional to the value of the second harmonic error. The value and direction of the first and the third harmonics can be identified by measuring the amplitude and the central angle of the standing wave under transverse vibration of the base [32]. Due to the presence of the first three harmonic errors, the resonator gyroscope is sensitive to external vibrations causing the deterioration in the bias drift and random walk.

Relationship between the Second Harmonic Error and the Frequency Spilt of the n = 1 Mode
The effective sensing element of a cylindrical resonator is the resonator shell. Therefore, the ring model is applied to describe the characteristics of a cylindrical resonator.
The cross-section of a thin circular ring with axial length L is shown in Figure 1. The mass point m i is attached at the central angle ϕ i . The imperfections of the resonator are modeled by the circumferentially varied density. The density of a resonator with imperfections can be expressed by Fourier series and only take the first three harmonic errors in the case: where ρ is the equivalent density of the imperfect resonator, ρ 0 is the density of the perfect resonator, ρ i is the i-th harmonic error of the imperfections, and ϕ is the central angle.
The natural frequency of the cylindrical resonator in the n-th mode is given by [33]: where ω n1 , ω n2 is the natural frequency of the n-th mode for an imperfect resonator, ω n0 is the natural frequency of a perfect resonator, m i is the attached point mass, ϕ i is the central angle of point mass m i , M 0 is the mass of the perfect resonator, and α n is the amplitude ratio of the radial and the tangential vibration. The natural frequency of the cylindrical resonator in the n-th mode is given by [33]: where ωn1, ωn2 is the natural frequency of the n-th mode for an imperfect resonator, ωn0 is the natural frequency of a perfect resonator, mi is the attached point mass, φi is the central angle of point mass mi, M0 is the mass of the perfect resonator, and αn is the amplitude ratio of the radial and the tangential vibration. The frequency split is defined as the frequency difference between the natural frequencies of the n-th mode. Therefore, the frequency split can be calculated from Equation (4) (5) where Δωn is the frequency split of the n-th mode in rad/s. As the number of discrete mass points increases, the mass defect can be expressed as an integration over the central angle φ:  (6), we can obtain the relationship between ρ2 and the frequency split of the n = 1 mode Δω1: Note that the frequency split of the n = 1 mode Δω1 is proportional to the second harmonic error ρ2. In other words, we can reduce the second harmonic error by frequency tuning of the n = 1 mode. Once the frequency split of the n = 1 mode decreases, the second harmonic error decreases accordingly. The frequency split of the n = 1 mode is readily measurable by the system developed for the measurement of the n = 2 mode, therefore, this new method avoids developing a special measurement system for the determination of the second harmonic error by vibration coupling. This novel method provides an effective and practical approach to identify the second harmonic error of the resonator which makes the balancing possible. The frequencies of the resonator in the FEM analysis and The frequency split is defined as the frequency difference between the natural frequencies of the n-th mode. Therefore, the frequency split can be calculated from Equation (4): where ∆ω n is the frequency split of the n-th mode in rad/s. As the number of discrete mass points increases, the mass defect can be expressed as an integration over the central angle ϕ: Substituting Equation (3) into Equation (6), we can obtain the relationship between ρ 2 and the frequency split of the n = 1 mode ∆ω 1 : Note that the frequency split of the n = 1 mode ∆ω 1 is proportional to the second harmonic error ρ 2 . In other words, we can reduce the second harmonic error by frequency tuning of the n = 1 mode. Once the frequency split of the n = 1 mode decreases, the second harmonic error decreases accordingly. The frequency split of the n = 1 mode is readily measurable by the system developed for the measurement of the n = 2 mode, therefore, this new method avoids developing a special measurement system for the determination of the second harmonic error by vibration coupling. This novel method provides an effective and practical approach to identify the second harmonic error of the resonator which makes the balancing possible. The frequencies of the resonator in the FEM analysis and our experiments were given in Hz, therefore, we calculate the frequency split of the n-th mode ∆f n in Hz, and ∆f n = ∆ω n /2π.

FEM Analysis of the First Three Harmonic Errors
It is widely recognized that the 2k-th harmonic error mainly contributes to the frequency split of the n = k mode. FEM analysis on our structured resonators is employed to investigate the quantitative relationship between the first three harmonic errors and the n-th (n = 1, 2, 3) mode. The half-section view of the cylindrical fused silica resonator is shown in Figure 2. The resonator shell is the main sensing element whose imperfections have a significant impact on the gyroscope performance. The suspension and the bottom plate is the elastic suspension, the stiffness of which can be adjusted by changing the thickness or the size of the eight holes on the bottom plate for different designs. The inner stem is used to mount the resonator on a base [12].

FEM Analysis of the First Three Harmonic Errors
It is widely recognized that the 2k-th harmonic error mainly contributes to the frequency split of the n = k mode. FEM analysis on our structured resonators is employed to investigate the quantitative relationship between the first three harmonic errors and the n-th (n = 1, 2, 3) mode. The half-section view of the cylindrical fused silica resonator is shown in Figure 2. The resonator shell is the main sensing element whose imperfections have a significant impact on the gyroscope performance. The suspension and the bottom plate is the elastic suspension, the stiffness of which can be adjusted by changing the thickness or the size of the eight holes on the bottom plate for different designs. The inner stem is used to mount the resonator on a base [12]. The parameters of the material and the structure of the resonator are listed in Table  1. The modes of vibration of a resonator without density errors is simulated, and the n = 1, 2 and 3 modes are shown in Figure 3a. Note that even the resonator without density errors has a frequency split around 0.1 Hz due to the mesh error of an irregular bottom plate. This frequency split is caused by the quality of the grids in the FEM model. Therefore, the zero-closed frequency splits which smaller than 0.1 Hz stems from the grids instead of the harmonic errors. The parameters of the material and the structure of the resonator are listed in Table 1. The modes of vibration of a resonator without density errors is simulated, and the n = 1, 2 and 3 modes are shown in Figure 3a. Note that even the resonator without density errors has a frequency split around 0.1 Hz due to the mesh error of an irregular bottom plate. This frequency split is caused by the quality of the grids in the FEM model. Therefore, the zero-closed frequency splits which smaller than 0.1 Hz stems from the grids instead of the harmonic errors.  According to Equation (3), the density of the resonator is set as a continuous function of central angle φ in the FEM model. To investigate the impact of the first three harmonic errors on frequency splits separately, three functions of nonuniform density are employed as follows: The degree of imperfections Δρi is ρi/ρ0. The circumferential distribution of imperfect density on resonator when Δρi = 0.05% is given in Figure 2. The distribution cycle of density with the first three harmonic errors is 2π, π, and 2π/3, respectively. According to Equation (3), the density of the resonator is set as a continuous function of central angle ϕ in the FEM model. To investigate the impact of the first three harmonic errors on frequency splits separately, three functions of nonuniform density are employed as follows: The degree of imperfections ∆ρ i is ρ i /ρ 0 . The circumferential distribution of imperfect density on resonator when ∆ρ i = 0.05% is given in Figure 2. The distribution cycle of density with the first three harmonic errors is 2π, π, and 2π/3, respectively.
The effects of the first three harmonics of density on the frequency split of the n = 1, 2, and 3 modes are investigated. As shown in Figure 4, the frequency split of the n = 1 and the n = 3 mode increases with the increasing first harmonic error, while the frequency split of the n = 2 mode remains small (around 0.0095 Hz) when the value of imperfect density is lower than 100 kg/m 3 and ∆ρ i = 4.5%. Moreover, the growth rate of the frequency split of the n = 1 mode (shown in grey line) is much faster than that of the n = 3 mode (shown in blue line). The first harmonic error has a larger impact on n = 1 mode than the other two, and the frequency split of n = 1 mode is approximately proportional to the square of the first harmonic error when the second harmonic error is more than 15 kg/m 3 , as shown in the curve fit of ∆f 1 . errors on frequency splits separately, three functions of nonuniform density are em as follows: The degree of imperfections Δρi is ρi/ρ0. The circumferential distribution of im density on resonator when Δρi = 0.05% is given in Figure 2. The distribution cycle sity with the first three harmonic errors is 2π, π, and 2π/3, respectively.
The effects of the first three harmonics of density on the frequency split of t 2, and 3 modes are investigated. As shown in Figure 4, the frequency split of the n the n = 3 mode increases with the increasing first harmonic error, while the freque of the n = 2 mode remains small (around 0.0095 Hz) when the value of imperfect is lower than 100 kg/m 3 and Δρi = 4.5%. Moreover, the growth rate of the frequen of the n = 1 mode (shown in grey line) is much faster than that of the n = 3 mode in blue line). The first harmonic error has a larger impact on n = 1 mode than t two, and the frequency split of n = 1 mode is approximately proportional to the s the first harmonic error when the second harmonic error is more than 15 kg/m 3 , a in the curve fit of Δf1.   Figure 5 shows that the frequency split of the n = 1 mode increases rapidly with the growth of the second harmonic error, while that of the n = 2 and the n = 3 mode increase only mildly. For example, the frequency split of n = 1 mode reaches 38.67 Hz when ρ 2 is 65 kg/m 3 , while the frequency split of the n = 2, 3 mode is only 1.19 Hz and 0.99 Hz, respectively. The frequency split of n = 1 mode is severely affected by ρ 2 compared with the other two modes, even when the defect is small. Moreover, ∆f 1 grows linearly with the increase of ρ 2 , as illustrated by the grey line and the least-square fit. Therefore, this simulation result is consistent with (5). Nearly 95% of the frequency split ∆f 1 comes from the second harmonic error when ∆ρ 2 is 2.95%. The frequency split ∆f 1 can be used for the evaluation of the second harmonic error. As we obtain the frequency split ∆f 1 , the value of the ρ 2 can be calculated from the fitted formula.
The simulation result shown in Figure 6 indicates that the frequency split of the n = 3 mode is proportional to the square of the third harmonic error. The influence of ρ 3 on the n = 1, 2 mode is rather small compared that on the n = 3 mode. The frequency split of the n = 1 mode keeps almost unchanged at 0.089 Hz. And the frequency split of n = 2 mode stay constant at 0.0096 Hz. Therefore, it is reasonable to consider that the third harmonic error has no impact on the frequency split of the n = 2, 3 mode. The frequency split of n = 3 mode is approximately proportional to the square of the third harmonic error when the second harmonic error is more than 40 kg/m 3 , as shown in the curve fit of ∆f 3 . other two modes, even when the defect is small. Moreover, Δf1 grows linearly with th increase of ρ2, as illustrated by the grey line and the least-square fit. Therefore, this simu lation result is consistent with (5). Nearly 95% of the frequency split Δf1 comes from th second harmonic error when Δρ2 is 2.95%. The frequency split Δf1 can be used for the eva uation of the second harmonic error. As we obtain the frequency split Δf1, the value of th ρ2 can be calculated from the fitted formula. The simulation result shown in Figure 6 indicates that the frequency split of the n = mode is proportional to the square of the third harmonic error. The influence of ρ3 on th n = 1, 2 mode is rather small compared that on the n = 3 mode. The frequency split of th n = 1 mode keeps almost unchanged at 0.089 Hz. And the frequency split of n = 2 mod stay constant at 0.0096 Hz. Therefore, it is reasonable to consider that the third harmon error has no impact on the frequency split of the n = 2, 3 mode. The frequency split of n 3 mode is approximately proportional to the square of the third harmonic error when th second harmonic error is more than 40 kg/m 3 , as shown in the curve fit of Δf3.   The simulation result shown in Figure 6 indicates that the frequency split of the n = 3 mode is proportional to the square of the third harmonic error. The influence of ρ3 on the n = 1, 2 mode is rather small compared that on the n = 3 mode. The frequency split of the n = 1 mode keeps almost unchanged at 0.089 Hz. And the frequency split of n = 2 mode stay constant at 0.0096 Hz. Therefore, it is reasonable to consider that the third harmonic error has no impact on the frequency split of the n = 2, 3 mode. The frequency split of n = 3 mode is approximately proportional to the square of the third harmonic error when the second harmonic error is more than 40 kg/m 3 , as shown in the curve fit of Δf3.    Figure 6 that the second harmonic error has the most significant impact on the frequency split of all three modes, especially the n = 1 mode. In fact, the second harmonic error accounts for the largest proportion of ∆f 1 . The ∆ρ i is usually small, no larger than 5%. As shown in Figure 6, the frequency split of the three modes is 19.352 Hz, 0.296 Hz and 0.592 Hz, respectively, when the second harmonic error is 32 kg/m 3 (where ∆ρ 2 is 1.45%). In contrast, the first and the third harmonic error have little influence on ∆f 1 . The frequency split ∆f 1 caused by the second harmonic error is 133 times and 215 times that of the first and third harmonic error when ∆ρ i is 1.45%. As the harmonic error reaches 65 kg/m 3 (where ∆ρ i is 2.95%), the frequency split ∆f 1 under harmonic errors ρ 1 , ρ 2 , and ρ 2 is 0.374 Hz, 38.672 Hz and 0.091 Hz, respectively. Hence, ρ 1 and ρ 3 make far less contribution to ∆f 1 compared with ρ 2 . Therefore, it is convenient to determine the value of the second harmonic error by measuring the frequency split of the n = 1 mode. harmonic error reaches 65 kg/m 3 (where Δρi is 2.95%), the frequency split Δf1 under harmonic errors ρ1, ρ2, and ρ2 is 0.374 Hz, 38.672 Hz and 0.091 Hz, respectively. Hence, ρ1 and ρ3 make far less contribution to Δf1 compared with ρ2. Therefore, it is convenient to determine the value of the second harmonic error by measuring the frequency split of the n = 1 mode.  Table 2 presents an overview of the qualitative impact of the harmonic errors of the imperfect density on the frequency splits of the three modes. The second harmonic error can be identified by the frequency split Δf1. However, there is no significant difference of the first three harmonic errors on the frequency split of the n = 2 and the n = 3 modes. The relations between the harmonic errors and the frequency split are also listed. There is a linear relation between the second harmonic error and the frequency split Δf1. The third harmonic error makes no significant impact on the frequency split of the n = 1 and the n = 2 mode. By contrast, the relation of the rest is quadratic or quadratic like.  Table 2 presents an overview of the qualitative impact of the harmonic errors of the imperfect density on the frequency splits of the three modes. The second harmonic error can be identified by the frequency split ∆f 1 . However, there is no significant difference of the first three harmonic errors on the frequency split of the n = 2 and the n = 3 modes. The relations between the harmonic errors and the frequency split are also listed. There is a linear relation between the second harmonic error and the frequency split ∆f 1 . The third harmonic error makes no significant impact on the frequency split of the n = 1 and the n = 2 mode. By contrast, the relation of the rest is quadratic or quadratic like. Table 2. Contribution of the first three harmonic errors on each mode when ρ i = 100 kg/m 3 and the relations between the frequency splits and the first three harmonic errors.

Mode
The

Chemical Balancing of the Second Harmonic Error
The chemical balancing procedure utilizes hydrofluoric acid to remove mass from fused silica resonators, which was first proposed by Basarab et al. [27,34] and extended in our previous work [30]. However, these studies have emphasized the chemical balancing on the fourth harmonic errors. Although Basarab et al. had presented brief examples for the treatment of other harmonic errors, the principles and properties of the chemical balancing on the second harmonic error need further investigation. In particular, the removal of the k-th harmonic error results in the generation of the ik-th (i = 2, 3, 4 . . . ) harmonics. This means the change in the second harmonic error also alters the fourth harmonic error, which may result in a larger frequency split of the n = 2 mode and dramatically degrade the performance of the resonator. Fortunately, by carefully selecting the chemical balancing parameters, it is possible to reduce the second harmonic error without introducing new harmonic errors.
Although the basic chemical trimming theory and its application on the fourth harmonic error balancing has been built and analyzed in our previous work [30], the properties of the chemical trimming theory in the second harmonic error case have not been illustrated and verified on cylindrical resonators. Therefore, we make efforts to figure out the relationship between the balancing mass and the balancing parameters.
As shown in Figure 8, the cylindrical fused silica resonator is partially soaked in the hydrofluoric acid at inclined angle α. Angle β is the included angle of the line segment OA and OB. R is the inner radius of the resonator shell, d is the thickness of the resonator shell, and h is the immersed depth of the resonator. Point O is the center of the resonator top surface. Point A is the lowest position of the immersed part, and point C is the intersection point of the line segment OA and the liquid level.
This means the change in the second harmonic error also alters the fourth harmonic error which may result in a larger frequency split of the n = 2 mode and dramatically degrad the performance of the resonator. Fortunately, by carefully selecting the chemical balanc ing parameters, it is possible to reduce the second harmonic error without introducing new harmonic errors.
Although the basic chemical trimming theory and its application on the fourth har monic error balancing has been built and analyzed in our previous work [30], the proper ties of the chemical trimming theory in the second harmonic error case have not been illustrated and verified on cylindrical resonators. Therefore, we make efforts to figure ou the relationship between the balancing mass and the balancing parameters.
As shown in Figure 8, the cylindrical fused silica resonator is partially soaked in th hydrofluoric acid at inclined angle α. Angle β is the included angle of the line segment OA and OB. R is the inner radius of the resonator shell, d is the thickness of the resonator shell and h is the immersed depth of the resonator. Point O is the center of the resonator top surface. Point A is the lowest position of the immersed part, and point C is the intersection point of the line segment OA and the liquid level. According to the multi-mode trimming theory [35][36][37], the second harmonic error could be balanced by removing mass from two positions starting from the heavy axis where the value of the second harmonic error reaches maximum. The mass element removed by chemical balancing from the two positions can be express as: where v is the etching rate of the hydrofluoric acid on fused silica in µm/min, C 0 (α, β), C 2 (α, β), C 4 (α, β) are the coefficients of Fourier series given in [30], and the coefficients higher than four is omitted in this case. As Equation (9) suggests, the removal of the second harmonic error could lead to the variation of the fourth harmonic error. However, it is possible to select the inclined angle α and the included angle β such that C 4 (α, β) = 0. Under this condition, the existing fourth harmonic error is immune to the balance of the second harmonic error. The desired parameters α, β lie on the intersecting lines of the zero-value surface in blue and the curved surface in red, as shown in Figure 8. There are three intersecting lines with numerous combinations of α and β. The larger angle β means the deeper immersed depth during chemical balancing, and the larger inclined angle α means that the resonator is more perpendicular to the liquid level. Taking both the experimental apparatus and the resonator structure into consideration, the angle α, β is finally set as 1.22 rad and 1.05 rad, which is marked green in Figure 9.
combinations of α and β. The larger angle β means the deeper immersed depth chemical balancing, and the larger inclined angle α means that the resonator is m pendicular to the liquid level. Taking both the experimental apparatus and the r structure into consideration, the angle α, β is finally set as 1.22 rad and 1.05 rad, marked green in Figure 9. In the practical balancing process, the resonator is etched at inclined angle α mersed depth h. According to balancing setup shown in Figure 7, the immersed d be expressed as: Taking the selected parameters α and β into Equation (10), the immersed d 2.25 mm in this case.

Determination of the Trimming Positions
One crucial problem to be solved is the determination of trimming positio balancing the second harmonic error. It has been demonstrated that the trimmi tions of the fourth harmonic error should be located at the principal axis of the lo ural frequency of the n = 2 mode [38]. However, the trimming positions of the n = has rarely been discussed in the open literature. In this section, simulation analys sented for investigation of the trimming positions of the n = 1 mode. In the practical balancing process, the resonator is etched at inclined angle α and immersed depth h. According to balancing setup shown in Figure 7, the immersed depth can be expressed as: Taking the selected parameters α and β into Equation (10), the immersed depth h is 2.25 mm in this case.

Determination of the Trimming Positions
One crucial problem to be solved is the determination of trimming positions when balancing the second harmonic error. It has been demonstrated that the trimming positions of the fourth harmonic error should be located at the principal axis of the lower natural frequency of the n = 2 mode [38]. However, the trimming positions of the n = 1 mode has rarely been discussed in the open literature. In this section, simulation analysis is presented for investigation of the trimming positions of the n = 1 mode.
As shown in Figure 10, angle σ is the circumferential location of the heavy axis, where the second harmonic error reaches the maximum. Angle ψ is the circumferential location of the low-frequency axis of the n = 1 mode, where the natural frequency reaches the minimum. Angle θ is the intersecting angle between the low-frequency axis and the heavy axis. We define the clockwise direction as the positive direction.
Micromachines 2021, 12, x FOR PEER REVIEW As shown in Figure 10, angle σ is the circumferential location of the heavy ax the second harmonic error reaches the maximum. Angle ψ is the circumferential of the low-frequency axis of the n = 1 mode, where the natural frequency reaches imum. Angle θ is the intersecting angle between the low-frequency axis and th axis. We define the clockwise direction as the positive direction. The relation between the angle ψ and σ is demonstrated by varying the ang then locating the low-frequency axis. Figure 11 illustrates that the low-frequency incides with the heavy axis. The angle of the low-frequency axis of the n = 1 mod portional to the angle of the heavy axis with a slope of 1.001. The intersecting The relation between the angle ψ and σ is demonstrated by varying the angle σ and then locating the low-frequency axis. Figure 11 illustrates that the low-frequency axis coincides with the heavy axis. The angle of the low-frequency axis of the n = 1 mode is proportional to the angle of the heavy axis with a slope of 1.001. The intersecting angle θ ranges from 0 to 1 degree when σ increases from 0 to 90 degrees. These slight variations of the intersecting angle θ mainly stem from simulation errors. Therefore, the defect mass could be removed from the direction of the low-frequency axis. Figure 10. The low-frequency axis of the n = 1 mode and the heavy axis of the second harmonic error.
The relation between the angle ψ and σ is demonstrated by varying the angle σ and then locating the low-frequency axis. Figure 11 illustrates that the low-frequency axis co incides with the heavy axis. The angle of the low-frequency axis of the n = 1 mode is pro portional to the angle of the heavy axis with a slope of 1.001. The intersecting angle θ ranges from 0 to 1 degree when σ increases from 0 to 90 degrees. These slight variations of the intersecting angle θ mainly stem from simulation errors. Therefore, the defect mass could be removed from the direction of the low-frequency axis. Figure 11. Relation between the heavy axis and the low-frequency axis of the n = 1 mode.

Results and Discussion
Chemical balancing experiments on two cylindrical fused silica resonators are performed. The frequency split and the direction of the low-frequency axis of the n = 1 mode are measured by the laser Doppler vibrometer (PSV-500). The resonator is excited by a loudspeaker. The axial vibration of the resonator shell is measured by the laser Doppler vibrometer. The natural frequency axis locates at the direction where the axial vibration reaches the maximum. Then, frequency sweeping is conducted in these directions. The direction of higher frequency is the high-frequency axis. The direction of lower frequency is the low-frequency axis. The frequency split is obtained by the frequency sweeping at the location which is the 22.5° (n = 2) or 45° (n = 1) apart from the low-frequency axis. The Figure 11. Relation between the heavy axis and the low-frequency axis of the n = 1 mode.

Results and Discussion
Chemical balancing experiments on two cylindrical fused silica resonators are performed. The frequency split and the direction of the low-frequency axis of the n = 1 mode are measured by the laser Doppler vibrometer (PSV-500). The resonator is excited by a loudspeaker. The axial vibration of the resonator shell is measured by the laser Doppler vibrometer. The natural frequency axis locates at the direction where the axial vibration reaches the maximum. Then, frequency sweeping is conducted in these directions. The direction of higher frequency is the high-frequency axis. The direction of lower frequency is the low-frequency axis. The frequency split is obtained by the frequency sweeping at the location which is the 22.5 • (n = 2) or 45 • (n = 1) apart from the low-frequency axis. The natural frequencies and frequency split of the n = 1 and n = 2 modes are measured respectively. The resonator is then transferred to the chemical balancing system. The inclined angle α and the immersed depth h are set as 1.22 rad and 2.25 mm according to theoretical analysis. The balancing time t depends on the etching rate v and the residual frequency split of the n = 1 mode. The resonator is lifted by the electronic platform after being etched for t/2 period, and is rotated 180 degrees for another t/2 period at the next trimming position.The detailed balancing process on resonator R01 is shown in Figure 12. Three rounds of chemical balancing experiments are performed on resonator R01. The frequency split ∆f 1 linearly decreases with the total balancing time t. The frequency split ∆f 2 ranges from 0.075 Hz to 0.098 Hz during the balancing process.
As shown in Figure 13a, the initial natural frequency of the n = 1 mode of resonator R01 is 3223.438 Hz and 3088.526 Hz, respectively, with an initial frequency split of up to 134.912 Hz, which implies the second harmonic error of resonator R01 is large. After several rounds of chemical balancing, the lower natural frequency increases to 3291.431 Hz and the higher natural frequency increases to 3309.424 Hz. The frequency split of the n = 1 mode decreases to 17.993 Hz with a drop of 86.6% compared to the original frequency split. Meanwhile, the vibration velocity at the low-frequency axis increases markedly from 0.433 µm/s to 11.54 µm/s, and the vibration velocity at the high-frequency axis increases from 1.327 µm/s to 5.16 µm/s. 134.912 Hz, which implies the second harmonic error of resonator R01 is large. After several rounds of chemical balancing, the lower natural frequency increases to 3291.431 Hz and the higher natural frequency increases to 3309.424 Hz. The frequency split of the n = 1 mode decreases to 17.993 Hz with a drop of 86.6% compared to the original frequency split. Meanwhile, the vibration velocity at the low-frequency axis increases markedly from 0.433 μm/s to 11.54 μm/s, and the vibration velocity at the high-frequency axis increases from 1.327 μm/s to 5.16 μm/s.  eral rounds of chemical balancing, the lower natural frequency increases to 3291.431 Hz and the higher natural frequency increases to 3309.424 Hz. The frequency split of the n = 1 mode decreases to 17.993 Hz with a drop of 86.6% compared to the original frequency split. Meanwhile, the vibration velocity at the low-frequency axis increases markedly from 0.433 μm/s to 11.54 μm/s, and the vibration velocity at the high-frequency axis increases from 1.327 μm/s to 5.16 μm/s. Experimental results on another same-structure resonator R02 are shown in Figure  13b. The initial natural frequency of resonator R02 is 2872.925 Hz and 2836.133 Hz, respectively, with an initial frequency split of 36.792 Hz. After one round of chemical balancing, the lower natural frequency increases to 2874.658 Hz and the higher natural frequency increases to 2867.261 Hz. The frequency split of the n = 1 mode reduces to 7.397 Hz with a drop of 79.8%. Meanwhile, the velocity of vibration at the high-frequency axis increases from 4.22 μm/s to 9.36 μm/s. The velocity of vibration at the low-frequency axis increases from 6.49 μm/s to 14.88 μm/s. Experimental results on another same-structure resonator R02 are shown in Figure 13b. The initial natural frequency of resonator R02 is 2872.925 Hz and 2836.133 Hz, respectively, with an initial frequency split of 36.792 Hz. After one round of chemical balancing, the lower natural frequency increases to 2874.658 Hz and the higher natural frequency increases to 2867.261 Hz. The frequency split of the n = 1 mode reduces to 7.397 Hz with a drop of 79.8%. Meanwhile, the velocity of vibration at the high-frequency axis increases from 4.22 µm/s to 9.36 µm/s. The velocity of vibration at the low-frequency axis increases from 6.49 µm/s to 14.88 µm/s.
The reason for the dramatic increase in the velocity of vibration could be that, as the second harmonic error decreases, the vibration of the resonator is also less coupled with the base, resulting in less energy dissipation from the base. Meanwhile, the etching rate v remains constant during the balancing experiments. Hence, the mass etched per unit time is the same, which means that the frequency split of the n = 1 mode is proportional to the second harmonic error of the defect mass. The linear relation between the frequency split of n = 1 mode and the etched mass obtained from balancing experiments is consistent with the theoretical calculation.
The natural frequency and frequency split of the two resonators before and after chemical balancing are shown in Figure 13c,d. The natural frequency of the n = 2 mode decreased after chemical balancing. It is noticed that, although we chose a set of balancing parameters for the second harmonic error that theoretically has no impact on the fourth harmonic error, ∆f 2 still has a slight change with the balancing of the second harmonic error. The fluctuation of the ∆f 2 mainly stems from the following two aspects: one is the inevitable error during the measurement of the low-frequency axis of the n = 1 mode. Meanwhile, the control of the inclined angle, immersed depth, and rotating angle may have a slight deviation from the ideal value. These possible sources of error could have affected the fourth harmonic error of the resonator. However, the maximum variation of the ∆f 2 during the balancing of is about 0.024 Hz, which indicates that the change the fourth harmonic error is quite small during the removal of the second harmonic error under the selected balancing parameters where inclined angle α is 1.22 rad and the immersed depth h is 2.25 mm.
The experimental results on resonator R01 and R02 illustrate that the second harmonic error can be reduced from the direction of the low-frequency axis, and the frequency split of the n = 1 mode can be a reference to the residual value of the second harmonic error. The velocity of the vibration is significantly increased which benefits from the decline of the second harmonic error. The performance of the resonator is further improved after balancing on the second harmonic error. The balancing method of the second harmonic error provides an effective and reliable guidance in the determination of the balancing parameters on resonators with different structures.

Conclusions
In this paper, a new method to identify and reduce the second harmonic error is proposed and verified. The results of the theoretical calculation clearly indicate that the frequency split of the n = 1 mode is proportional to the second harmonic error. The simulations have identified the effects of the first three harmonic error on the frequency splits. The second harmonic error makes the major contribution to the frequency split of the n = 1 mode. The second harmonic error can be calculated from the linear relation on the ∆f 1 . The chemical balancing parameters are specially selected to remove the second harmonic error without affecting other harmonics. The second harmonic error of resonator R01 and R02 is independently reduced by 86.6% and 79.8% when the inclined angle α is 1.22 rad and the immersed depth is 2.25 mm. Both the simulation and experimental result confirms the linear relation between the frequency split of n = 1 mode and the second harmonic error. This paper provides an effective and convenient approach to measure and reduce the second harmonic error, by which the performance of the cylindrical resonator gyroscope can be further improved.