Shupe效应非互易相移按匝传输的精确分析与计算

成茵瑛

（同济大学电子与信息工程学院，上海201804）

Shupe效应非互易相移按匝传输的精确分析与计算

成茵瑛

（同济大学电子与信息工程学院，上海201804）

为进一步提高干涉型光纤陀螺（IFOG）的Shupe效应非互易相移的计算精度，在柱坐标系中，根据热传递的原理和过程，以1匝光纤线圈为基本单位，在MATLAB平台上建立了随着轴向、径向和时间变化的IFOG的三维温度场分布模型。根据Shupe效应的原理，基于光波在光纤线圈中的实际传输过程，对热传递动态过程中的时间元素进行进一步分解，细化为每一匝传输的渡越时间，用插值法使温度场的热传递按匝进行，得到了不同时间点处每匝光纤线圈的温度分布。设计了按匝传输分析计算非互易相移的新方法。对四极对称绕法和交叉对称绕法进行了仿真数值计算，结果表明与原算法相比，新算法的Shupe效应非互易相移计算精度提高了3个量级，符合预期。

干涉型光纤陀螺；按匝传输；Shupe效应；三维温度场；四极对称绕法；交叉对称绕法；非互易相移；计算精度

0 引言

干涉型光纤陀螺是发展最早、技术最成熟、应用最广的光纤陀螺，其惯性级产品的研制和开发逐渐成熟。在国外，研制光纤陀螺的部门和机构主要集中在欧美与日本，水平在全世界领先，且产量最多。与国外研制水平相比，我国高温度稳定性、高精度级光纤陀螺的研制、实现与应用还有众多难题需要解决。其中，最难的是光纤陀螺的工程化，首先要解决是陀螺性能受温度影响的问题。当温度发生变化时，光纤陀螺的输出会发生一定漂移（温漂），显著的温漂会明显降低光纤陀螺输出的精度，阻碍了其发展及实际工程应用。

为提高光纤陀螺的温度稳定性，解决光纤陀螺的环境适应性问题，近年来国内对该问题进行了理论分析与研究，并提出了多种抑制温度效应的方法；部分从光纤敏感环缠绕方式、固胶技术、光学器件选取，以及光纤敏感环骨架选择等方面对陀螺结构进行优化改进以抑制陀螺的温度效应；部分则直接对光纤陀螺加温控或选用不同方法对温漂进行建模补偿［1］。传统的中高精度干涉型光纤陀螺虽然采用了对称绕法和环内热“短路”环外热“断路”的结构，并对温度漂移进行建模补偿，但光纤敏感环的热非互易误差仍是主要的误差源。上述方法都是从物理角度进行改善，缺少在度量方法上的改进。由文献报道可知：通常在建立光纤线圈的三维温度模型的基础上，计算Shupe效应非互易相移时，常从空间角度对公式作进一步分解，试图获得更精确的分析和更高精度的结果，即对原公式中的积分限做长度上的分割，而被积分项不变。本文从空间和时间两个角度，对公式作进一步分解和推导。对积分项，不选择导数形式等条件苛刻的极限形式，而是依据光波在光纤线圈中传播的微观详细过程，选择温度差值形式。由于光纤线圈中的温度场是时变的，本文首先建立光纤环的三维温度场分布模型，以模型的最小单位为每匝的温度。对经典的四极对称绕法和交叉对称绕法，研究一种基于渡越时间的按匝传输的Shupe效应非互易相移计算方法，对传统的计算方法精度进行改善，使结果更贴近实际变化过程，预期能使Shupe效应非互易相移结果的精度提高1个量级以上。在MATLAB仿真软件平台进行了初步仿真建模，根据仿真结果，对两种算法的性能进行比较。

1 光纤环的三维温度场

实际的光纤陀螺线圈框架呈圆柱形，以圆柱的底面圆心为原点建立柱坐标系。因温度在框架圆周角度方向对称，建模时仅考虑径向和轴向的温度变化，即温度场分布T（r，z，t）仅与半径坐标r、轴向坐标z和时间t三个变量有关［2］。由于骨架的对称性，研究中只考虑对称轴一半部分。添加温度条件的四极对称绕法和交叉对称绕法半剖面如图1所示。图1中：a为骨架内半径；b为骨架外半径；c为骨架高度；Tx，Tg，Te分别为线圈度、骨架和外部环境的初始温度；蓝色段表示0～L／2的光纤段，白色段表示L／2～L的光纤段；数字表示该层的匝序号。此处L为光纤总长。

根据热传递原理建模。考虑仿真中数据量过多时会降低仿真效率，假设光纤陀螺上绕制的光纤共16层，每层绕光纤20匝。光纤环温度模型中参量为：层数N＝16；每层匝数M＝20；线圈骨架半径a＝62mm；光纤半径r0＝0.081mm；热传导系数k＝0.32W／（m·K）；辐射换热角系数F＝1；斯蒂芬玻尔兹曼常数σ＝5.67×10－8W／（m2·K4）；光纤有效折射率n0＝1.46；光纤密度ρ＝1 600kg／m3；光纤比热容cc＝1 700J／（kg·K）；光纤表面发射率ε＝0.15；Tx＝273.15K；Tg＝293.15K；Te＝313.15K；L＝150m［3］。由此建立光纤线圈的热传导微分方程

式中：h为等效的表面传热系数，且h＝σεF；k为光纤敏感线圈的归一化热传导系数；α为热扩散系数，且α＝k／（ρcc）。

微分方程式（1）是一个方程与边界条件都为非齐次的方程。取进行辐射换热线性化。用贝塞尔函数、贝塞尔方程、傅里叶贝塞尔级数等可解得上述方程的静态解Ts和动态解Td分别为

式中：I0，K0分别为第一类和第二类修正贝塞尔函数；c为真空光速；

此处：

其中：J0，Y0分别为第一类和第二类零阶贝塞尔函数；J1，Y1分别为第一类和第二类一阶贝塞尔函数；

因此，总温度场分布可表示为

根据上述热传导微分方程组推导的结果，用MATLAB软件仿真所得时间t＝14s时柱坐标系中光纤陀螺的三维温度场的分布如图2所示。

2 按匝传输的Shupe效应非互易相移计算方法

Shupe效应是指光沿长度为L的光纤传播时，由于存在时变的温度场，两束反向传播的光波在不同时间经过这段光纤，因温度扰动而经历不同的相移［11］。一般情况下Shupe非互易相移可表示为

式中：β0＝2π／为自由空间的传播常数；ΔT（t，l）为沿光纤的温度分布的涨落部分；为光纤的温度系数；αh为光纤热膨胀系数。此处：为宽带光源的平均光波长。

假设在时间t处顺时针光波到达光纤线圈输出端，则顺时针波（CW）和逆时针波（CCW）在整个光纤线圈中的相位延迟分别为

式中：cm为光纤中的光传播速度，且cm＝c／n0。CW，CCW在整个光纤线圈中的相位差可表示为

由函数导数定义

可得相位差传统计算公式为

式中：lis，lie分别为第i匝光纤始端和末端距离光纤起点的坐标，且

式（16）要求Δt趋于0，才可转化为温度随时间的导数计算。但通常温度随时间的导数难以采集，而且仿真时Δt不可能趋于0，因此获得的温度随时间的导数不够精确。另从原理分析，在式（17）中，非互易相移大小不可能与光纤长度的二次方成正比。因此，为获得更精确的结果，本文计算Shupe效应引起的非互易相移时采用以匝为基本单位，细化为每一匝传输的渡越时间，通过插值方法使温度场的热传递按匝进行。这样，就可得每个时间点每一匝处的温度。具体算法如下。

三维温度场中计算Shupe非互易相移时热涨落对光传播的影响既基于折射率的温度依赖性，又依赖介质的热膨胀。但热膨胀系数随温度的变化率较折射率至少小1个量级，故可忽略。则相位大小可表示为

假设以每匝为基本单位展开，则一匝内温度涨落相同，式（18）可展开为

式中：Q为总匝数。

式（19）的关键是关于第i匝温度涨落ΔT（t，li）。可根据已建立的三维温度场分布模型，在总时间区间内进行插值，插值的最小单位为1匝的传输时间，这样可使温度分布能精确到每匝。以四极对称绕法中的第4层第13匝光纤线圈为例（物理模型如图1（a）所示），具体过程如下。

a）选定计算Shupe效应非互易相移的时间起点td，光在每匝线圈中传播的渡越时间为ti＝Δli／cm。分别以ti为时间增量点对已经建立的温度场T（r，z，t）在td后进行顺序插值，其中Δli＝lie－lis为第i匝光纤段的长度，顺次排序从光纤的一端到另一端，得到在每个tj＝ti＋td处的温度场分布。插值之前第4层第13匝光纤线圈从0.5～30s间的温度变化如图3所示。在本文建立的模型中，每匝的平均渡越时间为1.568 54ns，光波在长L的光纤中传播的渡越时间总时长为0.501 932μs。由于按匝传输的时间很短，温度变化很微小，为使表示更清楚，对图进行坐标平移和放大并取对数坐标，则Tsh＝（T－Tmin）×107，tsh＝（t－tmin）×108。此处：Tmin，ts分别为温度场温度最小值和计算时间起点。四极对称绕法光纤环按匝插值后，第4层第13匝光纤线圈分别以0.5，5.0，15.0，25.0s为时间起点，其后一段时间内（光波从光纤的一端传至另一端的总时长）的温度及其一阶导数分别如图4、5所示。

由图4、5可知：在热载荷刚影响光纤温度分布的初期（t＝0.5s时），温度上升非常快，上升趋势十分明显；当t＝5.0s时，上升速率略微降低，约10－1K／s；当t＝15.0s时，温度上升的趋势明显减慢，在光波从光纤的起始端到末端的传播时长内温度上升十分缓慢，速率约10－2K／s；当t＝25.0s时，温度趋于平稳，此时温度传递已完成，各匝的温度值基本确定。综上分析可得，模型的此变化趋势符合实际的热传递。

b）由于T（r，z，ti）是MN个时间点上的温度分布，本文需要的是特定匝在特定时间点上的温度值，用公式

进行转化变量。此处：p，q为温度场中径向和轴向的坐标序号；Tt为转化后只与变量r，z有关的温度分布场。四极对称绕法下，起始时刻为10s时，以匝渡越时间为时间增量热传递后的光纤环各匝在光波传播至此匝时的温度分布场如图6所示。

综合上述，由

求得第i匝温度涨落大小。此处：p，q为第i匝所在的位置，即第p层、第q匝。

3 经典四极绕法与交叉对称绕法实验

为验证新算法能提高非互易相移计算精度，在MATLAB软件平台中进行仿真实验。以建立的16层，每层20匝，总长度为L的光纤环温度场模型为例，四极对称和交叉对称绕法分别如图1（a）、（b）所示。分别用传统算法和按匝传输精确分析算法，计算四极对称绕法和交叉对称绕法两种绕法的Shupe效应引起的非互易相移。为更直观表现两种算法所得的非互易相移，取对数坐标，因零无对数，取近似值10－100代入。两种绕法的不同算法结果分别如图7、8所示。四极对称绕法的光纤在t＝0.5，30.0s时温度随径向的分布如图9所示。

分析上述实结果，可得以下结论。

a）对四极对称绕法，在温度开始传递的初期，四极对称绕法中各层的温度分布未能达到线性（如图9所示），因此关于光纤中点对称的两段光纤经历的温度变化并不相同，此时的Shupe效应非互易相移理论值并不为零。当温度传递结束，各层温度分布逐渐稳定且随径向呈现线性分布时（如图9所示），Shupe非互易相移理论值应为零。图7（a）中，相对理论值零，传统算法的精度可达10－11rad，而按匝传输的精确分析算法精度可达10－13rad，可见新算法精度可高3个量级。

b）对交叉对称绕法，由于温度分布的对称中心与光纤的对称中心重合。即在轴向上温度分布关于轴向中点对称，恰好对应光纤关于中点对称的两段，因此在温度场热传递初期和达到稳定后的温度分布，任何时刻交叉对称绕法下的Shupe效应非互易

相移理论值都应为零。图8中，相对理论值零，传统算法的精度最优可达10－24rad，按匝传输的精确分析算法的结果等于理论值0rad，可见新算法的精度更高，结果更准确。

c）由图7、8可知：在温度传递初期，四极对称绕法对Shupe效应非互易相移的抑制效果远不如交叉对称绕法。同一个时间点处，交叉对称绕法利用温度分布对称中心和光纤环长度对称中心的重合，大幅抑制了非互易相移，使其减小约10个量级。

4 结论

本文基于柱坐标系，用MATLAB仿真软件建立了IFOG的三维温度分布场，描述了在不同时间点上热传递过程。传统Shupe效应非互易相移计算方法是直接代入CW，CCW产生的相位差公式，涉及到温度场的时间导数等，但温度随时间的导数精度有限，且与抽样时间间隔有关。如时间间隔很大，就无法用温度的时间导数代替公式中的差值比项，从而导致非互易相移结果有偏差。新型计算方法中不需要温度场的时间导数值，仅要求获知某匝处的温度值，用两匝间的温度差代替温度的时间导数。由此本文提出了按每匝渡越时间进行温度插值的Shupe效应非互易相移计算方法，并在MATLAB仿真平台上用新型算法和传统算法分别对经典四极绕法和交叉对称绕法进行了仿真，根据某个确定的时间起点，以一匝为基本单位分别计算了两种绕法的Shupe效应产生的非互易相移。结果表明：与传统计算方法相比，按匝传输插值的Shupe效应非互易相移计算方法精度可有3个量级以上的改善。

本文的新算法提高了度量Shupe效应导致的非互易相移大小方法的精度，相较于前人文献中对光纤线圈长度分割的方法，进一步提出了时间维度上的分割。但本文研究也存在一些问题：首先，鉴于仿真条件的限制，每层光纤线圈只取了16层，与现实应用中的光纤陀螺的物理参数有差距，后续可进一步增加层数，并在真实的光纤陀螺实验环境中对其进行验证；其次，虽然本文有助于从更微观的角度了解Shupe效应原理，但未明确提出抑制的方法，后续研究可在本文算法的基础上，研究抑制Shupe效应非互易相移的方法。

［1］ 张桂才.光纤陀螺原理与技术［M］.北京：国防工业出版社，2008：156－178.

［2］ 王巍，向政，邢辉，等.光纤陀螺光纤敏感线圈三维温度瞬态响应模型［J］.中国惯性技术学报，2012，20（5）：605－601.

［3］ 李茂春，刘铁根，王佳佳，等.光纤陀螺光纤环Shupe效应的三维有限元模型［J］.天津大学学报，2008，41（5）：601－606.

［4］ 奚定平.贝塞尔函数［M］.北京：高等教育出版社，1998.

［5］ 陆全康.数学物理方法自学辅导［M］.上海：上海科学技术出版社，1989.

［6］ 俞昌铭.热传导及其数值分析［M］.北京：清华大学出版社，1981.

［7］ 因克罗普拉，德威特.传热基础［M］.北京：宇航出版社，1987.

［8］ 梁昆淼.数学物理方法［M］.北京：高等教育出版社，1998.

［9］ 方瑛，黄毅.数学物理方程与特殊函数［M］.北京：科学出版社，2012.

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

［11］ 冯文帅，王学锋，王巍.光纤环的每层匝数不同对光纤陀螺温度性能影响［J］.中国惯性技术学报，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

上海航天科技创新基金资助（SAST2015076）

成茵瑛（1991—），女，硕士生，主要研究方向为光纤陀螺的零偏误差优化。