成茵瑛

（同濟(jì)大學(xué)電子與信息工程學(xué)院，上海201804）

Shupe效應(yīng)非互易相移按匝傳輸?shù)木_分析與計(jì)算

成茵瑛

（同濟(jì)大學(xué)電子與信息工程學(xué)院，上海201804）

為進(jìn)一步提高干涉型光纖陀螺（IFOG）的Shupe效應(yīng)非互易相移的計(jì)算精度，在柱坐標(biāo)系中，根據(jù)熱傳遞的原理和過程，以1匝光纖線圈為基本單位，在MATLAB平臺上建立了隨著軸向、徑向和時間變化的IFOG的三維溫度場分布模型。根據(jù)Shupe效應(yīng)的原理，基于光波在光纖線圈中的實(shí)際傳輸過程，對熱傳遞動態(tài)過程中的時間元素進(jìn)行進(jìn)一步分解，細(xì)化為每一匝傳輸?shù)亩稍綍r間，用插值法使溫度場的熱傳遞按匝進(jìn)行，得到了不同時間點(diǎn)處每匝光纖線圈的溫度分布。設(shè)計(jì)了按匝傳輸分析計(jì)算非互易相移的新方法。對四極對稱繞法和交叉對稱繞法進(jìn)行了仿真數(shù)值計(jì)算，結(jié)果表明與原算法相比，新算法的Shupe效應(yīng)非互易相移計(jì)算精度提高了3個量級，符合預(yù)期。

干涉型光纖陀螺；按匝傳輸；Shupe效應(yīng)；三維溫度場；四極對稱繞法；交叉對稱繞法；非互易相移；計(jì)算精度

0 引言

干涉型光纖陀螺是發(fā)展最早、技術(shù)最成熟、應(yīng)用最廣的光纖陀螺，其慣性級產(chǎn)品的研制和開發(fā)逐漸成熟。在國外，研制光纖陀螺的部門和機(jī)構(gòu)主要集中在歐美與日本，水平在全世界領(lǐng)先，且產(chǎn)量最多。與國外研制水平相比，我國高溫度穩(wěn)定性、高精度級光纖陀螺的研制、實(shí)現(xiàn)與應(yīng)用還有眾多難題需要解決。其中，最難的是光纖陀螺的工程化，首先要解決是陀螺性能受溫度影響的問題。當(dāng)溫度發(fā)生變化時，光纖陀螺的輸出會發(fā)生一定漂移（溫漂），顯著的溫漂會明顯降低光纖陀螺輸出的精度，阻礙了其發(fā)展及實(shí)際工程應(yīng)用。

為提高光纖陀螺的溫度穩(wěn)定性，解決光纖陀螺的環(huán)境適應(yīng)性問題，近年來國內(nèi)對該問題進(jìn)行了理論分析與研究，并提出了多種抑制溫度效應(yīng)的方法；部分從光纖敏感環(huán)纏繞方式、固膠技術(shù)、光學(xué)器件選取，以及光纖敏感環(huán)骨架選擇等方面對陀螺結(jié)構(gòu)進(jìn)行優(yōu)化改進(jìn)以抑制陀螺的溫度效應(yīng)；部分則直接對光纖陀螺加溫控或選用不同方法對溫漂進(jìn)行建模補(bǔ)償［1］。傳統(tǒng)的中高精度干涉型光纖陀螺雖然采用了對稱繞法和環(huán)內(nèi)熱“短路”環(huán)外熱“斷路”的結(jié)構(gòu)，并對溫度漂移進(jìn)行建模補(bǔ)償，但光纖敏感環(huán)的熱非互易誤差仍是主要的誤差源。上述方法都是從物理角度進(jìn)行改善，缺少在度量方法上的改進(jìn)。由文獻(xiàn)報(bào)道可知：通常在建立光纖線圈的三維溫度模型的基礎(chǔ)上，計(jì)算Shupe效應(yīng)非互易相移時，常從空間角度對公式作進(jìn)一步分解，試圖獲得更精確的分析和更高精度的結(jié)果，即對原公式中的積分限做長度上的分割，而被積分項(xiàng)不變。本文從空間和時間兩個角度，對公式作進(jìn)一步分解和推導(dǎo)。對積分項(xiàng)，不選擇導(dǎo)數(shù)形式等條件苛刻的極限形式，而是依據(jù)光波在光纖線圈中傳播的微觀詳細(xì)過程，選擇溫度差值形式。由于光纖線圈中的溫度場是時變的，本文首先建立光纖環(huán)的三維溫度場分布模型，以模型的最小單位為每匝的溫度。對經(jīng)典的四極對稱繞法和交叉對稱繞法，研究一種基于渡越時間的按匝傳輸?shù)腟hupe效應(yīng)非互易相移計(jì)算方法，對傳統(tǒng)的計(jì)算方法精度進(jìn)行改善，使結(jié)果更貼近實(shí)際變化過程，預(yù)期能使Shupe效應(yīng)非互易相移結(jié)果的精度提高1個量級以上。在MATLAB仿真軟件平臺進(jìn)行了初步仿真建模，根據(jù)仿真結(jié)果，對兩種算法的性能進(jìn)行比較。

1 光纖環(huán)的三維溫度場

實(shí)際的光纖陀螺線圈框架呈圓柱形，以圓柱的底面圓心為原點(diǎn)建立柱坐標(biāo)系。因溫度在框架圓周角度方向?qū)ΨQ，建模時僅考慮徑向和軸向的溫度變化，即溫度場分布T（r，z，t）僅與半徑坐標(biāo)r、軸向坐標(biāo)z和時間t三個變量有關(guān)［2］。由于骨架的對稱性，研究中只考慮對稱軸一半部分。添加溫度條件的四極對稱繞法和交叉對稱繞法半剖面如圖1所示。圖1中：a為骨架內(nèi)半徑；b為骨架外半徑；c為骨架高度；Tx，Tg，Te分別為線圈度、骨架和外部環(huán)境的初始溫度；藍(lán)色段表示0～L／2的光纖段，白色段表示L／2～L的光纖段；數(shù)字表示該層的匝序號。此處L為光纖總長。

根據(jù)熱傳遞原理建模?？紤]仿真中數(shù)據(jù)量過多時會降低仿真效率，假設(shè)光纖陀螺上繞制的光纖共16層，每層繞光纖20匝。光纖環(huán)溫度模型中參量為：層數(shù)N＝16；每層匝數(shù)M＝20；線圈骨架半徑a＝62mm；光纖半徑r0＝0.081mm；熱傳導(dǎo)系數(shù)k＝0.32W／（m·K）；輻射換熱角系數(shù)F＝1；斯蒂芬玻爾茲曼常數(shù)σ＝5.67×10－8W／（m2·K4）；光纖有效折射率n0＝1.46；光纖密度ρ＝1 600kg／m3；光纖比熱容cc＝1 700J／（kg·K）；光纖表面發(fā)射率ε＝0.15；Tx＝273.15K；Tg＝293.15K；Te＝313.15K；L＝150m［3］。由此建立光纖線圈的熱傳導(dǎo)微分方程

式中：h為等效的表面?zhèn)鳠嵯禂?shù)，且h＝σεF；k為光纖敏感線圈的歸一化熱傳導(dǎo)系數(shù)；α為熱擴(kuò)散系數(shù)，且α＝k／（ρcc）。

微分方程式（1）是一個方程與邊界條件都為非齊次的方程。取進(jìn)行輻射換熱線性化。用貝塞爾函數(shù)、貝塞爾方程、傅里葉貝塞爾級數(shù)等可解得上述方程的靜態(tài)解Ts和動態(tài)解Td分別為

式中：I0，K0分別為第一類和第二類修正貝塞爾函數(shù)；c為真空光速；

此處：

其中：J0，Y0分別為第一類和第二類零階貝塞爾函數(shù)；J1，Y1分別為第一類和第二類一階貝塞爾函數(shù)；

因此，總溫度場分布可表示為

根據(jù)上述熱傳導(dǎo)微分方程組推導(dǎo)的結(jié)果，用MATLAB軟件仿真所得時間t＝14s時柱坐標(biāo)系中光纖陀螺的三維溫度場的分布如圖2所示。

2 按匝傳輸?shù)腟hupe效應(yīng)非互易相移計(jì)算方法

Shupe效應(yīng)是指光沿長度為L的光纖傳播時，由于存在時變的溫度場，兩束反向傳播的光波在不同時間經(jīng)過這段光纖，因溫度擾動而經(jīng)歷不同的相移［11］。一般情況下Shupe非互易相移可表示為

式中：β0＝2π／為自由空間的傳播常數(shù)；ΔT（t，l）為沿光纖的溫度分布的漲落部分；為光纖的溫度系數(shù)；αh為光纖熱膨脹系數(shù)。此處：為寬帶光源的平均光波長。

假設(shè)在時間t處順時針光波到達(dá)光纖線圈輸出端，則順時針波（CW）和逆時針波（CCW）在整個光纖線圈中的相位延遲分別為

式中：cm為光纖中的光傳播速度，且cm＝c／n0。CW，CCW在整個光纖線圈中的相位差可表示為

由函數(shù)導(dǎo)數(shù)定義

可得相位差傳統(tǒng)計(jì)算公式為

式中：lis，lie分別為第i匝光纖始端和末端距離光纖起點(diǎn)的坐標(biāo)，且

式（16）要求Δt趨于0，才可轉(zhuǎn)化為溫度隨時間的導(dǎo)數(shù)計(jì)算。但通常溫度隨時間的導(dǎo)數(shù)難以采集，而且仿真時Δt不可能趨于0，因此獲得的溫度隨時間的導(dǎo)數(shù)不夠精確。另從原理分析，在式（17）中，非互易相移大小不可能與光纖長度的二次方成正比。因此，為獲得更精確的結(jié)果，本文計(jì)算Shupe效應(yīng)引起的非互易相移時采用以匝為基本單位，細(xì)化為每一匝傳輸?shù)亩稍綍r間，通過插值方法使溫度場的熱傳遞按匝進(jìn)行。這樣，就可得每個時間點(diǎn)每一匝處的溫度。具體算法如下。

三維溫度場中計(jì)算Shupe非互易相移時熱漲落對光傳播的影響既基于折射率的溫度依賴性，又依賴介質(zhì)的熱膨脹。但熱膨脹系數(shù)隨溫度的變化率較折射率至少小1個量級，故可忽略。則相位大小可表示為

假設(shè)以每匝為基本單位展開，則一匝內(nèi)溫度漲落相同，式（18）可展開為

式中：Q為總匝數(shù)。

式（19）的關(guān)鍵是關(guān)于第i匝溫度漲落ΔT（t，li）?？筛鶕?jù)已建立的三維溫度場分布模型，在總時間區(qū)間內(nèi)進(jìn)行插值，插值的最小單位為1匝的傳輸時間，這樣可使溫度分布能精確到每匝。以四極對稱繞法中的第4層第13匝光纖線圈為例（物理模型如圖1（a）所示），具體過程如下。

a）選定計(jì)算Shupe效應(yīng)非互易相移的時間起點(diǎn)td，光在每匝線圈中傳播的渡越時間為ti＝Δli／cm。分別以ti為時間增量點(diǎn)對已經(jīng)建立的溫度場T（r，z，t）在td后進(jìn)行順序插值，其中Δli＝lie－lis為第i匝光纖段的長度，順次排序從光纖的一端到另一端，得到在每個tj＝ti＋td處的溫度場分布。插值之前第4層第13匝光纖線圈從0.5～30s間的溫度變化如圖3所示。在本文建立的模型中，每匝的平均渡越時間為1.568 54ns，光波在長L的光纖中傳播的渡越時間總時長為0.501 932μs。由于按匝傳輸?shù)臅r間很短，溫度變化很微小，為使表示更清楚，對圖進(jìn)行坐標(biāo)平移和放大并取對數(shù)坐標(biāo)，則Tsh＝（T－Tmin）×107，tsh＝（t－tmin）×108。此處：Tmin，ts分別為溫度場溫度最小值和計(jì)算時間起點(diǎn)。四極對稱繞法光纖環(huán)按匝插值后，第4層第13匝光纖線圈分別以0.5，5.0，15.0，25.0s為時間起點(diǎn)，其后一段時間內(nèi)（光波從光纖的一端傳至另一端的總時長）的溫度及其一階導(dǎo)數(shù)分別如圖4、5所示。

由圖4、5可知：在熱載荷剛影響光纖溫度分布的初期（t＝0.5s時），溫度上升非?？欤仙厔菔置黠@；當(dāng)t＝5.0s時，上升速率略微降低，約10－1K／s；當(dāng)t＝15.0s時，溫度上升的趨勢明顯減慢，在光波從光纖的起始端到末端的傳播時長內(nèi)溫度上升十分緩慢，速率約10－2K／s；當(dāng)t＝25.0s時，溫度趨于平穩(wěn)，此時溫度傳遞已完成，各匝的溫度值基本確定。綜上分析可得，模型的此變化趨勢符合實(shí)際的熱傳遞。

b）由于T（r，z，ti）是MN個時間點(diǎn)上的溫度分布，本文需要的是特定匝在特定時間點(diǎn)上的溫度值，用公式

進(jìn)行轉(zhuǎn)化變量。此處：p，q為溫度場中徑向和軸向的坐標(biāo)序號；Tt為轉(zhuǎn)化后只與變量r，z有關(guān)的溫度分布場。四極對稱繞法下，起始時刻為10s時，以匝渡越時間為時間增量熱傳遞后的光纖環(huán)各匝在光波傳播至此匝時的溫度分布場如圖6所示。

綜合上述，由

求得第i匝溫度漲落大小。此處：p，q為第i匝所在的位置，即第p層、第q匝。

3 經(jīng)典四極繞法與交叉對稱繞法實(shí)驗(yàn)

為驗(yàn)證新算法能提高非互易相移計(jì)算精度，在MATLAB軟件平臺中進(jìn)行仿真實(shí)驗(yàn)。以建立的16層，每層20匝，總長度為L的光纖環(huán)溫度場模型為例，四極對稱和交叉對稱繞法分別如圖1（a）、（b）所示。分別用傳統(tǒng)算法和按匝傳輸精確分析算法，計(jì)算四極對稱繞法和交叉對稱繞法兩種繞法的Shupe效應(yīng)引起的非互易相移。為更直觀表現(xiàn)兩種算法所得的非互易相移，取對數(shù)坐標(biāo)，因零無對數(shù)，取近似值10－100代入。兩種繞法的不同算法結(jié)果分別如圖7、8所示。四極對稱繞法的光纖在t＝0.5，30.0s時溫度隨徑向的分布如圖9所示。

分析上述實(shí)結(jié)果，可得以下結(jié)論。

a）對四極對稱繞法，在溫度開始傳遞的初期，四極對稱繞法中各層的溫度分布未能達(dá)到線性（如圖9所示），因此關(guān)于光纖中點(diǎn)對稱的兩段光纖經(jīng)歷的溫度變化并不相同，此時的Shupe效應(yīng)非互易相移理論值并不為零。當(dāng)溫度傳遞結(jié)束，各層溫度分布逐漸穩(wěn)定且隨徑向呈現(xiàn)線性分布時（如圖9所示），Shupe非互易相移理論值應(yīng)為零。圖7（a）中，相對理論值零，傳統(tǒng)算法的精度可達(dá)10－11rad，而按匝傳輸?shù)木_分析算法精度可達(dá)10－13rad，可見新算法精度可高3個量級。

b）對交叉對稱繞法，由于溫度分布的對稱中心與光纖的對稱中心重合。即在軸向上溫度分布關(guān)于軸向中點(diǎn)對稱，恰好對應(yīng)光纖關(guān)于中點(diǎn)對稱的兩段，因此在溫度場熱傳遞初期和達(dá)到穩(wěn)定后的溫度分布，任何時刻交叉對稱繞法下的Shupe效應(yīng)非互易

相移理論值都應(yīng)為零。圖8中，相對理論值零，傳統(tǒng)算法的精度最優(yōu)可達(dá)10－24rad，按匝傳輸?shù)木_分析算法的結(jié)果等于理論值0rad，可見新算法的精度更高，結(jié)果更準(zhǔn)確。

c）由圖7、8可知：在溫度傳遞初期，四極對稱繞法對Shupe效應(yīng)非互易相移的抑制效果遠(yuǎn)不如交叉對稱繞法。同一個時間點(diǎn)處，交叉對稱繞法利用溫度分布對稱中心和光纖環(huán)長度對稱中心的重合，大幅抑制了非互易相移，使其減小約10個量級。

4 結(jié)論

本文基于柱坐標(biāo)系，用MATLAB仿真軟件建立了IFOG的三維溫度分布場，描述了在不同時間點(diǎn)上熱傳遞過程。傳統(tǒng)Shupe效應(yīng)非互易相移計(jì)算方法是直接代入CW，CCW產(chǎn)生的相位差公式，涉及到溫度場的時間導(dǎo)數(shù)等，但溫度隨時間的導(dǎo)數(shù)精度有限，且與抽樣時間間隔有關(guān)。如時間間隔很大，就無法用溫度的時間導(dǎo)數(shù)代替公式中的差值比項(xiàng)，從而導(dǎo)致非互易相移結(jié)果有偏差。新型計(jì)算方法中不需要溫度場的時間導(dǎo)數(shù)值，僅要求獲知某匝處的溫度值，用兩匝間的溫度差代替溫度的時間導(dǎo)數(shù)。由此本文提出了按每匝渡越時間進(jìn)行溫度插值的Shupe效應(yīng)非互易相移計(jì)算方法，并在MATLAB仿真平臺上用新型算法和傳統(tǒng)算法分別對經(jīng)典四極繞法和交叉對稱繞法進(jìn)行了仿真，根據(jù)某個確定的時間起點(diǎn)，以一匝為基本單位分別計(jì)算了兩種繞法的Shupe效應(yīng)產(chǎn)生的非互易相移。結(jié)果表明：與傳統(tǒng)計(jì)算方法相比，按匝傳輸插值的Shupe效應(yīng)非互易相移計(jì)算方法精度可有3個量級以上的改善。

本文的新算法提高了度量Shupe效應(yīng)導(dǎo)致的非互易相移大小方法的精度，相較于前人文獻(xiàn)中對光纖線圈長度分割的方法，進(jìn)一步提出了時間維度上的分割。但本文研究也存在一些問題：首先，鑒于仿真條件的限制，每層光纖線圈只取了16層，與現(xiàn)實(shí)應(yīng)用中的光纖陀螺的物理參數(shù)有差距，后續(xù)可進(jìn)一步增加層數(shù)，并在真實(shí)的光纖陀螺實(shí)驗(yàn)環(huán)境中對其進(jìn)行驗(yàn)證；其次，雖然本文有助于從更微觀的角度了解Shupe效應(yīng)原理，但未明確提出抑制的方法，后續(xù)研究可在本文算法的基礎(chǔ)上，研究抑制Shupe效應(yīng)非互易相移的方法。

［1］ 張桂才.光纖陀螺原理與技術(shù)［M］.北京：國防工業(yè)出版社，2008：156－178.

［2］ 王巍，向政，邢輝，等.光纖陀螺光纖敏感線圈三維溫度瞬態(tài)響應(yīng)模型［J］.中國慣性技術(shù)學(xué)報(bào)，2012，20（5）：605－601.

［3］ 李茂春，劉鐵根，王佳佳，等.光纖陀螺光纖環(huán)Shupe效應(yīng)的三維有限元模型［J］.天津大學(xué)學(xué)報(bào)，2008，41（5）：601－606.

［4］ 奚定平.貝塞爾函數(shù)［M］.北京：高等教育出版社，1998.

［5］ 陸全康.數(shù)學(xué)物理方法自學(xué)輔導(dǎo)［M］.上海：上?？茖W(xué)技術(shù)出版社，1989.

［6］ 俞昌銘.熱傳導(dǎo)及其數(shù)值分析［M］.北京：清華大學(xué)出版社，1981.

［7］ 因克羅普拉，德威特.傳熱基礎(chǔ)［M］.北京：宇航出版社，1987.

［8］ 梁昆淼.數(shù)學(xué)物理方法［M］.北京：高等教育出版社，1998.

［9］ 方瑛，黃毅.數(shù)學(xué)物理方程與特殊函數(shù)［M］.北京：科學(xué)出版社，2012.

［10］ SHUPE D.Thermally induced nonreciprocity in the fiber－optic interferometer［J］.Applied Optics，1980，19（5）：654－655.

［11］ 馮文帥，王學(xué)鋒，王巍.光纖環(huán)的每層匝數(shù)不同對光纖陀螺溫度性能影響［J］.中國慣性技術(shù)學(xué)報(bào)，2011，19（4）：487－493.

Accurate Analysis and Calculation of Nonreciprocal Phase Shift Caused by Shupe Effect According to Transmission by Circles

CHENG Yin－ying
（School of Electronics and Information Engineering，Tongji University，Shanghai 201804，China）

To further improve the computation accuracy of non－reciprocal phase shift caused by the Shupe effect of interferometric fiber optic gyroscope（IFOG），a three－dimensional temperature field distribution model of IFOG was established under the cylindrical coordinate system which was based on the basic unit of circle and varied by axis direction，radial direction and time according to the principles and process of heat transmission.Based on specific microscopic process in light transmission in the optical fiber coil，a novel calculating method of non－reciprocal phase shift was introduced according to the principle of Shupe effects，which was defined by the transmission time of each circle and the heat transmission of the temperature field was calculated for each circle by interpolation methed.Time element in dynamic process of heat transmission was decomposed to obtain the temperature distribution of each circle of fiber coil at different time point.Then some numerical calculation experiments were simulated for quadrupole symmetry winding method and crossing symmetry winding method.The computation accuracy of non－reciprocal phase shift caused by the Shupe effect was increased by three orders of magnitude compared with the original algorithm.The results conformed to the expectations.

interferometric fiber optic gyroscope；transmission by circles；Shupe effect；three－dimensional temperature field distribution；quadrupole symmetry winding method；crossing symmetry winding method；nonreciprocal phase shift；calculation accuracy

V241.5

A

10.19328／j.cnki.1006－1630.2017.01.016

1006－1630（2017）01－0097－08

2016－06－21；

2016－07－10

上海航天科技創(chuàng)新基金資助（SAST2015076）

成茵瑛（1991—），女，碩士生，主要研究方向?yàn)楣饫w陀螺的零偏誤差優(yōu)化。