拉格朗日插值方法和等参单元逆变换方法在热固耦合中的应用

金赛英，陈铁锋

软件应用

拉格朗日插值方法和等参单元逆变换方法在热固耦合中的应用

金赛英，陈铁锋

（中国航发商用航空发动机有限责任公司，上海200241）

基于拉格朗日插值方法，结合等参单元逆变换方法，将目标插值节点的全局坐标变化成等参单元的局部坐标，获得拉格朗日插值所需的形函数，从而得到目标插值节点所对应的热固耦合结果。对该方法进行了验证并程序化，适用于航空发动机二维整机有限元分析模型进行了热固耦合，算例表明该方法具有较高的精度，工程适用性好，能用于类似的多种航空发动机结构进行热固耦合。

拉格朗日插值方法；等参单元逆变换；航空发动机；热固耦合

在工程设计上一般需要热固耦合加载温度载荷，目前常用的热固耦合方法采用就近插值方法或者商用软件自带的插值方法，在解决温度梯度大、热模型和结构模型不一致时存在工程误差大的特点。针对航空发动复杂精细的结构特征，常规的热固耦合方法，无法满足篦齿和轮筒等截面变化梯度比较大的复杂位置的温度载荷精度要求，也无法满足获得精细结构轮廓的冷热态换算计算的温度载荷精度要求。目前国内在耦合方法方面开展了大量的研究工作，但关于热固耦合的两种物理结构网格上的参数转换问题，主要针对插值方法以及单元间结构参数的有效变换等问题开展了研究[1-5]，Mutri V and Valliappan S和钱向东等人提出了基于泰勒展开的等参单元逆变换方法[6-8]，朱以文等开展了比较精细的等参单元换算方法研究[9-10]，其提出的等参有限元逆变换高效算法工程适用性非常好。本文对某型号二维整机有限元模型的温度场进行了研究，基于拉格朗日-权函数插值方法[11]，采用等参单元逆变换思想，将目标插值节点的全局坐标变成等参单元的局部坐标，从而获得插值形函数，最终得到目标插值节点的温度载荷。插值结果表明基于拉格朗日-权函数插值方法和采用等参单元逆变换方法得到的转子温度分布准确度高，工程适用性好。航空发动机大部分结构都需要承受温度载荷，本方法能广泛应用于转子、机匣、轴类等航空发动机结构，并已形成了相应的工程软件。

1 拉格朗日插值方法及形函数

1.1 拉格朗日插值方法

采用拉格朗日形函数插值方法，需要考虑单元类型以及阶次的组合以及不同单元类型和阶次组合下的形函数类别。对于热固耦合问题，基于热模型的单元、节点布局，获得所有结构模型节点的相应形参，从而得到相应所有目标节点的插值载荷。如图1圆块为热模型，方块为结构模型，结构模型上的一个待插值节点p处在热模型的单元ijk内。如公式（1）所示Li、Lj、Lk为目标节点p在三角形单元△ijk内对应的形函数，LOAD（i）、LOAD（j）、LOAD（k）为相应的温度模型单元节点的温度。

图1 热固耦合模型

形函数是定义于单元内部的、坐标的连续函数，它应满足一下条件：

a）在节点i，形函数Ni=1；在其他节点，Ni=0；

b）能保证用它定义的未知量（u，v或x，y）在相邻单元之间的连续性；

c）它包含任意线性项，以便用它定义的单元位移可满足常应变条件；

d）应满足下列等式∑Ni.

因此（1）式中，Li+Lj+Lk=1.

1.2 单元形函数

常用的二维热模型和二维结构模型的离散有限元网格有一次三角形单元、二次三角形单元、一次四边形单元和二次四边形单元。如图2所示为四种单元类型的母单元。

图2 二维母单元（从左至右分别为：一次三角形单元，二次三角形单元，一次四边形单元，二次四边形单元）

a）对于一次三角形单元（3节点），基于面积计算得到形函数

本文研究的热固耦合模型基于二维结构，考虑单元类型所要用到形函数为以上四种。对于三角形单元，只要定位目标插值节点在热模型上的精确位置，就可通过单元内的形函数积分获得其插值结果。而对于四边形单元，除了定位其在热模型上的精确位置，还需获得在等参单元对应的ξ和η值。通常的四边形单元不是规则的母单元，需要等参变换为规则的母单元，每个母单元都具有一个独立的局部坐标系，这里目标插值节点只有全局坐标系下的坐标值，基于拉格朗日形函数法，需要目标插值接节点所在热模型单元对应的母单元的局部坐标系下的坐标值。因此要开展等参单元逆变换，获得目标插值节点的局部坐标值，即ξ和η.

2 等参单元逆变换方法

如图3所示为四边形等参单元的全局坐标系在全局坐标系和局部坐标系下的示意图，在四个节点的局部坐标分别为（-1，-1），（1，-1），（1，1），（-1，1），在局部坐标系下单元内任意位置的坐标为满足

图3 四边形等参单元坐标模型

在全局坐标系下（x为横坐标，y为纵坐标），单元内部任意一点的坐标可以用如（8）式表达：

式中，x0和y0为单元内任意位置在全局坐标系下的坐标值，ξ和η为对应的在局部坐标系下的坐标值。

本文所采用的等参单元逆变换，将结构模型的目标插值节点映射到热模型相应的单元中，获得相应的局部坐标值。如（8）式所示，已知xi和yi（i =0，1，2，3，4），求ξ和η值。理论上通过一组联立方程，可以数值求解两个未知量。但是之前的研究表明，这种数值求解存在较大的误差，特别是在单元交界位置存在无解的情况。因此必须建立一种有效的等参单元逆变换的方法，在保证求解精度的同时能完成所有目标节点的变换。本文所采用的变换方法如下公式所示：

3 算法验证及程序

3.1 拉格朗日算法验证

为了验证拉格朗日算法的准确性，建立了用于验证的热模型和有限元模型。如图4所示，热分析模型耦合了方块状有限元模型的网格，在热分析过程中将有限元模型对应的网格节点温度也获得。输出温度时将有限元部分的网格温度信息不输出，只输出如图4所示的圆块状分布。根据如图4所输出的圆块状温度，需要通过插值分析获得如图5所示的方块状有限元模型的温度。

图4 热分析模型

图5 有限元分析模型

从图6可见，方块状的温度分布在与热模型叠加时的温度情况，就近插值方法不能复原热模型的温度场分布规律，而本方法能达到与初始热模型完全吻合的程度。

图6 热固耦合结果

通过对照有限元模型各个节点的温度结果，如表1所示。在节点3、4、5、13、14、17、18、19等节点处的误差最大，这些节点的共同特性是位于单元内部，离最近的热模型节点较远。采用拉格朗日插值方法能克服网格节点之间的局限性，通过形函数将目标节点与热模型的有限元模型连接起来，作为一种同等的离散节点，从而实现目标节点的精确插值，本方法插值误差最高为6.4%，绝大多数目标耦合点能达到零误差。常用的就近插值方法有7个点的插值误差超过10%，在温度梯度大的点高达53.8%.

表1 热固耦合插值结果

3.2 程序算法

航空发动机结构复杂，插值方法需要很强大的通用性，因此在处理算法的时候应该考虑到各个算法模块化，同时热固耦合的输入输出规范化。

具体的程序方法实现如图7所示。

图7 程序流程图

4 算例

采用基于拉格朗日和等参单元逆变换方法及其完善的耦合程序应用于某型二维整机有限元模型温度场加载，如图8为高压涡轮一级盘盘心位置热固耦合情况，图9所示为低压涡轮一级盘轮缘及鼓筒位置热固耦合情况，图10为涡轮部件的整个部件热固耦合结果。

通过对某型二维整机有限元模型的表明，本方法插值结果精确较高，程序稳定性好，可以用于解决整机的转子、轮盘和轴类的热固耦合问题。

图8 高压涡轮一级盘盘心拉格朗日插值热固耦合

图9 低压涡轮一级盘轮缘及鼓筒处拉格朗日插值热固耦合

图10 高低压涡轮转静子有限元热固耦合结果

5 结束语

采用拉格朗日形函数思想，结合等参单元逆变换方法，实现了在热模型离散单元内部插值，结论如下：

（1）本文针对热固耦合的离散单元模型，将结构节点理解为热模型离散单元连续体的某个点，通过采用拉格朗日形函数方法，将所有的结构节点归一化到所对应的热模型的某个单元内。同时，采用等参单元逆变换方法，将目标节点从全局坐标系变换到局部坐标系下，最终通过单元内插值，得到所有目标插值节点的热耦合值。

（2）针对热模型为一次以及二次的三角形单元和四边形单元混合的离散模型，采用本方法形成了通用的工程软件，并成功用于某型二维整机有限元模型的温度场插值，插值结果非常精确，具有很强的工程适用性，可以广泛应用于类似结构的强度、寿命评估，以及需要获得精确结构冷、热态结构的冷热态换算。

（3）目前考虑的单元类型为二维平面内的单元，基于本文提出的方法，完善三维形函数库和三维空间的等参单元逆变换，可用于空间壳体模型（如叶片表面单元）以及三维模型的热固耦合问题，并扩展至空间流固耦合问题。

[1]刘振宇，傅云，谭建荣，基于异构网格耦合的产品多物理场有限元数据集成与可视化仿真[J].机械工程学报，2010，46（7）：114-121.

[2]方锡武，刘振宇，谭建荣，等.异构有限元网格多场信息的等效集成方法[J].浙江大学学报，2014，48（6）：301-302

[3]沈毅，三维曲面插值在发动机涡轮叶片有限元温度场计算中的应用[J].航空发动机，2002（2）：43－45.

[4]包劲青，杨强，官福海.线性等参单元逆变换的二分法及其修正[J].计算力学学报，2010，27（5）：770-774.

[5]李春光，郑宏，葛修润，等.六面体单元等参逆变换的一种迭代解法[J].岩土力学，2004，25（7）：1050-1052.

[6]Mutri V and Valliappan S numerical inverse isoparametric mapping in remeshing and nodal quantity continuing[J].com puter&structure，1986，22：1011－1012.

[7]Mutri V，Wang Y and and Valliappan S numerical inverse isoparametric mapping in 3D FEM[J].computer&structure，1988，29：611－622.

[8]钱向东，任文青.一种高效的等参有限元逆变换算法[J].计算力学学报，1998（4）：437-441.

[9]朱以文，李伟，蔡元奇.基于解析性质的等参有限元逆变换高效算法[J].武汉大学学报，2002，35（2）：62-65.

[10]徐燕萍，项阳，刘泉声，等.参逆变换插值法的研究及其应用[J].岩土力学，2001，22（2）：226-228.

[11]王勖成.有限单元法[M].北京：清华大学出版社，2003：468－472.

Lagrange Interpolation Method and Isoparametric Element Inverse Transformation Method Application in Thermosetting Coupling

JIN Sai－ying，CHEN Tie－feng
（Chinese Hangfa Commercial Aircraft Engine Co.，Ltd.，Shanghai 200241，China）

Based on the lagrangian interpolation method and inverse is oparametric element mapping method，to get the local coordinate of a target point in the aera of a isoparametric element by changing it’s global coordinate，then obtain the desired parameters of the lagrangian interpolation shape function，so the corresponding thermalstruture coupling result of the target mapping point be arrived.In this paper，the above method be validated and be programed，and applicable to thermal-structure coupling of the two-dimensional aero-engine rotor structure，the result show that this method has high accuracy and good adaptability in engineering utilize，it can be used in a variety of similar structure of aero-engine to thermal-structure coupling.

lagrangian interpolation method；inverse isoparametric element mapping；aero-engine thermal-structure coupling

TP301.6

B

1672-545X（2017）07-0232-05

2017-04-09

金赛英（1986-），女，江西吉安人，工学硕士，工程师，研究方向：航空发动机结构强度。