光纖陀螺儀(FOG)中光纖尾纖的裝配主要依靠工人手工完成，裝配應(yīng)力控制憑借人工經(jīng)驗，一致性難以保障。針對光纖裝配應(yīng)力定量控制的難題，提出一種基于離散微分幾何理論的光纖動力學(xué)建模仿真方法。首先，基于離散微分幾何理論建立了光纖運動學(xué)框架，并由彈性勢能推導(dǎo)了包含拉伸、彎曲和扭轉(zhuǎn)變形的光纖動力學(xué)方程以及應(yīng)力計算方法。其次，提出了基于Newmark隱式算法的光纖動力學(xué)方程數(shù)值求解方法，并設(shè)計了光纖在裝配過程中與周圍物體的碰撞檢測與響應(yīng)策略。隨后，進行了光纖圓周布設(shè)和彎扭變形的仿真測試。結(jié)果表明，光纖控制端與固定端之間的自由長度越短，光纖圓周布設(shè)的位置越精確;光纖截面半徑越大，在同樣的扭轉(zhuǎn)角下光纖產(chǎn)生的彎扭變形程度越大。最后，設(shè)計了一種光纖陀螺儀裝配布局，并通過動力學(xué)仿真得到了光纖裝配過程的應(yīng)力分布及最大應(yīng)力變化，為陀螺儀光路低應(yīng)力裝配與路徑規(guī)劃提供了理論支撐。

　　關(guān)鍵詞

　　光纖陀螺儀;光纖尾纖;動力學(xué)模型;離散微分幾何;應(yīng)力計算

　　論文《光纖陀螺儀光纖尾纖裝配動力學(xué)建模仿真與應(yīng)力計算方法》發(fā)表在《計算機集成制造系統(tǒng)》，版權(quán)歸《計算機集成制造系統(tǒng)》所有。本文來自網(wǎng)絡(luò)平臺，僅供參考。

　　0 引言

　　光纖陀螺儀是一種基于薩格納克效應(yīng)的角速度傳感器，是慣性導(dǎo)航系統(tǒng)的核心部件之一，它以全固態(tài)、高精度、高穩(wěn)定和低成本等優(yōu)點在航空航天、船舶及軍事等領(lǐng)域得到廣泛應(yīng)用[1]。光路裝配是光纖陀螺儀裝配的核心環(huán)節(jié)，光路裝配對象包括5大光學(xué)器件、光纖、熔接點等，光路裝配質(zhì)量的好壞直接關(guān)系到光纖陀螺儀的功能與輸出性能參數(shù)。其中光纖尾纖裝配是光路裝配的重要組成步驟，光纖尾纖的裝配質(zhì)量對于中高精度光纖陀螺儀最終性能的影響不容忽視。隨著中高精度光纖陀螺儀應(yīng)用范圍越來越廣泛，對光纖尾纖裝配質(zhì)量和一致性控制的研究顯得日益重要。

　　光纖尾纖裝配質(zhì)量對于光纖陀螺儀精度和性能的影響主要體現(xiàn)在光纖彎曲應(yīng)力、扭曲應(yīng)力造成的損耗及光波偏振態(tài)變化上[2-3]。在目前光纖陀螺儀光路的應(yīng)力研究中，學(xué)者們大多是針對光纖環(huán)中的光纖，通過應(yīng)力分析儀來進行光纖應(yīng)力的分析和測試[4-6]，較少涉及光路裝配過程中光纖尾纖應(yīng)力狀態(tài)的評估。而光纖尾纖在光纖陀螺儀中所受到的扭曲應(yīng)力和彎曲應(yīng)力與光路裝配路徑有很大關(guān)系。

　　另一方面，光纖陀螺儀中光纖本身由纖芯和包層組成，涉及到復(fù)雜的材料和力學(xué)特性，目前研究者主要通過將光纖簡化為彈性變形進行建模分析[7-8]。而近年來，在彈性線體的變形預(yù)測中，通過建立包含重力、拉壓、彎曲和扭轉(zhuǎn)等特性的物理模型來進行動態(tài)仿真已成為主要研究方法[9-11]。Cosserat彈性桿模型[12-13]包含了彈性細(xì)桿的軸向線應(yīng)變及彎曲剪應(yīng)變等物理參數(shù)，具有更精確的仿真計算，但涉及到求解非線性偏微分方程難以獲得精確解、離散化后參數(shù)過多等問題。考慮到模型仿真的真實性和實時性，離散微分幾何理論被用于處理連續(xù)線體的離散化建模中[10,14]，通過保留表征實際(光滑)物理系統(tǒng)的關(guān)鍵幾何結(jié)構(gòu)來對連續(xù)線體進行離散化，在保證連續(xù)系統(tǒng)的魯棒性的同時極大地簡化了計算量[15]。其中，BERGOU等[16]通過一種離散處理的自適應(yīng)框架曲線，建立了彈性細(xì)桿的離散幾何模型，并模擬了彈性細(xì)桿的彎曲、扭轉(zhuǎn)以及打結(jié)現(xiàn)象。LV等[17]利用離散微分幾何理論建立了柔性線纜的離散彈性細(xì)桿動力學(xué)模型，并采用中心線-角度的表示方法，提出了離散彈性細(xì)桿模型的動力學(xué)求解方法，從而對柔性線纜進行實時動態(tài)變形模擬。目前，離散微分幾何理論已經(jīng)廣泛地應(yīng)用于具有復(fù)雜的線體或類線體的場景建模中，特別是在計算機圖形學(xué)中[18-20]。

　　因此，本文提出一種基于離散微分幾何理論的光纖動力學(xué)建模仿真方法，研究光纖陀螺儀裝配過程中光纖變形與應(yīng)力變化，用于建立光纖尾纖應(yīng)力與布設(shè)路徑之間的關(guān)系，以規(guī)劃最佳光纖布設(shè)路徑。該方法采用自適應(yīng)參考框架與中心線結(jié)合的方式，描述光纖裝配過程中包括拉伸、彎曲和扭轉(zhuǎn)變形的運動學(xué)位姿，并通過動力學(xué)模型快速計算光纖裝配過程中的應(yīng)力變化。

　　1 光纖動力學(xué)模型

　　1.1 運動學(xué)框架

　　在光纖陀螺儀中，光纖尾纖的幾何表示是建立光纖動力學(xué)模型的關(guān)鍵問題。本文假設(shè)光纖由一條中心線和一系列均勻各向同性的圓形截面組成，其中，中心線穿過每個圓形截面的幾何中心，中心線是一條空間光滑曲線且至少具有二階可微性，中心線上的截面可以相對移動和旋轉(zhuǎn)。

　　在全局坐標(biāo)系中(Oξηζ)，光纖的中心線C可以用向量r(s, t)表示，其中s是以中心線的一端為原點建立的弧坐標(biāo)系，t為時間。則中心線的運動學(xué)可由其運動速度$dot{r}(s, t)=partial r / partial t$與變形梯度即中心線的切向量$r'(s, t)=partial r / partial s$表示。為描述光纖的彎扭狀態(tài)，定義一個與剛性橫截面固連的材料框架{d?(s, t), d?(s, t), d?(s, t)}，3個單位向量d?(k=1,2,3)兩兩互相垂直，且d?=d?×d?是橫截面的單位法向量。

　　對于細(xì)長光纖，通常忽略其剪切變形，這時光纖橫截面始終與中心線正交，即橫截面單位法向量d?與中心線切線方向的單位矢量t重合：

　　[d_{3}=t=frac{r'}{left|r' ight|}]

　　則光纖的拉伸應(yīng)變ε?可表示為：

　　[varepsilon_{r}=r' cdot d_{3}=left|r' ight|]

　　而根據(jù)無窮小旋轉(zhuǎn)的定義，材料框架{d?, d?, d?}在空間上沿弧長的變化d?'(k=1,2,3)和在時間上的變化$dottl7lrfpl_{k}$(k=1,2,3)為：

　　[d_{k}'=omega × d_{k}, dottl7lrfpl_{k}=Omega × d_{k}]

　　其中：ω表示光纖的彎扭應(yīng)變即彎扭度矢量，Ω表示光纖的彎扭速度，即角速度矢量：

　　[omega=d_{3} × d_{3}'+left(d_{1}' cdot d_{2} ight) d_{3}]

　　[Omega=d_{3} × dottl7lrfpl_{3}+left(dottl7lrfpl_{1} cdot d_{2} ight) d_{3}]

　　在不考慮剪切應(yīng)變的情況下，彎扭度矢量ω與角速度矢量Ω可表示為：

　　[omega=t × t'+left(d_{1}' cdot d_{2} ight) t=K b+omega_{3} t]

　　[Omega=t × dot{t}+left(dottl7lrfpl_{1} cdot d_{2} ight) t=t × dot{t}+Omega_{3} t]

　　其中：Kb=t×t'為中心線副法線曲率，K=|t'|為曲率大小，b=t×t'/|t'|為副法線方向的單位矢量，ω?=d?'·d?為扭轉(zhuǎn)應(yīng)變即ε?，Ω?=$dottl7lrfpl_{1}·d_{2}$為扭轉(zhuǎn)角速度。則光纖的彎扭應(yīng)變ε?、ε?可表示為：

　　[varepsilon_{b}=|omega-(omega cdot t) t|=K, varepsilon_{t}=|omega cdot t|=omega_{3}]

　　為了表示光纖材料框架的位置和方向，需要光纖的中心線坐標(biāo)和材料框架的扭轉(zhuǎn)度。對于材料框架的表示，有四元數(shù)[21-22]和自適應(yīng)框架[14,16]兩種典型的方案。中心線與四元數(shù)組合的方案需要7個變量，其中包含冗余變量和多余約束，而中心線與自適應(yīng)框架組合的方案只需要4個變量。因此，采用中心線和自適應(yīng)框架組合的方案來表示材料框架的位置和方向。中心線的自適應(yīng)框架{h?, h?, h?}通常可以用Bishop框架即空間平行框架B={b?, b?, t}和時間平行框架A={a?, a?, t}表示[14]，其中空間平行框架的彎扭度矢量ω_B和時間平行框架的角速度矢量Ω_A分別表示為：

　　[omega_{B}=K b, Omega_{A}=t × dot{t}]

　　進一步將光纖離散為N+2個節(jié)點(0,1,?,N+1)和N+1個圓柱段(0,1,?,N)。初始時刻離散節(jié)點均勻分布，節(jié)點坐標(biāo)為x?={x?, y?, z?}，光纖第i條線段的切向量為s?=x???−x?，對應(yīng)單位切向量為t?=s?/|s?|，每個光纖線段上材料框架{d??, d??, d??}與時間平行框架A={a??, a??, t?}之間的夾角為θ。離散化后光纖位姿的廣義坐標(biāo)系表示為X，大小為4N+7：

　　[X=left{x_{0}, heta^{0}, x_{1}, heta^{1}, cdots, x_{N}, heta^{N}, x_{N+1} ight}^{T}]

　　離散化后光纖第i條線段的拉伸應(yīng)變ε??可表示為：

　　[varepsilon_{e}^{i}=frac{left|s^{i} ight|}{left|ar{s}^{i} ight|}]

　　其中|s?|和|$ar{s}$?|分別為第i條線段的當(dāng)前長度和初始長度。

　　第i個節(jié)點的彎扭應(yīng)變ε??、ε??可表示為：

　　[varepsilon_{i}^{b}=K_{i}=left|t_{i}' ight|=frac{kappa_{i}}{ar{l}_{i}}, varepsilon_{i}^{t}=frac{alpha_{i}}{ar{l}_{i}}]

　　其中：$ar{l}_{i}$=(|$ar{s}$??¹|+|$ar{s}$?|)/2為第i個節(jié)點的Voronoi區(qū)域長度，即第i-1條線段和第i條線段長度和的一半;α?為第i個節(jié)點的離散扭轉(zhuǎn)角;K?為第i個節(jié)點的副法線曲率;κ?為第i個節(jié)點的離散積分曲率，定義為：

　　[kappa_{i}=2 tan frac{varphi_{i}}{2},(kappa b)_{i}=frac{2 t^{i-1} × t^{i}}{1+t^{i-1} cdot t^{i}}]

　　其中φ?為相鄰線段i-1和i之間的夾角。

　　為了計算離散扭轉(zhuǎn)角α?，引入了平行傳遞的概念，即向量c?到c?的最小旋轉(zhuǎn)為：

　　[P_{c_{1}}^{c_{2}}=Rleft(c_{1} × c_{2}, angleleft(c_{1}, c_{2} ight) ight)]

　　具體而言，平行傳遞指沿著c?×c?的方向使向量c?旋轉(zhuǎn)角度∠(c?, c?)。空間平行框架只取決于中心線位置，這也意味著空間平行框架之間的平行傳遞為：

　　[b_{k}^{i}left(t_{n} ight)=P_{t^{i}}^{t^{i-1}}left(t_{n} ight) b_{k}^{i-1}left(t_{n} ight), k=1,2]

　　將材料框架{d?, d?, d?}和空間平行框架{b?, b?, t}從第i-1條線段平行傳遞到第i條線段得到離散扭轉(zhuǎn)角α?=θ_B?−θ_B??¹，其中θ_B??¹和θ_B?分別為相鄰的第i-1條線段和第i條線段上的材料框架與空間平行框架之間的夾角。類似地，同一線段上的時間平行框架也具有平行傳遞關(guān)系：

　　[a_{k}^{i}left(t_{n} ight)=P_{tleft(t_{n-1} ight)}^{tleft(t_{n} ight)} a_{k}^{i}left(t_{n-1} ight), k=1,2]

　　將時間平行框架的向量a???¹進行空間上的平行傳遞得到向量P???¹a???¹，使得向量P???¹a???¹和b??夾角為β??¹，即與第i-1條線段上的空間平行框架和時間平行框架之間的夾角保持一致。最后第i個節(jié)點上的離散扭角計算為：

　　[alpha_{i}= heta_{B}^{i}- heta_{B}^{i-1}= heta^{i}- heta^{i-1}+alpha_{r j}^{i}]

　　其中α???為向量P???¹a???¹和a??之間的夾角。

　　1.2 動力學(xué)模型

　　光纖陀螺儀中的光纖在裝配過程中變形大但應(yīng)變小，可近似視作純彈性線。基于1.1節(jié)中的光纖運動學(xué)框架，在不考慮剪切變形時，光纖的彈性勢能E?包括拉壓彈性勢能E?(d)、彎曲彈性勢能E?(d)和扭轉(zhuǎn)彈性勢能E?(d)，離散化形式為：

　　[

　　left{

　　egin{array}{c}

　　E_{e(d)}=frac{1}{2} sum_{i=0}^{N} frac{E_{e} A^{i}}{left|ar{s}^{i} ight|}left(left|s^{i} ight|-left|ar{s}^{i} ight| ight)^{2}

　　E_{b(d)}=frac{1}{2} sum_{i=1}^{N} frac{E_{b} I_{i}}{ar{l}_{i}}left(kappa_{i}-ar{kappa}_{i} ight)^{2}

　　E_{t(d)}=frac{1}{2} sum_{i=1}^{N} frac{G I_{p, i}}{ar{l}_{i}}left(alpha_{i}-ar{alpha}_{i} ight)^{2}

　　end{array}

　　 ight.

　　]

　　其中：E?、E?和G分別為光纖的拉伸模量、彎曲模量和剪切模量(本研究中E?=E?=E，E為楊氏模量，G=E/2(1+ν)，ν為泊松比);A?為第i條線段上的橫截面積;I?和I?,?分別為第i個節(jié)點上的截面慣性矩和極慣性矩;$ar{varepsilon}_{i}^{b}$、$ar{varepsilon}_{i}^{t}$、$ar{kappa}_{i}$和$ar{alpha}_{i}$分別為第i個節(jié)點上Voronoi區(qū)域在初始時刻的彎曲應(yīng)變、扭轉(zhuǎn)應(yīng)變、離散積分曲率和離散扭轉(zhuǎn)角。

　　節(jié)點力是總彈性勢能對位置矢量的偏導(dǎo)[23]，而線段上的扭矩是彈性勢能對角度的偏導(dǎo)。第i個節(jié)點上的彈性力F??(F??=F??+F??+F??)包括拉伸力F??、彎曲力F??和扭轉(zhuǎn)力F??，而扭矩M?可以由扭轉(zhuǎn)彈性勢能E?(d)對角度θ?的偏導(dǎo)得到。

　　因此，光纖在第i個節(jié)點受到的軸向應(yīng)力σ??和切向應(yīng)力σ??可由節(jié)點彈性力F??得到：

　　[sigma_{i}^{n}=frac{left(F_{i}^tl7lrfpl ight)^{n}}{A}, sigma_{i}^{t}=frac{left(F_{i}^tl7lrfpl ight)^{t}}{A}]

　　其中(F??)?和(F??)?分別為第i個節(jié)點的軸向彈性力與切向彈性力，定義第i個節(jié)點的軸向矢量(即沿中心線的切向矢量)為相鄰線段矢量的角平分線矢量t??=t??¹+t?，節(jié)點處的截面面積A=πr²。

　　廣義彈性力F?的列向量可表示為：

　　[F_{e}=left{cdots, F_{i, x}^tl7lrfpl, F_{i, y}^tl7lrfpl, F_{i, z}^tl7lrfpl, M_{i}, cdots ight}^{T}]

　　此外，光纖還受到廣義外力F???的作用(例如重力和來自周圍環(huán)境的接觸力及接觸力矩等)，以及廣義速度阻尼力F?。根據(jù)牛頓定律，總的動力學(xué)方程為：

　　[M cdot a=F_{e}+F_{n}+F_{e x t}]

　　其中：M為節(jié)點質(zhì)量和轉(zhuǎn)動慣量組成的廣義質(zhì)量矩陣，a為廣義加速度的列向量。

　　2 動力學(xué)算法及碰撞策略

　　2.1 動力學(xué)算法

　　為保證動力學(xué)模型求解的高效性與良好的魯棒性，本文采用Newmark多步時間積分隱式算法[24]求解動力學(xué)系統(tǒng)。首先，將動力學(xué)方程改寫為二階矩陣形式：

　　[M ddot{U}left(t_{n+1} ight)+D dot{U}left(t_{n+1} ight)+K Uleft(t_{n+1} ight)=F_{e x t}left(t_{n} ight)]

　　其中：U(t???)、$dot{U}(t_{n+1})$和$ddot{U}(t_{n+1})$分別為下一時刻t???的未知廣義位移、廣義速度和廣義加速度;F???(t?)為當(dāng)前時刻t?的已知廣義外力;廣義位移U包括節(jié)點位移和線段扭轉(zhuǎn)角：

　　[U=left{u_{1}, eta^{1}, u_{2}, eta^{2}, cdots, u_{N}, eta^{N}, u_{N+1} ight}]

　　其中u?={u?, v?, w?}為第i個節(jié)點的位移，β?為第i條線段的扭轉(zhuǎn)角。因此可以得到廣義位移U與廣義坐標(biāo)X之間的關(guān)系。D為速度阻尼剛度矩陣(即對角線上元素為阻尼系數(shù)的對角矩陣)，K為彈性剛度矩陣，可由廣義彈性力F?線性化得到：

　　[K=frac{partial F_{e}}{partial U}, F_{e}=K U+b_{e}]

　　其中b?(t?)=F?(t?)−K U(t?)為常數(shù)項，最后需要移動到F???(t?)中。

　　在確定好動力學(xué)方程的系數(shù)矩陣后，進行Newmark算法。首先，計算t???時刻的廣義位移和廣義速度的預(yù)測值：

　　[ ilde{U}left(t_{n+1} ight)=Uleft(t_{n} ight)+dot{U}left(t_{n} ight) Delta t_{n+1}+ddot{U}left(t_{n} ight)left(frac{1}{2}-lambda_{2} ight) Delta t_{n+1}^{2}]

　　[ ilde{dot{U}}left(t_{n+1} ight)=dot{U}left(t_{n} ight)+ddot{U}left(t_{n} ight)left(1-lambda_{1} ight) Delta t_{n+1}]

　　其中取權(quán)重系數(shù)λ?=1/2和λ?=1/4，得到梯形規(guī)則下時間步長的無條件穩(wěn)定性。最后得到t???時刻的廣義加速度$ddot{U}(t_{n+1})$、廣義速度$dot{U}(t_{n+1})$、廣義位移U(t???)為：

　　[

　　egin{gathered}

　　ddot{U}left(t_{n+1} ight)=left(M+D lambda_{1} Delta t_{n+1}+K lambda_{2} Delta t_{n+1}^{2} ight)^{-1} cdot

　　left(F_{e x t}left(t_{n} ight)-D ilde{dot{U}}left(t_{n+1} ight)-K ilde{U}left(t_{n+1} ight) ight)

　　end{gathered}

　　]

　　[dot{U}left(t_{n+1} ight)= ilde{dot{U}}left(t_{n+1} ight)+ddot{U}left(t_{n+1} ight) lambda_{1} Delta t_{n+1}]

　　[Uleft(t_{n+1} ight)= ilde{U}left(t_{n+1} ight)+ddot{U}left(t_{n+1} ight) lambda_{2} Delta t_{n+1}^{2}]

　　2.2 碰撞檢測與響應(yīng)

　　在光纖陀螺儀的光纖裝配過程中存在著光纖與其他物體之間的接觸碰撞，為了保證動力學(xué)模型仿真的真實性，需要進行光纖的碰撞檢測與響應(yīng)策略研究。

　　在本研究中，采用楊嘯東等[25]提出的分層碰撞檢測方法，即第一層為包圍光纖整體的最小圓柱體包圍盒，第二層為包圍光纖線段的最小球體包圍盒，第三層為包圍光纖線段的最小圓柱體包圍盒。當(dāng)檢測到第一層發(fā)生碰撞時，執(zhí)行第二層碰撞檢測;當(dāng)檢測到第二層發(fā)生碰撞時，執(zhí)行第三層碰撞檢測。此外，對于相鄰的光纖線段認(rèn)為不存在碰撞情況。

　　對于光纖與周圍物體之間的碰撞響應(yīng)(如光纖裝配過程與陀螺儀結(jié)構(gòu)面的接觸碰撞)，采用一種接觸力的響應(yīng)策略。當(dāng)光纖第i個節(jié)點發(fā)生碰撞時，節(jié)點會受到碰撞物體的接觸力F?，包括法向支撐力F??和切向摩擦力F??：

　　[F_{i}=F_{i}^{n}+F_{i}^{f}]

　　其中支撐力F??使得節(jié)點沿接觸面法向的受力和速度變?yōu)榱?假設(shè)碰撞過程為完全非彈性碰撞)，則碰撞后的法向支撐力F??和節(jié)點切向速度v????為：

　　[F_{i}^{n}=-left(F_{i}^{g}+F_{i}^tl7lrfpl ight)^{n}]

　　[v_{i}^{t a n}=v_{i}-v_{i}^{n}]

　　其中：F??為節(jié)點重力，F(xiàn)??為節(jié)點彈性力，v?為節(jié)點碰撞前的速度，v??為節(jié)點法向速度。對于切向摩擦力F??，若節(jié)點碰撞后速度v????不為零或碰撞后節(jié)點所受切向力F????大于最大靜摩擦力F?????(其中：F?????=k?F??，k?為靜摩擦系數(shù))，則節(jié)點受到切向的滑動摩擦力：

　　[F_{i}^{f}=F_{i}^{f d}=k_tl7lrfplleft|F_{i}^{n} ight| frac{v_{i}^{t a n}}{left|v_{i}^{t a n} ight|}]

　　其中k_d為動摩擦系數(shù)。若節(jié)點碰撞后速度v????為零且碰撞后節(jié)點所受切向力F????小于最大靜摩擦力F?????，則節(jié)點將保持靜止，對應(yīng)的切向摩擦力為：

　　[F_{i}^{f}=F_{i}^{f s}=-F_{i}^{t a n}]

　　3 光纖動力學(xué)仿真

　　本章主要基于光纖動力學(xué)模型的MATLAB代碼進行數(shù)值仿真分析，首先檢驗了光纖動力學(xué)模型的收斂性，然后基于該模型進行了光纖圓周布設(shè)及彎扭的仿真分析，最后設(shè)計了一種光纖陀螺儀裝配結(jié)構(gòu)，并通過光纖的動力學(xué)仿真來分析光纖沿軸向與切向的應(yīng)力分布及最大應(yīng)力變化。光纖這類單元尺寸較小、剛度較大的柔性體在有實時性要求的動力學(xué)仿真情況下容易出現(xiàn)穩(wěn)定性問題，而裝配仿真屬于實時性要求較高的場合。考慮到本章主要探究裝配過程中布設(shè)路徑對光纖應(yīng)力的影響，且彈性范圍內(nèi)材料的彈性模量只影響光纖應(yīng)力的大小，不影響光纖的應(yīng)力與布設(shè)路徑之間的相對關(guān)系，因此本節(jié)采用彈性模量較低的塑料光纖的模量進行仿真分析，本章所有仿真案例中涉及的主要參數(shù)如表1所示。

　　表1 光纖動力學(xué)仿真參數(shù)

　　| 參數(shù)名稱 | 參數(shù)值 |

　　| 密度ρ_f/(kg/m³) | 0.32×10³ |

　　| 移動速度v_f/(m/s) | 0.01 |

　　| 截面半徑r_f/m | 1×10?³ |

　　| 楊氏模量E_f/Pa | 1×10? |

　　| 泊松比λ_f | 0.3 |

　　| 阻力系數(shù) | 1×10?? |

　　| 靜摩擦系數(shù) | 0.3 |

　　| 動摩擦系數(shù) | 0.2 |

　　3.1 仿真收斂性檢驗

　　在進行光纖動力學(xué)仿真分析前，需要對動力學(xué)模型的網(wǎng)格尺寸(光纖離散化的每條線段長度)與時間步長進行收斂性檢驗。以純彎曲仿真為例，研究不同網(wǎng)格尺寸與時間步長下仿真結(jié)果的形狀變化與應(yīng)力分布。具體過程為：一根長度L_f=0.4m的光纖懸空并水平放置，然后控制光纖兩端以固定的移動速度相向移動，仿真時間為10s，其他相關(guān)的仿真參數(shù)采用表1中的參數(shù)值。

　　對于網(wǎng)格尺寸的收斂性研究，固定時間步長Δt為0.001s，網(wǎng)格尺寸分別取0.04m、0.02m和0.01m，其他仿真參數(shù)保持一致。結(jié)果顯示，隨著網(wǎng)格尺寸的等比減小，中心線形狀差異越來越小，且應(yīng)力分布的變化逐漸收斂。考慮到計算精度與效率，選取0.01m的網(wǎng)格尺寸作為后續(xù)仿真參數(shù)。

　　對于時間步長的收斂性研究，固定網(wǎng)格尺寸為0.01m，時間步長Δt分別取0.01s、0.001s和0.0001s，其他仿真參數(shù)保持一致。結(jié)果表明，時間步長對中心線形狀及應(yīng)力分布的影響很小，考慮到計算精度與效率，選取0.001s的時間步長作為后續(xù)的仿真參數(shù)。

　　3.2 仿真案例測試

　　光纖裝配過程中光纖自由部分的長度(即光纖控制端與光纖固定端之間的距離)會影響裝配過程中光纖尾纖布設(shè)的位置精度。為探究不同自由長度下光纖裝配精度的差異，進行光纖圓周布設(shè)的仿真分析。XY平面表示光纖陀螺儀中的一個用于固定光纖的結(jié)構(gòu)平面，光纖固定端定義為固定在XY平面上的光纖端部，光纖控制端定義為沿著半圓周軌跡移動的光纖端部，布設(shè)過程中光纖不斷下落至平面上并保持光纖控制端距XY平面的高度H不變。分別取高度H=0.04m、0.008m和0.002m，進行3種不同高度下的圓周布設(shè)仿真。其中，光纖控制端的水平移動速度v_xy=0.01m/s，豎直下落速度v_z=0.01m/s，圓周半徑R=0.1m，其它相關(guān)的仿真參數(shù)采用表1中的參數(shù)值。結(jié)果顯示，圓周布設(shè)的高度H越小，光纖圓周布設(shè)越來越接近目標(biāo)圓周曲線。因此，為保證光纖尾纖在裝配過程中具有較好的布設(shè)位置精度，光纖控制端與固定端的距離應(yīng)盡可能短，以便更好地控制光纖布設(shè)。

　　為研究光纖的彎扭特性，進行彎扭仿真測試分析。結(jié)果顯示了光纖在不同的圓周切向扭轉(zhuǎn)角度時的空間構(gòu)型的變化，同一光纖在固定端和控制端相對扭轉(zhuǎn)變形一周與四周的仿真模擬中，光纖長度L?取為0.4m，光纖截面半徑r_f取為1×10?³m，其他相關(guān)的仿真參數(shù)采用表1中的參數(shù)值，仿真過程分別為光纖兩端相對扭轉(zhuǎn)一周(總扭轉(zhuǎn)角度為360°)和四周(總扭轉(zhuǎn)角度為1440°)后再逐漸沿著水平直線分別靠近移動0.15m和0.06m(每端移動分別為0.15m和0.06m，兩端總靠近位移分別為0.3m和0.12m)。

　　在同一扭轉(zhuǎn)角度下，不同半徑的光纖會產(chǎn)生不同的空間變形效果。在光纖兩端相對扭轉(zhuǎn)四周、相向位移為0.12m的情況下，光纖半徑r_f=0.002m時可以實現(xiàn)光纖空間螺旋變形兩周的效果，但通過仿真發(fā)現(xiàn)當(dāng)光纖半徑r_f=0.001m時不能實現(xiàn)光纖空間螺旋變形兩周的效果。進一步分析不同光纖半徑下的拉伸力F?、彎曲力F?和扭轉(zhuǎn)力F?作用，發(fā)現(xiàn)在光纖半徑r_f=0.001m時，沿光纖中心線分布的拉伸力F?占主導(dǎo)作用，而彎曲力F?和扭轉(zhuǎn)力F?幾乎為零，從而使得光纖無法產(chǎn)生較大的扭轉(zhuǎn)變形;在光纖半徑r_f=0.002m時，沿光纖中心線分布的扭轉(zhuǎn)力F?占主導(dǎo)作用，其次是彎曲力F?與拉伸力F?，這也使得光纖受到扭轉(zhuǎn)力的作用較大，從而產(chǎn)生較大的扭轉(zhuǎn)變形。因此，在光纖長度和扭轉(zhuǎn)角度一定的情況下，光纖半徑越大，所產(chǎn)生的扭轉(zhuǎn)力也越大，光纖越容易扭轉(zhuǎn)變形。

　　3.3 光纖陀螺儀光路應(yīng)力分析

　　在光纖陀螺儀裝配過程中，光纖的盤纖過程是非常關(guān)鍵的環(huán)節(jié)。由于裝配過程的復(fù)雜性，其裝配工藝參數(shù)及盤纖路徑的設(shè)計目前仍主要基于人工手動裝配的經(jīng)驗積累，缺乏定量評估光纖尾纖應(yīng)力和光纖布設(shè)路徑的方法。本節(jié)基于光纖動力學(xué)模型對光纖裝配過程進行動力學(xué)仿真，研究裝配過程中光路的應(yīng)力分布及變化。

　　一種光纖陀螺儀光路裝配中探測器的光纖尾纖盤纖過程的切向應(yīng)力云圖，基于C++語言編程開發(fā)可視化操作界面進行仿真操作與顯示，并采用三維造型引擎ACIS和三維顯示交互工具包HOOPS建立三維環(huán)境。光纖裝配過程的變形及光纖沿軸向和切向的應(yīng)力分布首先由光纖動力學(xué)模型的MATLAB代碼仿真計算得到，然后通過C++與MATLAB的數(shù)據(jù)交互實現(xiàn)光纖在可視化操作界面上的裝配過程與應(yīng)力云圖顯示。其中，探測器及周圍結(jié)構(gòu)為Solidworks中的三維模型，其凹槽內(nèi)徑D=0.12m，凹槽深度H?=0.02m，主要仿真參數(shù)采用表1的參數(shù)值，對探測器的光纖尾纖盤纖過程進行裝配仿真，盤纖過程從探測器管殼端的光纖開始，在盤纖過程中光纖控制端高度H?=0.04m保持不變。

　　為進一步分析光纖陀螺儀中光路應(yīng)力分布及變化，對探測器光纖尾纖作軸向應(yīng)力和切向應(yīng)力沿中心線長度分布的變化曲線。結(jié)果顯示，光纖最大的軸向應(yīng)力σ?和切向應(yīng)力σ?均出現(xiàn)在光纖尾纖的端部。其中，光纖最大切向應(yīng)力遠(yuǎn)遠(yuǎn)大于最大軸向應(yīng)力，且光纖固定端的切向應(yīng)力遠(yuǎn)大于控制端。這是由于盤纖過程中光纖主要發(fā)生彎扭變形，而光纖受到的拉伸作用很小，且相較于光纖固定端，控制端的光纖彎扭程度相對較小。此外，由于盤纖過程近似為圓周軌跡運動且每一圈的圓周半徑近似相等，光纖中間部分的軸向和切向應(yīng)力較小且變化較均勻。

　　在光纖裝配過程中，光路的最大應(yīng)力也在隨時間變化。結(jié)果表明，最大軸向應(yīng)力與最大切向應(yīng)力在最開始盤纖過程中變化較大，隨后逐漸穩(wěn)定，但最大軸向應(yīng)力相較于最大切向應(yīng)力波動較大。這是因為光纖在剛開始盤纖時固定端附近的光纖彎扭變形較大，穩(wěn)定布設(shè)后光纖彎扭變形較均勻，而光纖拉伸變形在光纖固定端與控制端附近變形均較大，所以在光纖穩(wěn)定布設(shè)過程中最大軸向應(yīng)力受光纖控制端影響較大。此外，整個盤纖過程中切向應(yīng)力的最大值出現(xiàn)在盤纖最開始的時刻(即仿真時間t=0s)，對應(yīng)于初始狀態(tài)的切向應(yīng)力云圖，表明此時光纖彎扭變形程度最大;而軸向應(yīng)力的最大值出現(xiàn)在盤纖開始一段時間后(即仿真時間t=1.2s)，對應(yīng)于該時刻的軸向應(yīng)力云圖，表明此時光纖拉伸變形程度最大。這也為未來進一步定量設(shè)計最小應(yīng)力的光纖盤纖路徑提供了借鑒意義。

　　4 結(jié)束語

　　本文針對光纖陀螺儀中光纖尾纖裝配應(yīng)力定量控制的難題，以光纖陀螺儀的光纖尾纖為研究對象，對其裝配過程中的變形及應(yīng)力變化進行了動力學(xué)建模與仿真，主要研究結(jié)論如下：

　　1. 基于離散微分幾何理論建立了光纖運動學(xué)框架，該方法能在簡化光纖自由度、降低計算量的同時，保證計算結(jié)果具有良好的魯棒性與精度。隨后，基于該運動學(xué)框架由彈性勢能推導(dǎo)了光纖動力學(xué)方程，得到了描述光纖裝配過程中拉伸、彎曲和扭轉(zhuǎn)變形以及光纖應(yīng)力分布的動力學(xué)模型。

　　2. 基于Newmark隱式算法提出了一種求解光纖動力學(xué)模型的數(shù)值方法，并設(shè)計了光纖裝配過程中的碰撞檢測與響應(yīng)策略，實現(xiàn)了光纖裝配仿真的高效穩(wěn)定計算與物理意義上的接觸碰撞響應(yīng)效果。

　　3. 通過對光纖圓周布設(shè)的仿真發(fā)現(xiàn)，在光纖裝配過程中，光纖控制端與固定端之間自由長度越短，光纖布設(shè)的位置越精確;在光纖彎扭仿真過程中，同樣的光纖長度與扭轉(zhuǎn)角度下，光纖半徑越大使得扭轉(zhuǎn)力作用越大，光纖扭轉(zhuǎn)變形的程度越大。

　　4. 基于光纖動力學(xué)模型對一種光纖陀螺儀的盤纖過程進行了動力學(xué)仿真，得到了光纖盤纖過程中光路應(yīng)力分布與最大應(yīng)力變化曲線，確定了光纖軸向應(yīng)力與切向應(yīng)力最大時的光纖布設(shè)路徑，實現(xiàn)了光纖布設(shè)路徑下的軸向應(yīng)力與切向應(yīng)力的定量分析。

　　本文通過建立光纖動力學(xué)模型，對光纖陀螺儀裝配過程中光路應(yīng)力分布與變化進行了研究，為光纖陀螺儀低應(yīng)力裝配的定量控制提供了理論支撐。未來也將進一步研究小尺寸大剛度下光纖裝配過程的動力學(xué)仿真，并開發(fā)光纖裝配仿真系統(tǒng)，實現(xiàn)虛擬環(huán)境下光纖低應(yīng)力裝配的實時仿真與路徑規(guī)劃。

　　參考文獻

　　[1] WANG Wei. Interferometric fiber optic gyroscope technology[M]. Beijing: China Aerospace Press, 2010 (in Chinese). [王巍. 干涉型光纖陀螺儀技術(shù)[M]. 北京: 中國宇航出版社, 2010.]

　　[2] WU Changxin, YANG Mingwei, YANG Yuanhong, et al. Effect of fiber bend on interferometric fiber optical gyroscope performance[J]. Acta Optica Sinica, 2014, 34(3): 87-91 (in Chinese). [吳長莘, 楊明偉, 楊遠(yuǎn)洪, 等. 光纖彎曲對干涉式光纖陀螺性能的影響[J]. 光學(xué)學(xué)報, 2014, 34(3): 87-91.]

　　[3] WU Weifeng. Research on the bending loss of panda polarization maintaining fiber[D]. Nanjing: Nanjing University of Science and Technology, 2009 (in Chinese). [吳衛(wèi)鋒. 熊貓型保偏光纖彎曲損耗的研究[D]. 南京: 南京理工大學(xué), 2009.]

　　[4] YANG Yuanhong, YI Xiaosu, MENG Zhaokui. Experimental study on strain distribution in fiber coil used in fiber optic gyroscope[J]. Piezoelectrics and Acoustooptics, 2005, 27(2): 98-101 (in Chinese). [楊遠(yuǎn)洪, 伊小素, 孟照魁. 光纖陀螺用光纖環(huán)的應(yīng)力分布實驗研究[J]. 壓電與聲光, 2005, 27(2): 98-101.]

　　[5] YIN Qiqi, LI Jing, ZHANG Zhihua. Experimental study on the stress in fiber sensing coil[J]. Navigation and Control, 2016, 15(2): 79, 98-103 (in Chinese). [尹其其, 李晶, 張智華. 光纖環(huán)內(nèi)部應(yīng)力試驗研究[J]. 導(dǎo)航與控制, 2016, 15(2): 79, 98-103.]

　　[6] QIU Hongfang. Analysis on influence of inflection point in optical fiber loop on performance of optical fiber gyroscope[J]. Tianjin Science & Technology, 2023, 50(5): 23-25 (in Chinese). [邱紅芳. 光纖環(huán)圈拐點對光纖陀螺性能影響的分析[J]. 天津科技, 2023, 50(5): 23-25.]

　　[7] YANG Xuewei, DU Xinzheng, LIN Heng, et al. Thermal-stress effect of package optical fiber coil adhesive on optical fiber[J]. Journal of Chinese Inertial Technology, 2003(6): 95-98 (in Chinese). [楊學(xué)圍, 杜新政, 林恒, 等. 封裝光纖線圈的膠粘劑對光纖產(chǎn)生的熱應(yīng)力影響[J]. 中國慣性技術(shù)學(xué)報, 2003(6): 95-98.]

　　[8] TANG Xu, ZHU Yazhou, SHEN Zhongxiang, et al. Research on application of optical fiber sensing material in ship structural health monitoring[J]. Ship Science and Technology, 2020, 42(23): 181-185 (in Chinese). [唐旭, 朱亞洲, 沈中祥, 等. 光纖傳感材料在船舶結(jié)構(gòu)健康監(jiān)測中的應(yīng)用[J]. 艦船科學(xué)技術(shù), 2020, 42(23): 181-185.]

　　[9] LV N J, LIU J H, XIA H X, et al. A review of techniques for modeling flexible cables[J]. Computer-Aided Design, 2020, 122: 102826. DOI: 10.1016/j.cad.2020.102826.

　　[10] JAWED M K, NOVELIA A, O'REILLY O M. A primer on the kinematics of discrete elastic rods[M]. Berlin, Germany: Springer-Verlag, 2018.

　　[11] ZHANG Y, CHEN D S, QIAN H L. Computational method for the deformation mechanism of non-prestressed cable net structures based on the vector form intrinsic finite element method[J]. Engineering Structures, 2021, 231: 111788. DOI: 10.1016/j.engstruct.2020.111788.

　　[12] LIN Haili, LIU Jianhua, TANG Chengtong, et al. Physical characteristic oriented modeling for flexible cable harness based on Cosserat elastic rod theory[J]. Journal of Graphics, 2016, 37(1): 34-42 (in Chinese). [林海立, 劉檢華, 唐承統(tǒng), 等. 基于Cosserat彈性桿理論的柔性線纜物理建模方法[J]. 圖學(xué)學(xué)報, 2016, 37(1): 34-42.]

　　[13] LIU L, YAN Z, HE Y M, et al. On the configuration evolution of soft filaments under combined tension and torsion[J]. Acta Mechanica Sinica, 2023, 39(11): 222498. DOI: 10.1007/s10409-023-22498-x.

　　[14] BERGOU M, AUDOLY B, VOUGA E, et al. Discrete viscous threads[J]. ACM Transactions on Graphics, 2010, 29(4): 116. DOI: 10.1145/1778765.1778853.

　　[15] DESBRUN M, POLTHIER K, SCHRÖDER P, et al. Discrete differential geometry: An applied introduction[EB/OL]. (2006-05-11) [2023-11-29]. http://geometry.caltech.edu/pubs/GSD06.pdf.

　　[16] BERGOU M, WARDETZKY M, ROBINSON S, et al. Discrete elastic rods[J]. ACM Transactions on Graphics, 2008, 27(3): 63. DOI: 10.1145/1399504.1360662.

　　[17] LV N J, LIU J H, JIA Y Y. Dynamic modeling and control of deformable linear objects for single-arm and dual-arm robot manipulations[J]. IEEE Transactions on Robotics, 2022, 38(4): 2341-2353.

　　[18] AUDOLY B, CLAUVELIN N, BRUN P T, et al. A discrete geometric approach for simulating the dynamics of thin viscous threads[J]. Journal of Computational Physics, 2013, 253: 18-49.

　　[19] SPERL G, SÁNCHEZ-BANDERAS R M, LI M, et al. Estimation of yarn-level simulation models for production fabrics[J]. ACM Transactions on Graphics (TOG), 2022, 41(4): 65. DOI: 10.1145/3528223.3530167.

　　[20] ARMANINI C, BOYER F, MATHEW A T, et al. Soft robots modeling: A structured overview[J]. IEEE Transactions on Robotics, 2023, 39(3): 1728-1748.

　　[21] GREGOIRE M, SCHÖMER E. Interactive simulation of one-dimensional flexible parts[J]. Computer-Aided Design, 2007, 39(8): 694-707.

　　[22] LANG H, LIN N J, ARNOLD M. Multi-body dynamics simulation of geometrically exact Cosserat rods[J]. Multibody System Dynamics, 2011, 25(3): 285-312.

　　[23] LOOCK A, SCHÖMER E, STADTWALDI. A virtual environment for interactive assembly simulation: From rigid bodies to deformable cables[EB/OL].[2023-11-29]. https://www.staff.uni-mainz.de/schoemer/publications/SCI01.pdf.

　　[24] NEWMARK N M. A method of computation for structural dynamics[J]. Journal of the Engineering Mechanics Division, 1959, 85(3): 67-94.

　　[25] YANG Xiaodong, LIU Jianhua, ZHAO Yingfeng, et al. Rigid and flexible hybrid assembly process simulation technology for complex mechatronic products[J]. Computer Integrated Manufacturing Systems, 2019, 25(2): 340-349 (in Chinese). [楊嘯東, 劉檢華, 趙瑛峰, 等. 面向復(fù)雜機電產(chǎn)品的剛?cè)峄旌涎b配過程仿真技術(shù)[J]. 計算機集成制造系統(tǒng), 2019, 25(2): 340-349.]