一種半球諧振陀螺諧振子動力學建模方法

徐澤遠，伊國興,，魏振楠，趙萬良

1. 哈爾濱工業大學 航天學院，哈爾濱 150001 2. 上海航天控制技術研究所，上海 201109

半球諧振陀螺(HRG)是一種新型的高精度、高穩定性、長壽命的固體振動陀螺[1-4]，非常適合長時間工作場合的使用，在航空航天領域的應用越來越廣泛[5-6]。與機械陀螺、光纖陀螺和激光陀螺等相比[7-9]，它結構簡單，無損耗部件(如機械轉子和光源等)，無需后期維護；不需預熱，啟動時間短；功耗低，體積小，重量輕；具有很強的抗沖擊能力，能承受大的機動過載；抗輻照能力強，失效因素少；諧振子物理特性穩定，陀螺具有很高的可靠性和超長的壽命，連續工作15年的可靠度高達0.99，這些優點使得半球諧振陀螺在慣性技術領域具有廣闊的應用前景，為此發展半球諧振陀螺技術對于中國導航技術的快速發展具有十分重要的意義。

目前，半球諧振陀螺誤差機理分析與抑制技術仍是制約其發展的關鍵因素之一。而建立半球諧振陀螺諧振子的動力學模型為研究其誤差機理問題提供了力學基礎。所以，在半球諧振陀螺的理論研究和實際制造過程中，其誤差機理分析的基礎是含有各種誤差源的半球殼諧振子動力學模型，這些誤差源包括諧振子的密度不均勻、厚度不均勻、品質因數不均勻等加工工藝誤差和溫度、加速度等環境因素的影響[10-16]，還可用于研究電極的加工誤差及其不同控制方式對陀螺精度和性能的影響[17-19]。

對于諧振子的動力學建模一直在不斷地完善，目前的諧振子動力學模型主要對諧振子的唇緣處進行動力學建模，將半球殼諧振子簡化成環形諧振子，解算的動力學參數偏差較大，這是因為只是單純地分析了變形最大的諧振子邊緣處的振動特性，而諧振子是整體振動，必須對整個諧振子的振動特性進行分析。另外，有的研究雖然是對整個諧振子進行建模，但是在推導諧振子中面物理方程過程中，對中面內力作了簡化，直接影響所建立動力學方程的精度；在外載荷加載過程中，將外載荷進行簡化，以致不能準確地反映其實際作用效果，建立的動力學模型精度較低。

針對以上問題，本文在現有的諧振子動力學模型基礎上進行了改進，基于彈性力學薄殼理論，提出了一種完善的諧振子動力學建模方法，結合諧振子的整體動力學分析，推導了更完整的諧振子動力學方程，利用布勃諾夫-伽遼金法建立了更精確的諧振子2階諧振狀態動力學模型。諧振子動力學模型參數的計算值與測試數據結果一致，證明了本文提出的諧振子動力學建模方法的正確性。建立的諧振子動力學模型能更全面地用于半球諧振陀螺誤差機理問題研究。

1 半球諧振陀螺簡介

1.1 半球諧振陀螺工作原理

力反饋式半球諧振陀螺的組成結構主要由外基座、諧振子和內基座三件套構成，如圖1(a)所示。外基座的內表面均勻分布激勵電極，用于諧振子的振幅、速率、正交以及頻率控制。陀螺內基座的外表面均布檢測電極，用于諧振子諧振頻率、振幅、振型角和正交漂移等振動狀態的檢測。下面簡單介紹核心部件諧振子的工作原理。

諧振子如圖1(b)所示，它由高純熔融石英制成，并在其表面噴鍍金屬薄層，用以構成電極的一個電極板。諧振子工作于2階諧振狀態下，其半球殼唇緣處振動駐波為四波腹形式，如圖2(a)所示。振動過程中，唇緣變形達到最大時為橢圓形，橢圓的長軸稱為波腹軸。四波腹狀態下存在2個相互正交的波腹軸。2個橢圓存在4個交點，即為波節點。波腹軸位置用于諧振子的振幅控制和檢測振幅信息，波節點位置用于諧振子的速率控制和檢測振型角信息。

當諧振子敏感軸方向有外界角速度Ω輸入時，諧振子的振型會相對于初始位置發生進動，進動角度?稱為振型角，如圖2(b)所示。理論和實驗都可以證明，在2階諧振狀態下，?與Ω的關系為

(1)

式中：K為諧振子的比例系數。

式(1)說明，只要測量出諧振子的2階振型角，就可以得到諧振子相對于慣性空間轉過的絕對角度。

基于以上分析可知，對諧振子進行動力學分析與精確建模是十分必要的。表1為諧振子的幾何參數(如圖1(b)所示)與物理參數。

表1 諧振子的幾何參數與物理參數Table 1 Geometry and physical parameters of resonator

1.2 諧振子的半球殼模型

在建模過程中，將半徑為R、厚度為h的傘形諧振子結構簡化為半球形進行建模，以半球殼中面的曲率線為θ、φ坐標線，如圖3所示。圖中： 以諧振子球心為原點，定義諧振子坐標系Oxryrzr。

2 諧振子動力學特性分析

2.1 薄殼的彈性力學幾何方程

薄殼的厚度h遠小于殼體中面的最小曲率半徑R，諧振子符合薄殼的條件。在薄殼的彈性力學幾何方程基礎上建立半球殼諧振子動力學模型。為建立符合諧振子振動特性的薄殼彈性力學幾何方程，在基希霍夫-李雅夫假設基礎上作出如下基本假設：

1) 垂直于中面的正應變可以忽略不計。

2) 變形前任何垂直于中面的法線在變形后仍然垂直于中面，而且中面法線及其垂直線段之間的直角保持不變，即該方向的切應變為零。

3) 與中面平行截面上的正應力遠小于其垂直面上的正應力，因而它對應變的影響可以忽略不計。

4) 薄殼上所有加載的面力均可轉化為作用于中面的載荷。

諧振子中面內任意一點P0的位移在θ、φ和γ坐標方向的分量分別用p1、p2和p3表示，沿坐標軸方向的正應變用τ1、τ2和τ3表示，切應變用τ23、τ31和τ12表示。點P0發生位移后到達Q0點，Q0點的坐標為θ+dθ、φ+dφ和γ+dγ。建立中面應變與中面位移關系的幾何方程，在任意一點P0處取一個體積微元，體積微元的所有邊都沿著坐標線θ、φ和γ的方向，如圖4所示。

當坐標改變時，θ線的弧長增量與θ坐標增量的比值稱為θ方向的拉梅系數H1[20]。同理定義φ方向的拉梅系數H2,γ方向的拉梅系數H3。其中γ為直線坐標，直線的拉梅系數為1。

(2)

式中：k1和k2分別為P0點沿θ和φ方向的曲率，且

(3)

其中：x、y和z為直角坐標，對應下文中的半球殼諧振子坐標系。

將應變分量用位移分量來表示，首先計算薄殼中面的正應變，以P0P1的正應變τ1為例。邊P0P3與P1Q2的交角為dη13，邊P0P1在θγ面的曲率半徑為R13,P0P1在θφ面的曲率半徑為R12，且

(4)

(5)

(6)

應用式(4)～式(6)可得正應變為

(7)

其次來考慮切應變，以直角∠P1P0P2的切應變τ12為例。此項切應變是由P0P1和P0P2在θφ面內相向的轉角相加而成，推導過程同上，可得

(8)

由以上推導過程，同理可得彈性力學幾何方程的6個表達式為

(9a)

(9b)

(9c)

(9d)

(9e)

(9f)

建立諧振子動力學方程的過程中以半球殼的中面位移、中面應變、內力以及中面載荷作為討論對象。下面將以彈性力學幾何方程式(9)為基礎推導半球殼諧振子的變形幾何方程、物理方程和平衡微分方程。

2.2 半球殼諧振子變形的幾何方程

p3=w

(10)

令面上各點沿θ和φ方向的位移分別為u和v，即

(11)

根據2.1節中的第2)個基本假設，采用正交曲線坐標系，有τ31=0和τ23=0。根據式(9d)和式(9e)，利用式(2)，以w代替p3，并對γ從0到γ進行積分，注意w不隨γ變化，求解p1和p2，簡化以后，與式(10)聯立可得：

(12)

在半球殼中，γ的最大絕對值是h/2，可見k1γ和k2γ的最大絕對值分別為k1h/2和k2h/2，與1相比，是很小的數值。在文獻[20]中，計算精確到1階小量，1+k1γ和1+k2γ可以用式(13)所示的展開式來代替，達到提高建模精度的目的。

(13)

然后，將式(2)、式(12)、式(13)代入式(9a)、式(9b)和式(9f)，化簡得

(14)

式中：

(15)

式(15)表明中面應變與中面位移之間關系的方程就是半球殼諧振子變形的幾何方程。

2.3 半球殼諧振子的物理方程

在θ面上(在θ為常量的橫截面上)，作用于中面單位寬度上的法向力用T1表示，剪力用T12表示；同樣，在φ面上(在φ為常量的橫截面上)，法向力為T2，剪力為T21；以上4個力稱為半球殼諧振子的中面內力，如圖5(a)所示。在θ面上，作用于單位寬度上的彎矩用M1表示，扭矩用M12表示，橫向剪力用S1表示；在φ面上，彎矩用M2表示，扭矩用M21表示，橫向剪力用S2表示；以上6個力稱為彎曲內力，如圖5(b)所示。

通過分析計算得出中面內力及彎曲內力為

(16)

2.4 半球殼諧振子的平衡微分方程

下面建立表述半球殼諧振子中面的內力與所受載荷之間關系的平衡微分方程。考慮任意微分殼體OO1O3O2的平衡，如圖6所示。在圖6中，把中面內力和橫向剪力畫在一個圖中，如圖6(a)所示；而把彎矩和扭矩畫在一個圖中，如圖6(b)所示，圖中的X、Y、Z是單位中面面積范圍內的載荷。

首先分析了各力在OO1θ軸上的投影，從而建立平衡微分方程∑(T+S)θ=0，將相關的投影分量相加，令總和等于零，再進行微分，得到平衡微分方程∑(T+S)θ=0的投影方程。同理可得∑(T+S)φ=0和∑(T+S)γ=0的平衡微分方程。將所有的力對OO1θ、OO2φ和Oγ求矩，得到∑Mθ=0、∑Mφ=0和∑Mγ=0的平衡微分方程，建立的平衡微分方程為

(17a)

(17b)

ABZ=0

(17c)

(17d)

(17e)

(17f)

3 諧振子動力學方程

3.1 諧振子的半球殼模型

1) 諧振子坐標系

以諧振子球心為原點，定義諧振子坐標系Oxryrzr，如圖7所示。諧振子上點Q的矢徑q可以表示為

q=xxr+yyr+zzr

(18)

式中：xr、yr和zr為諧振子坐標系的單位矢量；矢徑q的坐標值x、y和z表達式為

(19)

2) 諧振子中面一點沿θ和φ方向的拉梅系數為

(20)

3) 任意一點沿θ和φ方向的曲率半徑為

R1=R2=R

(21)

4) 當諧振子中面發生形變時，Q(θ,φ,R)變形之后到達Q′，產生的位移矢量W在局部坐標系(t1t2n系)中的表達式為

W=ut1+vt2+wn

(22)

5)xryrzr系與t1t2n系的轉換關系為

(23)

3.2 諧振子動力學方程

由于求解式(15)～式(17)這3個偏微分方程組的難度很大，故而采用諧振子中面不可拉伸的假設[23]，諧振子中面上的正應變和切應變等于零，滿足

ε1=ε2=ε12=0

(24)

為了分析半球諧振陀螺諧振子的振動特性，必須要分析其位移分量，因此以位移分量表示諧振子動力學方程。將諧振子的半球殼模型參數式(18)～式(23)及中面不可拉伸的假設式(24)代入式(15)～式(17)后化簡，并將式(15)～式(17)聯立可得諧振子動力學方程為

(25)

由式(25)可知，該式能用于建立含有中面半徑R不均勻、厚度h不均勻的動力學方程。式(25)是在現有動力學方程基礎上進行改進建立的更精確的諧振子動力學方程，體現在以下4個方面：

1) 更全面的動力學分析

在建立物理方程時考慮了曲率改變量和扭率改變量，對物理方程中內力的分析更加完善，有利于建立更精確的動力學方程。

在建立諧振子動力學方程方面，現有的動力學方程在諧振子唇緣處進行動力學建模，將半球殼諧振子簡化成環形諧振子，不能全面地描述諧振子的振動特性。為解決這一簡化問題，動力學方程式(25)描述了整個半球殼的振動特性，可以得到半球殼任意θ∈[0,π/2]的動力學方程。

2) 更完整的半球殼諧振子動力學建模

在加載外載荷方面，動力學方程式(25)沒有將外載荷近似簡化加載在諧振子唇緣處，而是在建模過程中，保留了外載荷的實際作用效果，有利于反映諧振子的實際振動特性，可以加載求解得到更精確的諧振子2階諧振狀態動力學模型。

3) 有利于建立更準確的激勵電極模型

因為激勵電極是分布在θ∈[75°,85°]、φ∈[0°,360°]上的一個獨立區間，而不是諧振子邊緣上的一個點，進行激振力的加載時是作用在整個激勵電極上的，而不是作用在諧振子的邊緣上，可以更準確地進行激振力的加載，對于后續控制系統激勵電極的建模很有意義。

4) 更方便求解外載荷的影響

動力學方程式(25)是線性偏微分方程組，外載荷分析滿足疊加原理，可對不同形式的外載荷進行單獨分析。

4 諧振子2階諧振狀態動力學模型

4.1 布勃諾夫-伽遼金法求解

諧振子處于工作狀態時，不考慮外界加速度和溫度的影響，作用在諧振子上的慣性載荷為

(26)

考慮加速度和溫度等外界環境因素的影響時，對諧振子的影響主要通過式(26)體現，由此可以分析環境因素通過影響諧振子從而導致的陀螺漂移機理。

諧振子在2階諧振狀態時，環向波數為n的變形方程滿足瑞利里茨函數為

(27)

式中：振型方程W2=-[(2+cosθ)tan2θ]/2,U2=V2=(sinθtan2θ)/2；ω2為2階諧振角頻率；t為時間；p(t)和q(t)為與時間相關的系數。

由文獻[23]可知，當諧振子進入2階穩態振動時，p(t)和q(t)近似滿足

(28)

式中：a和b為與諧振子振幅相關的系數。

本文利用布勃諾夫-伽遼金法求解近似解析解，過程為

L=ut1+vt2+wn=[(U2cos 2φ)t1+(V2sin 2φ)

t2+(W2cos 2φ)n]p(t)+[(U2sin 2φ)t1-

(V2cos 2φ)t2+(W2sin 2φ)n]q(t)=

L1(x)p(t)+L2(x)q(t)

(29)

將式(26)～式(28)代入動力學方程式(25)中，并用Ax、Ay和Az分別表示動力學方程組的3個式子，則半球殼諧振子邊緣的運動方程M可表示為

M=Axt1+Ayt2+Azn=0

(30)

由布勃諾夫-伽遼金法方法得：

(31)

式中：Vhs為半球殼區域，化為球坐標，可得半球殼諧振子的2階固有振型動力學方程為

(32)

求解結果為

(33)

(34)

其中：

(35)

(36)

(37)

(38)

建模過程中必須考慮到實際諧振子中存在材料內摩擦，模型應該是一個具有衰減環節的2階系統，因此在模型中加入阻尼比ξ2，將方程寫成標準形式為

(39)

式中：ξ2=1/(Qfω2)。

考慮諧振子的材料內摩擦，由式(39)可得帶有品質因數的p(t)和q(t)表達式為

(40)

由式(32)、式(39)和式(40)可知，能建立含有密度ρ不均勻、品質因數Qf不均勻的動力學方程，以此為基礎分析加工工藝誤差對陀螺的影響。

4.2 測試驗證

將表1中參數代入文獻[21]中的動力學模型，求解出的比例系數和2階諧振頻率為K≈0.313,f2≈2 578.0 Hz。將表1中參數代入文獻[22]中得到比例系數為K≈0.277。文獻[23]中動力學模型求解得到的比例系數和2階諧振頻率為K≈0.277,f2≈10 269 Hz。

與文獻[21-23]中的比例系數和2階諧振頻率的結果相比，本文比例系數、2階諧振頻率的解算結果與實際測試結果基本一致。

5 結 論

針對半球諧振陀螺諧振子動力學建模的實際應用問題，基于彈性力學薄殼理論，提出了完善的諧振子動力學建模方法，為其他振動類陀螺的動力學建模提供了理論基礎。

1) 在彈性力學幾何方程的基礎上，推導了諧振子的變形幾何方程、物理方程以及平衡微分方程，建立了精確的諧振子動力學方程，給出了動力學方程建立過程中相關問題的處理方法。

2) 與現有的動力學方程相比，建立的諧振子動力學方程更準確、更有優勢，主要體現在：① 更全面的動力學分析；② 更完整的半球殼諧振子動力學建模；③ 有利于建立更準確的激勵電極模型；④ 更方便求解外載荷的影響。

3) 通過半球諧振陀螺諧振子的參數計算和實際測試結果的對比驗證，結果一致，證明建立了準確的諧振子動力學模型。

4) 建立的諧振子動力學方程為以后研究諧振子的中面半徑不均勻、密度不均勻、厚度不均勻、品質因數不均勻等加工工藝誤差以及加速度、溫度等環境因素導致的陀螺誤差機理問題提供了動力學基礎。

參 考 文 獻

[1] ROZELLE D M. The hemispherical resonator gyro: From wineglass to the planets[C]∥Proceeding 19th AAS/AIAA Space Flight Mechanics Meeting. Reston, VA: AIAA, 2009: 1157-1178.

[2] XU Z Y, YI G X, QI Z Y, et al. Structural optimization research on hemispherical resonator gyro based on finite element analysis[C]∥The 35th Chinese Control Conference. Beijing: Technical Committee on Control Theory, Chinese Association of Automation, 2016: 5737-5742.

[3] MATTHEWS A, RYBAK F J. Comparison of hemispherical resonator gyro and optical gyros[J]. IEEE Aerospace and Electronic Systems Magazine, 1992, 7(5): 40-46.

[4] MEYER A D, ROZELLE D M. Milli-HRG inertial navigation system[C]∥The IEEE/ION Position, Location and Navigation Symposium (PLANS’12). Piscataway, NJ: IEEE Press, 2012: 24-29.

[5] ZHBANOV Y K. Vibration of a hemispherical shell gyro excited by an electrostatic field[J]. Mechanics of Solids, 2008, 43(3): 328-332.

[6] ZHBANOV Y K. Self-tuning control loop for suppression of quadrature in a hemispherical resonator gyro[J]. Gyroscopy and Navigation, 2007(2): 37-42.

[7] 王錦瑜, 馮培德. 激光陀螺速率偏頗系統的分析研究[J]. 航空學報, 2001, 22(1): 46-50.

WANG J Y, FENG P D. Research on rate-bias system of laser gyro[J]. Acta Aeronautica et Astronautica Sinica, 2001, 22(1): 46-50 (in Chinese).

[8] 李金濤, 房建成. 高精度光纖IMU的磁屏蔽方法及實驗研究[J]. 航空學報, 2011, 32(11): 2106-2116.

LI J T, FANG J C. Magnetic shielding method and experiment study of inertial measurement unit based on high precision fiber-optic gyroscope[J]. Acta Aeronautica et Astronautica Sinica, 2011, 32(11): 2106-2116 (in Chinese).

[9] 魏玉淼, 董永貴, 李昊. 微機械陀螺非線性特性的自由衰減振蕩測量方法[J]. 儀器儀表學報, 2016, 37(11): 2465-2472.

WEI Y M, DONG Y G, LI H. Free damped oscillation measurement method for the nonlinear features of micromechanical gyroscopes[J]. Chinese Journal of Scientific Instrument, 2016, 37(11): 2465-2472 (in Chinese).

[10] 李巍, 任順清, 王常虹. 半球諧振陀螺諧振子品質因數不均勻引起的誤差分析[J]. 航空學報, 2013, 34(1): 121-129.

LI W, REN S Q, WANG C H. Analysis for impact of resonator’sQ-factor nonuniformity on the error of hemispherical resonator gyro[J]. Acta Aeronautica et Astronautica Sinica, 2013, 34(1): 121-129 (in Chinese).

[11] SHATALOV M Y, JOUBERT S V, COETZEE C E. The influence of mass imperfections on the evolution of standing waves in slowly rotating spherical bodies[J]. Journal of Sound and Vibration, 2009, 330(1): 127-135.

[12] ZHBANOV Y K. Amplitude control contour in a hemispherical resonator gyro with automatic compensation for difference inQ-factors[J]. Mechanics of Solids, 2008, 43(3): 328-332.

[13] SHATALOV M Y, COETZEE C E. Dynamics of rotating and vibrating thin hemispherical shell with mass and damping imperfections and parametrically driven by discrete electrodes[J]. Gyroscopy and Navigation, 2011, 2(1): 27-33.

[14] 伊國興, 謝陽光, 王常虹, 等. 加速度對半球諧振陀螺振動檢測系統影響分析[J]. 中國慣性技術學報, 2013, 21(5): 676-681.

YI G X, XIE Y G, WANG C H, et al. Analysis of acceleration influence on HRG vibration detection system[J]. Journal of Chinese Inertial Technology, 2013, 21(5): 676-681 (in Chinese).

[15] FREITAG S, BEER M, GRAF W, et al. Lifetime prediction using accelerated test data and neural networks[J]. Computers & Structures, 2009, 87(19-20): 1187-1194.

[16] WANG X, WU W Q, FANG Z, et al. Temperature drift compensation for hemispherical resonator gyro based on natural frequency[J]. Sensors, 2012, 12(5): 6434-6446.

[17] SONG J W, SONG H M, LEE Y J, et al. Design of oscillation control loop with coarse-precision mode transition for solid-state resonant gyroscope[J]. IEEE Sensors Journal, 2016, 16(6): 1730-1742.

[18] WANG X, WU W Q, LUO B, et al. Force to rebalance control of HRG and suppression of its errors on the basis of FPGA[J]. Sensors, 2011, 11(12): 11761-11773.

[19] ZHURAVLEV V P. Hemispherical resonator gyro withmdata electrodes andncontrol electrodes[J]. Mechanics, 2015, 50(4): 375-378.

[20] 徐芝綸. 彈性力學下冊[M]. 5版. 北京: 高等教育出版社, 2016: 100-200.

XU Z L. Elasticity[M]. 5th ed. Beijing: Higher Education Press, 2016: 100-200 (in Chinese).

[21] 趙洪波, 任順清, 李巍. 半球諧振子動力學方程的建立及固有頻率的計算[J]. 哈爾濱工業大學學報, 2010, 42(11): 1702-1706.

ZHAO H B, REN S Q, LI W. Establishment of dynamics equation of HRG resonator and calculation of natural frequency[J]. Journal of Harbin Institute of Technology, 2010, 42(11): 1702-1706 (in Chinese).

[22] 謝陽光. 基于半球諧振陀螺的姿態測量系統研究[D]. 哈爾濱: 哈爾濱工業大學, 2013: 13-36.

XIE Y G. Research on attitude measurement systems based on hemispherical resonator gyro[D]. Harbin: Harbin Institute of Technology, 2013: 13-36 (in Chinese).

[23] BA 馬特維耶夫, H И 利帕特尼科夫, AB 阿廖欣, 等. 固體波動陀螺[M]. 北京: 國防工業出版社, 2009: 1-43.

MATVEEV V A, BASARAB M A, ALEKIN A V, et al. Solid state wave gyro[M]. Beijing: National Defense Industry Press, 2009: 1-43 (in Chinese).