REVIEW OF THERMAL-DYNAMICAL ANALYSIS METHODS FOR LARGE SPACE STRUCTURES
-
摘要: 近年来, 各种大型空间结构逐渐在我国航天工程中得到应用, 相应的热诱发振动问题也日益受到重视. 在此背景下, 有必要进一步梳理热诱发振动的机理和分析设计中的关键问题. 本文将结合作者的研究工作对此问题进行全面介绍, 并主要强调在分析复杂工程结构的热诱发振动问题时需要注意的特殊问题. 本文首先介绍了可以高效地分析空间薄壁杆件结构(包含开口和闭口薄壁杆件)在辐射换热条件下的瞬态温度场的Fourier有限元方法; 随后介绍了热诱发振动的线性和非线性分析方法, 强调了热−动力学耦合效应. 为了对复杂空间结构产生热诱发振动的必要条件给出解析表达式, 本文将瞬态温度场与振动位移场统一在模态空间中进行分析, 从而得到评价热诱发振动剧烈程度的一般性条件. 在此基础上, 本文进一步讨论了热诱发振动的运动稳定性问题, 以悬臂杆件的热颤振准则为例揭示了热诱发振动发散的物理机理, 并给出了评定复杂工程结构热颤振的分析方法. 论文最后概要地指出了在热诱发振动的地面试验和抑制方法中需要注意的问题, 并对将来的研究工作进行了展望.Abstract: In recent years, various large space structures are gradually implemented in the aerospace industry of China. Thus, the corresponding thermally induced vibration problems are drawn more and more attentions. Under this background, it is necessary to clarify the underling mechanism of the thermally induced vibration phenomenon and the corresponding critical issues in the analysis and design. Based on the research work of the authors, this article gives a comprehensive review of the related problems and mainly focuses on some special aspects in the thermally induced vibration analysis of complex engineering structures, which are compose of many thin-walled bars. Firstly, this article introduces a Fourier finite element that decomposes the temperature into the average part and the perturbation part. In this way, the thermal conduction equation under thermal radiation can be decoupled into the corresponding two parts due to the orthogonal property of the Fourier series. Thus, the transient temperature field of closed-section or open-section thin-walled bars can be efficiently analyzed. Based on this kind of element, both linear and nonlinear methods for the thermally induced vibration analysis are presented with the emphasis on the thermal-dynamic coupling effect. In order to give the analytical form of the necessary condition of the thermally induced vibration, this paper analyzes the properties of the transient temperature and the oscillation displacement in the mode space, and thus it obtains a general criterion to evaluate the intensity of the thermally induced vibration. Based on these work, the dynamic stability of the thermally induced vibration is further discussed by not only the mechanism reflected in the thermal flutter criterion of a cantilever bar, but also the thermal flutter analysis of complex engineering structures. Finally, the conclusion part briefly addresses some important factors in the underground testing and the method of suppressing the thermally induced responses. Some research topics need further investigating in the future are also envisaged.
-
引 言
航天器中有很多可展开附件, 比如太阳翼、天线、桅杆等等, 这些结构大多由薄壁杆件构成. 随着我国航天事业的进步, 这类结构逐渐向大型化和柔性化方向发展. 在轨展开之后, 它们在太空失重环境中就主要受到热载荷的作用, 可能产生影响结构功能的热致响应. 比如准静态热变形会改变天线的形状, 导致其聚焦能力下降. 而本文主要关注发生在低轨航天器上的热诱发振动(thermally induced vibration)现象. 因为低轨的地球半影区很窄, 所以航天器只需要几秒钟就能进出地球阴影, 从而使这类轻柔结构受到急剧变化的太阳热流作用. 对于这类低热容和低刚度的结构, 快速变化的热流会迅速改变其截面温差, 并进一步激发振动响应, 甚至产生不稳定的热诱发振动, 即热颤振(thermal flutter)现象.
热诱发振动现象在人类早期的航天活动中就已经出现. 美国的第一颗卫星Explorer I于1958年发射入轨后很快就发生了姿态偏转, 其原因就是柔性天线的热诱发振动耗散了能量[1]. 加拿大于1962年发射的第一颗卫星Alouette 1也出现了类似的事故[2]. 随着航天事业的发展, 热诱发振动现象也随着各种大型柔性结构的使用而陆续被观测到[3], 甚至出现过热颤振现象. 其中最著名的是1990年发射的哈勃太空望远镜的柔性太阳翼弯扭耦合的热颤振事故[4], 直接影响了望远镜的成像精度. 后来人们只能借助航天飞机将其太阳翼更换成半刚性结构, 才使望远镜能正常工作. 经过多次在轨维护后, 目前的哈勃太空望远镜使用的是刚性太阳翼[5]. 哈勃热颤振事故对国际航天工程产生了深远的影响[6], 对我国相关航天器的研制也具有重要的借鉴价值.
其实在人类开展航天活动之前, Boley[7]已从理论上预测了热诱发振动现象, 并引入了一个无量纲参数B来衡量发生弯曲型热诱发振动时的动态位移dmax和稳态位移dst之间的比值, 即热诱发振动的剧烈程度[8]
$$ \begin{array}{*{20}{c}} {\dfrac{{{d_{\max }}}}{{{d_{{\text{st}}}}}} = 1 + \dfrac{1}{{\sqrt {1 + {B^2}} }},}&{B = \omega \tau } \end{array} $$ (1) 式中, ω是结构的第一阶固有振动圆频率, 代表了机械振动的快慢程度; τ称为结构的热特征时间, 代表了截面温差随时间变化的快慢程度. 该公式简洁明了地解释了为什么轻柔结构容易出现热诱发振动. 参数B也因此被称为“Boley参数”[9]. 不过式(1)只适用于突加热流所产生的杆件弯曲型热诱发振动, 本文的第3节将给出复杂的空间薄壁杆系结构任意形式热诱发振动剧烈程度的度量公式, 并讨论它在渐变热流作用下的准确性.
热颤振机理的研究最早见于Augusti[10]对于一根同时受到轴力和热流作用的杆模型的分析工作. 他指出热颤振是一种热−对结构[11]强耦合效应所产生的自激振动现象. Donohue和Frisch[11]也通过讨论开口薄壁杆件的扭转振动解释了OV1-10卫星桅杆的扭转热颤振现象. Yu[12]曾试图建立悬臂杆弯曲热颤振的判断准则: 当变形前的悬臂杆指向太阳时就会发生热颤振. 但Graham[13]发现Yu[12]论文中的边界条件有误, 并将该热颤振准则修正为: 当变形前的悬臂杆背离太阳时就会发生热颤振. 该准则也被Thornton和Kim[4]用来解释哈勃太空望远镜的热颤振事故. 不过, 人们在实验中发现当热流垂直于悬臂杆入射时也会发生热颤振[3, 14], 并不符合Graham准则. 造成理论预测和实验结果不一致的原因是: 稳定性是非线性问题, 必须针对变形后的平衡状态进行分析, 而Graham热颤振准则却是基于变形前构型建立的. 于是文献[15]将Graham热颤振准则修正为: 当变形后的悬臂杆背离太阳时就会发生热颤振. 本文的第4节将对此准则进行具体的介绍, 并进一步推广至复杂工程结构的热颤振分析.
明确了物理机理后, 人们更关心如何分析复杂工程结构的热诱发振动和热颤振响应. 由于实际航天结构大多由开口和闭口薄壁杆件构成复杂的空间杆系, 需要一体地分析该复杂结构中的传热学与动力学及其互相耦合问题, 这种分析工作往往需要借助于有限元程序才能开展. Mason[16]于1968年就开始用有限元程序计算简支杆和板的热诱发振动响应. 但他的有限元模型中直接采用了温度场的解析解, 回避了辐射换热条件下的瞬态温度场分析这个高度非线性的难题. 由于考虑辐射换热的温度计算工作量很大, 用当时的计算机很难分析复杂结构的热诱发振动问题. 该状况一直持续到20世纪90年代, 特别是哈勃太空望远镜的热颤振事故之后, 人们才更加重视相关的研究. Namburu和Tamma[17]直接采用三维有限元进行热诱发振动分析, 导致计算工作量非常大. 为了提高计算效率, 为文献[18]将闭口薄壁杆的横截面温度场分解为平均温度和摄动温度两部分, 从而发展出一种谱单元来分析杆件横截面的稳态温度场, 并可求出结构的准静态热变形. 不过这种单元的平均温度和摄动温度相互耦合, 计算工作量也不小. Thornton和Kim[4]采用Ritz法分析哈勃太空望远镜太阳翼的热诱发振动问题时, 将薄壁圆管杆件横截面的温度场表示成平均温度加余弦摄动温度的形式, 简化了热传导控制方程. 而文献[19-20]进一步将任意形式(闭口或单支开口)薄壁杆件横截面的温度场展开成一般的Fourier级数形式, 发展出一种全新的一维薄壁杆件单元. 该单元的每一个节点包含一个平均温度自由度和若干摄动温度自由度. 利用Fourier级数的正交性, 可将瞬态热传导方程解耦为仅包含平均温度的非线性方程组和同时包含平均温度和摄动温度的线性方程组. 因为有限元模型中的平均温度自由度少, 因此可以很快地求解该非线性方程组, 然后再高效地求解含有摄动温度自由度的线性方程组, 从而提高了计算效率. 另外, 采用该瞬态温度单元可以与结构动力学分析的杆件单元共用节点, 因此非常有利于进行热−结构的耦合分析[20-21]. 本文的第1节将详细介绍Fourier温度单元, 并分析了平均温度和摄动温度时变特性的差别. 本文的第2节将介绍基于Fourier温度单元的热诱发振动分析方法, 主要强调在开口薄壁杆件、复杂结构中瞬态热传导−动力学耦合和几何非线性分析方面需要注意的问题.
到20世纪末, 国际上对热诱发振动的研究已经比较完整[22], 但是分析精度尚不能完全满足蓬勃发展的太空望远镜等具有高指向精度要求的航天器的需要[6,23]. 近年来, NASA还进一步开展了在轨热诱发振动实验[24], 相关成果已经应用于2021年的“双小行星重定向测试 (Double Asteroid Redirection Test, DART)”任务[25]. 我国的相关研究始于1997年的国家高技术研究发展计划. 和其他国家的研究相比, 虽然起步很晚, 但是一直立足国情致力于发展能解决复杂工程问题的高效分析方法, 并在热颤振准则方面做出了有特色的工作. 该研究工作断断续续地持续了二十多年, 最近几年才因我国航天工程的迫切需求而引起重视, 并通过地面测试验证了所发展的理论和方法[26-29], 也在工程中得到了应用. 在本文的后续小节中, 将详细地总结这个领域的研究成果, 并展望未来的工作, 以期为我国航天事业的发展做一些微薄的贡献. 限于篇幅, 本文只强调解决这种强非线性问题(辐射换热、几何大变形、热−动力学耦合)的思路, 而略去了很多具体的数学推导. 有兴趣的读者请在参考文献中寻找更多的细节.
1. 辐射换热下薄壁杆件的瞬态温度场分析
航天器可展附件通常包含由许多薄壁杆件组成的复杂构件. 它们在太空中受到太阳和地球热流照射, 同时向空间辐射散热. 研究这类结构的热诱发振动问题, 就必须发展高效的瞬态温度场分析算法, 以及能与之协调工作的结构动力学分析算法. 虽然有限元方法特别适合分析这种复杂结构, 但是传统温度杆件单元的节点自由度只有平均温度, 无法刻画杆件横截面内的温度分布, 从而无法算出对热诱发振动至关重要的热弯矩和热双力矩. 采用温度壳单元或温度实体单元虽然可以解决上述问题[17,30-31], 但是需要在薄壁杆件横截面内划分很多单元, 从而在求解辐射换热这种强非线性方程时会花费较多的计算工作量. 另外, 这些温度单元和结构单元往往不能共用一套有限元网格, 不利于分析复杂工程结构的热−动力学耦合响应. 针对这些问题, 提出了一种Fourier薄壁杆件温度有限单元, 通过增加节点摄动温度的方法来刻画杆件横截面的温度分布, 并利用Fourier级数的正交性来提高计算效率. 该温度单元可以和结构动力学分析的杆单元共用节点, 因此采用一套有限元网格就可以高效而准确地计算出复杂工程结构的热诱发振动响应.
1.1 薄壁杆件温度单元[19-20]
如图1所示, 薄壁杆件单元的几何尺寸中包含轴线长度l、横截面周长Sh和壁厚ts三个尺度. 由于壁厚远小于其他两个尺度, 所以可以忽略温度沿壁厚的变化. 于是在分析随时间t变化的温度场T时只考虑它在杆件轴线和截面周向的变化, 即
$T = T(\theta ,\xi ,t)$ , 其中$\xi = x/l$ (0 ≤ x ≤ l)和$\theta = 2{\text{π}} s/{S_h}$ (0 ≤ s ≤ Sh) 分别是单元轴向和周向的参数坐标.![]()
图 1 薄壁杆件温度单元及其单元局部坐标系Figure 1. thin-walled bar temperature element and its local coordinate system忽略薄壁杆件内壁的辐射换热, 仅考虑其外壁向外辐射散热. 根据单位长度内“体内存储的热量 + 体内传导出去的热量 + 表面辐射出去的热量 = 表面吸收的热量”可在单元局部坐标系中建立控制方程
$$ c\rho \frac{{\partial T}}{{\partial t}} - {k_\xi }\frac{{{\partial ^2}T}}{{{l^2}\partial {\xi ^2}}} - {k_s}\frac{{4{{\text{π}} ^2}{\partial ^2}T}}{{S_h^2\partial {\theta ^2}}} + \frac{{\varepsilon \sigma }}{{{t_s}}}({T^4} - T_{bk}^4) = \frac{{{\alpha _s}\bar q}}{{{t_s}}} $$ (2) 式中, c是比热容; ρ是密度; kξ和ks分别是杆件轴向和横截面周向的导热系数; αs和ε分别是杆件表面的热流吸收率和辐射率; σ = 5.67 × 10−8 W/(m2·K4)是斯蒂芬−玻尔兹曼常数; Tbk是环境背景温度. 在轨环境中的太阳热流非常均匀, 不沿杆件轴线发生变化(除非其他结构的遮挡产生了阴影, 但分析过程并不发生实质性的改变), 所以有效热流
$\bar q$ 为$$ \left.\begin{split} & \bar q\left( {\theta ,t} \right) = - \delta {\boldsymbol{q}}\left( t \right) \cdot {\boldsymbol{n}}\left( \theta \right) \\ &\delta = \left\{ {\begin{array}{*{20}{c}} {1,{\text{ }}\;\; - {\text{π}} {\text{/2 < }}\psi {\text{ < }}{\text{π}} {\text{/2}}} \\ {0,{\text{ }}\;\;{\text{π}} {\text{/2}} \leqslant \psi \leqslant {\text{3}}{\text{π}} {\text{/2}}} \end{array}} \right. \end{split}\right\}$$ (3) 式中, n是外法线方向矢量;
${\boldsymbol{q}} = - q{\boldsymbol{\nu }}$ 是入射热流矢量(方向矢量${\boldsymbol{\nu }}$ 从杆件表面向外为正方向); 吸收热流角$ \psi (\theta ) $ 为${\boldsymbol{\nu }}$ 和n之间的夹角, 它沿薄壁杆件横截面周向的变化规律取决于杆截面形状. 所以特定形状薄壁杆件所吸收的有效热流$\bar q$ 取决于各时刻入射热流${\boldsymbol{q}}(t)$ 与该杆件横截面各点的相对位置.假设单元温度沿轴向线性分布, 则可采用两节点Lagrange形函数Ni(ξ) (i = 1, 2)在轴向插值. 杆件横截面内的温度场可分解为平均温度
$\bar T$ 和摄动温度$\tilde T$ 两部分, 并将$\tilde T$ 沿周向进行Fourier级数展开$$ \begin{split} T(\theta ,\xi ,t) =& \bar T(\xi ,t) + \tilde T(\xi ,\theta ,t) \approx\\ & \sum\limits_{i = 1}^2 {\left[ {{{\bar T}_i}(t) + \sum\limits_{n = 1}^M {\tilde T_i^n(t){N^n}(\theta )} } \right]} {N_i}(\xi ) \end{split} $$ (4) 因为闭口截面的温度场在周向呈周期性分布, 而开口的断面几乎为绝热边界(即
$\partial T/\partial \theta = 0$ ), 所以上式中的周向插值函数${N^n}(\theta )$ 为$$ \tag{5a} 闭口截面 :{N}^{n}(\theta ) = \left\{ \begin{array}{l}\mathrm{sin}\left(n\theta \right),\;\;n = 1,\;3,\;5\cdots \\ \mathrm{cos}\left(n\theta \right),\;\;n = 2,\;4,\;6\cdots \end{array}\right. $$ $$ \tag{5b} 开口截面: {N^n}(\theta ) = \cos \left( {n\theta /2} \right)$$ 根据式(4), 每个单元节点包含1个平均温度和M个摄动温度, 其中摄动温度反映了沿杆件横截面的温差. 因为杆件横截面的温差(100量级)远小于其平均温度(102量级), 所以方程(2)中
$$ {T^4} \approx {\bar T^4} + 4{\tilde T^3}\tilde T $$ (6) 再注意到式(5)中的三角函数均满足正交性条件[19-21], 于是可以通过方程(2)和(3)得到将平均温度和摄动温度解耦的单元有限元方程.
单元方程经过组集后[19-21]得到的整体有限元方程仍然是解耦的, 其中的平均温度方程是非线性的
$$\boldsymbol C \dot{\bar{{\boldsymbol{T}}}} + [\bar{{\boldsymbol{K}}} + \bar{{\boldsymbol{R}}}(\bar{{\boldsymbol{T}}})] \bar{{\boldsymbol{T}}}=\bar{{\boldsymbol{Q}}} $$ (7) 根据该方程解出平均温度后, 再求解如下线性方程
$$ \boldsymbol{C} \dot{\tilde{{\boldsymbol{T}}}} + [\tilde{\boldsymbol{K}} + \tilde{\boldsymbol{R}}(\overline{\boldsymbol{T}})] \tilde{\boldsymbol{T}}=\tilde{\boldsymbol{Q}} $$ (8) 就可以得到摄动温度.
基于上述理论框架, 可以针对不同截面形状的薄壁杆件开发出具体的单元[32-34], 甚至可以考虑入射热流被部分遮挡[33]等实际工程状况. 在某些特殊的热流入射情况下, 该方法还可以推广至不等壁厚的薄壁杆[34].
1.2 热特征时间
特定变形模式的热诱发振动的剧烈程度和其对应温度随时间变化的快慢相关. 通常沿杆件轴线方向的入射热流分布是比较均匀的, 所以只需要分析杆件横截面内的平均温度和摄动温度随时间变化的特性就能揭示热诱发振动的机理. 忽略方程(2)中沿杆件轴线方向的变化项后, 再注意到式(6), 可以得到解耦的平均温度和摄动温度的齐次方程
$$\tag{9a} c\rho \frac{{\partial \bar T}}{{\partial t}} + \frac{{\varepsilon \sigma }}{{{t_s}}}{\bar T^4} = 0 $$ $$\tag{9b} c\rho \frac{{\partial \tilde T}}{{\partial t}} - \frac{{4{{\text{π}} ^2}{k_s}}}{{S_h^2}}\frac{{{\partial ^2}\tilde T}}{{\partial {\theta ^2}}} + \frac{{\varepsilon \sigma }}{{{t_s}}}{\bar T^3}\tilde T = 0 $$ 式(9)通解为
$ {T_\infty } + \left( {{T_0} - {T_\infty }} \right)\exp \left( { - t/\tau } \right) $ , 其中τ就是式(1)中反映温度变化快慢的热特征时间,$ {T_0} $ 是初始温度,${T_\infty }$ 是t→∞时的稳态温度. 根据式(9a)可估计出平均温度的热特征时间为$\bar \tau \approx c\rho {t_s}/\left( {\varepsilon \sigma {{\hat T}^3}} \right)$ , 其中$\hat T$ 的数量级和$ \bar T $ 相当. 从数量级分析可知式(9b)中${k_s}/S_h^2 \gg \varepsilon \sigma {\bar T^3}/{t_s}$ , 于是根据式(5)和式(9b)可估计出摄动温度的热特征时间为$\tilde \tau \approx c\rho S_h^2/\left( {4{{\text{π}} ^2}{k_s}} \right)$ .显然
$\tilde \tau \ll \bar \tau $ , 即摄动温度(和横截面温差相关)相对于平均温度变化得非常迅速. 根据式(1), 可见弯曲型热诱发振动主要由横截面内快速变化的温差引起; 而轴向伸缩型热诱发振动可以忽略不计.上述分析表明1.1节中利用Fourier级数正交性得到解耦的有限单元方程的方法不仅是一种数学技巧, 而且符合热传导问题的物理特性. 这种方法避免了同时求解慢变量(平均温度)和快变量(摄动温度), 既保证了求解精度又提高了求解效率. 具体来说, 方程(7) 仅包含了少量的慢变量成分, 是“浓缩”后的纯非线性部分. 该方程可以用Galerkin两点差分格式进行时间积分, 在每一时间步内用Newton-Raphson法求解. 而方程(8)虽然包含了大量的快变量成分, 但由于
${\bar{\boldsymbol T}}$ 已经从方程(7)中解出, 所以该方程是线性的, 可以高效求解.考虑到
$4{{\text{π}} ^2}{k_s}/S_h^2 \gg \varepsilon \sigma {\bar T^3}/{t_s}$ , 故方程(8)中的$ {\boldsymbol{\tilde K}} \gg {\boldsymbol{\tilde R}} $ . 再注意到平均温度变化得非常缓慢, 因此可以近似地认为辐射换热矩阵不随时间变化, 即$ {\boldsymbol{\tilde R}}\left( {{\bar{\boldsymbol T}}\left( t \right)} \right) \approx {\boldsymbol{\tilde R}}\left( {{\bar{\boldsymbol T}}\left( {{t_0}} \right)} \right) $ . 经这样处理后, 方程(8)也可以用模态叠加法求解[35-37], 即用模态空间中的前J阶模态来表达摄动温度的齐次成分${\tilde {\boldsymbol{T}}^s} = \displaystyle\sum\limits_{j = 1}^J {{{\boldsymbol{\phi}} _{\;Tj}}\;{{\rm{e}}^{ - {\tau _j}t}}}$ , 其中$ {{\boldsymbol{\phi}} _{\;Tj}} $ 和$ {{\boldsymbol{\tau }}_j} $ 分别是方程(8)的第j阶热特征主模态和热特征值(热特征时间). 对于复杂的工程结构, 可能出现很多特征值相近的密集模态, 从而给特征方程的求解造成困难. 此时可以参考Nour-Omid[36]的建议, 采用Lanczos方法[37-38]来提高求解效率.1.3 温度单元验证算例[33]
算例1 如图2(a)所示, 一根两端绝热、长度为100 mm的闭口梯形薄壁杆受到被电池片部分遮挡的太阳热流q(1)和地球热流q(2)的照射并向外层空间辐射散热. 在初始温度为290 K, 环境温度为0 K的情况下, 分别用ANSYS软件的四节点薄壳温度单元和Fourier杆件温度单元计算其温度场. 在杆件横截面周向均匀划分了21个薄壳单元, 而Fourier杆件单元中取M = 6. 如图2(b)和图2(c)所示, 在此复杂情况下Fourier杆件单元和薄壳单元的计算结果非常吻合. 另外, 同一截面上薄壳有限元模型一共有21个温度自由度, 而Fourier杆件有限元模型只需7个广义温度自由度, 从而减少了计算工作量.
算例2 如图3(a)所示, 将文献[4]中哈勃太空望远镜的伸展杆由闭口圆管改为开口薄壁圆管, 它受太阳热流q的照射并向外层空间辐射散热. 初始温度为290 K, 环境温度为0 K, 分别用ANSYS软件的薄壳温度单元(沿截面周向划分24个单元)和Fourier杆件温度单元(M = 3)计算其温度场. 图3(b)和图3(c)表明Fourier有限元的结果有足够的精度. 图3(d)给出各谐摄动温度的时变历程, 和图3(b)相比, 可以明显地看出摄动温度的变化速度远大于平均温度的变化速度, 验证了本文1.2节的分析结果.
图3(c)显示当入射热流方向与开口薄壁杆件的几何对称面不重合时, 沿杆件横截面周向温度分布不对称, 因此不仅造成弯曲变形, 还引起非对称的截面翘曲变形, 从而产生热双力矩使杆件发生扭转.
2. 热−动力学耦合有限元分析
Fourier温度杆件单元可以和常规的结构动力学杆件单元协调工作, 从而便于分析热诱发振动问题. 热−动力学耦合效应分为弱耦合和强耦合两种情况. 弱耦合只关心瞬态温度场所产生的结构动态热变形, 而强耦合还需考虑结构热变形改变了它所吸收的有效热流后对温度场造成的影响. 分析热颤振问题时, 必须考虑热−动力学强耦合效应[10].
柔性可展开附件的结构动力学分析应该考虑几何大变形的影响. 对于刚性太阳翼等比较刚硬的结构, 则完全可以采用线性方法进行分析.
2.1 热−动力学耦合有限元分析的一般框架
最一般的情况, 需同时考虑辐射换热、几何非线性和热−动力学强耦合三种非线性问题. 此时应采用迭代法逐步求解[20-21,33]. 如图4示, 设t时刻的所有场变量都为已知, 则可在全局坐标系(X, Y, Z)中求解
$ t + \Delta t $ 时刻的节点温度$ {}^{t + \Delta t}{{\boldsymbol{T}}^g} $ 和节点位移$ {}^{t + \Delta t}{{\boldsymbol{u}}^g} $ .首先用第1节的方法求解
$ {}^{t + \Delta t}{{\boldsymbol{T}}^g} $ . 为表述方便, 可认为是求解方程(7)和(8)的合并形式$$ {}^t{\boldsymbol{C}}{}^{t + \Delta t}{{\dot{\boldsymbol T}}^g} + {}^t{{\boldsymbol{K}}_T}{}^{t + \Delta t}{{\boldsymbol{T}}^g} + {}^t{\boldsymbol{R}}{}^{t + \Delta t}{{\boldsymbol{T}}^g} = {}^{t + \Delta t}{\boldsymbol{Q}}\left( {{}^{t + \Delta t}{\boldsymbol{q}},{}^t{{\boldsymbol{u}}^g}} \right) $$ (10) 上式中的热流载荷向量是先在t时刻单元局部坐标系(
$ {}^tx $ ,$ {}^ty $ ,$ {}^tz $ )中根据式(3)计算, 然后再转换到全局坐标系(X, Y, Z)中的.在假设薄壁杆件的横截面是刚周边的前提下, 定义在单元局部坐标系中的杆件外表面法向矢量n沿杆件横截面周向的分布只和杆件横截面形状有关[37-38], 不随杆件变形而改变. 但是, 由于太阳热流方向是恒定的, 定义在单元局部坐标系中的入射热流方向矢量
$ {}^t{\boldsymbol{\nu }} $ 将随着杆件的变形$ {}^t{{\boldsymbol{u}}^g} $ 而改变(具体可见2.2节).根据
$ {}^{t + \Delta t}{{\boldsymbol{T}}^g} $ 可计算结构所受的热载荷向量$ {}^{t + \Delta t}{\boldsymbol{F}} $ (包括了由平均温度产生的热轴力, 摄动温度产生的热弯矩和热双力矩), 然后通过如下更新拉格朗日格式的非线性动力学方程求解$ {}^{t + \Delta t}{{\boldsymbol{u}}^g} $ $$ {}^t{\boldsymbol{M}}{}^{t + \Delta t}{{\ddot{\boldsymbol u}}^g} + {}^t{\boldsymbol{D}}{}^{t + \Delta t}{{\dot{\boldsymbol u}}^g} + \left[ {{}^t{{{\boldsymbol{ K}}}_L} + {}^t{{\boldsymbol{K}}_{NL}}\left( {{}^t{{\boldsymbol{u}}^g},{}^t{{\boldsymbol{T}}^g}} \right)} \right]{{\boldsymbol{u}}^g}= $$ $$ {}^{t + \Delta t}{\boldsymbol{F}}\left( {{}^{t + \Delta t}{{\boldsymbol{T}}^g}} \right) - {}^t{\boldsymbol{f}}\left( {{}^t{{\boldsymbol{u}}^g},{}^t{{\boldsymbol{T}}^g}} \right) $$ (11) 其中,
$ {}^t{\boldsymbol{f}} $ 是内力向量;$ {}^t{\boldsymbol{M}} $ 和$ {}^t{\boldsymbol{D}} $ 分别是质量矩阵和阻尼矩阵. 由于热诱发振动是典型的小应变大转动问题, 所以$ {}^t{{\boldsymbol{K}}_L} $ 就是常规的线弹性刚度矩阵. 而几何刚度矩阵$ {}^t{{\boldsymbol{K}}_{NL}}\left( {{}^t{{\boldsymbol{u}}^g},{}^t{{\boldsymbol{T}}^g}} \right) $ 包含了初应力的贡献, 也考虑了温度变化的影响. 方程(11)可用Newmark法进行时间积分, 并在每一时间步内用Newton-Raphson法迭代求解.结构变形后其所受热流会发生变化, 从而导致温度改变, 并进一步影响热变形. 这种热−变形的强耦合效应在热传导方程(10)中体现为
$ {}^{t + \Delta t}{\boldsymbol{Q}} $ 和位移$ {}^t{{\boldsymbol{u}}^g} $ 相关, 动力学方程(11)中$ {}^{t + \Delta t}{\boldsymbol{F}} $ 和温度$ {}^{t + \Delta t}{{\boldsymbol{T}}^g} $ 相关. 为保证解的正确性, 需要迭代求解方程(10)和(11), 直至$ {}^{t + \Delta t}{\boldsymbol{Q}} $ 的变化量小于预设值. 迭代过程中, 这些方程中的矩阵和向量都需要根据构型的变化进行更新.在求解过程中还需要注意以下几个问题.
(1) Newton-Raphson法的切线预测步骤所用的斜率是比较灵活的, 所以在KNL中只考虑杆件轴力引起的非线性效应即已体现关键物理特征[33].
(2) Newton-Raphson法的路径修正步骤需要正确地计算不平衡力. 若直接用积分方法计算内力向量, 会因构型转换而非常繁琐. 而根据刚体转动准则[39], 考虑到Cauchy应力在随体坐标系下不会因为刚体转动而改变, 就可以直接将k−1构型单元局部坐标系中的内力
$ {}^{k - 1}{{\boldsymbol{f}}^e} $ (而不需要转换到k构型)加上表达在k构型局部坐标系中的增量内力$ \Delta {{\boldsymbol{f}}^e} $ , 得到k构型单元局部坐标系中的内力$$ {}^k{{\boldsymbol{f}}^e} = {}^{k - 1}{{\boldsymbol{f}}^e} + \Delta {{\boldsymbol{f}}^e} $$ (12) 根据文献[39], 若单元几何刚度阵
$ {}^{k - 1}{\boldsymbol{K}}_{NL}^e $ 符合刚体转动准则($ {}^{k - 1}{\boldsymbol{K}}_{NL}^e $ 乘以刚体位移等于不平衡力), 那么式(12)中的内力增量就是$ \Delta {{\boldsymbol{f}}^e} = \left( {{}^{k - 1}{{\boldsymbol{ K}}}_L^e + {}^{k - 1}{\boldsymbol{K}}_{NL}^e} \right){\boldsymbol{u}}_n^{e\left( k \right)} $ . 不过$ {\boldsymbol{u}}_n^{e\left( k \right)} $ 是去掉从k−1到k构型的位移增量$ {{\boldsymbol{u}}^{e\left( k \right)}} $ 中的刚体位移后的纯变形部分, 在实际计算中并不容易得到. 考虑到本问题的$ {}^{k - 1}{\boldsymbol{K}}_{NL}^e $ 只计及了轴力的影响, 而且通常$ {}^{k - 1}{{{ {\boldsymbol{K}}}}}_L^e \gg {}^{k - 1}{\boldsymbol{K}}_{NL}^e $ , 所以近似取$ \Delta {{\boldsymbol{f}}^e} \approx {}^{k - 1}{{{\boldsymbol K}}}_L^e{{\boldsymbol{u}}^{e\left( k \right)}} $ 也能算出很好的结果.(3) 对于开口薄壁杆需要考虑截面翘曲位移所产生的约束扭转效应, 以及杆件横截面形心和扭心之间的坐标转换关系[33-34, 40-41].
(4) 应采用2.2节介绍的Rodrigues公式来描述杆单元的空间大转动关系, 以保证求解精度.
2.2 大转动情况下的构型转换关系[33]
方程(10)和(11)中的矩阵和向量是首先在某个构型的单元局部坐标系中建立后, 再转换到全局坐标系(X, Y, Z)中的. 这些坐标转换会直接影响求解精度, 特别是热−变形强耦合分析的准确性.
如图4所示, 以位移为例, 因为初始构型单元局部坐标系(0x, 0y, 0z)是已知的, 所以可以用常规的坐标转换矩阵
${\boldsymbol{G}}_{{e}}^g$ 将初始构型中的向量$ {}^0{{\boldsymbol{u}}^e} $ 表达为全局坐标系(X, Y, Z)中的向量$ {}^0{{\boldsymbol{u}}^g} $ $$ {}^0{{\boldsymbol{u}}^g} = {\boldsymbol{G}}_{{e}}^g{}^0{{\boldsymbol{u}}^e} $$ (13) 将初始构型的单元局部坐标系(0x, 0y, 0z)进行刚体转动后就可以得到t时刻构型的单元局部坐标系(
$ {}^tx $ ,$ {}^ty $ ,$ {}^tz $ ), 从而可以用相应的正交矩阵$ {\boldsymbol{G}}_0^t $ 建立$ {}^0{{\boldsymbol{u}}^e} $ 和$ {}^t{{\boldsymbol{u}}^e} $ 之间的关系. 但是对于热诱发振动这种小应变大转动问题, 若用Euler转动公式来计算($ {}^tx $ ,$ {}^ty $ ,$ {}^tz $ )不但失去了保内积的性质, 而且所得结果还依赖于绕坐标轴的转动顺序, 因此必须使用精确的正交矩阵进行计算. 如图5所示, 此时认为(0x, 0y, 0z)是先绕正交矩阵$ {\boldsymbol{G}}_0^{\tilde t} $ 的轴方向矢量$ {}^t{{\boldsymbol{e}}_\beta } $ 旋转$ {}^t\beta $ 角度到达($ {}^tx $ ,$ {}^t\tilde y $ ,$ {}^t\tilde z $ )坐标系后, 再绕$ {}^tx $ 轴旋转$ {}^t\gamma $ 角度到达t时刻单元局部坐标系($ {}^tx $ ,$ {}^ty $ ,$ {}^tz $ ). 相应的转换关系为$$ {}^t{{\boldsymbol{u}}^e} = {\boldsymbol{G}}_0^t{}^0{{\boldsymbol{u}}^e} = {\boldsymbol{G}}_{\tilde t}^t\left( {{}^t\gamma } \right){\boldsymbol{G}}_0^{\tilde t}\left( {{}^t\beta } \right){}^0{{\boldsymbol{u}}^e} $$ (14) 式中,
$ {\boldsymbol{G}}_{\tilde t}^t\left( {{}^t\gamma } \right) $ 只涉及到从($ {}^t\tilde y $ ,$ {}^t\tilde z $ )到($ {}^ty $ ,$ {}^tz $ )的平面内转动, 可以直接用Euler转动公式; 而大转动所对应的精确正交矩阵$ {\boldsymbol{G}}_0^{\tilde t}\left( {{}^t\beta } \right) $ 为[42]$$ {\boldsymbol{G}}_0^{\tilde t} = {\left[ {{\boldsymbol{I}} + \frac{{\sin {}^t\beta }}{{{}^t\beta }}{\boldsymbol{S}}\left( {{}^t{{\beta }}} \right) + \frac{{\left( {1 - \cos {}^t\beta } \right)}}{{{}^t{\beta ^2}}}{\boldsymbol{S}}\left( {{}^t{{\beta }}} \right){\boldsymbol{S}}\left( {{}^t{{\beta }}} \right)} \right]^{{\rm{T}}} } $$ (15) 上式又称为Rodrigues公式.
$ {}^{\rm{T}}{\boldsymbol{\beta }}={\left[{}^{\rm{T}}\beta {}_{1}, {}^{\rm{T}}\beta {}_{2}, {}^{\rm{T}}\beta {}_{3}\right]}^{\text{T}} $ 是定义在(0x, 0y, 0z)中的旋转角向量, 其模为$ {}^t\beta $ , 各分量和位移增量有关.${\boldsymbol{S}}\left( {{}^t{{\beta }}} \right)$ 是反对称矩阵.式(15)是精确的, 和绕坐标轴的转动顺序无关. 当
$ {}^t\beta $ 很小时, 式(15)就退化为Euler转动公式. 根据式(13)和式(14), 就可以实现在全局坐标系和任意时刻单元局部坐标系之间的转换.2.3 热−动力学强耦合算例[21, 34]
采用上述全耦合方法分析图6所示哈勃太空望远镜的热颤振事故. 文献[4]所建立的传热学和力学模型是闭口薄壁圆杆, 只能模拟热诱发弯曲振动现象, 无法回答“热为什么会引起扭转振动”的问题. 实际上该太阳翼的Bi-Stem支承梁是开口薄壁杆件, 在轨展开后其开口方向是随机的. 为了更加接近这种实际情况, 在有限元模型[24,40-41]中采用开口C型圆杆来模拟支承梁, 并假设左右梁的开口方向有
$\theta = {20^{ \circ }}$ 的偏差. 太阳翼中部的柔性太阳毯采用等效的闭口薄壁圆杆单元网格模拟. 模型中还考虑了在轨展开后太阳毯受到的29.5 N的预张力, 以及左右梁受到的14.75 N预压力. 模拟时取太阳热流为1350 W/m2, 沿$\psi = {80^{\circ} }$ 角入射.从图7可以看出, 热诱发振动包含了弯曲(0.097 Hz)与扭转(0.0040 Hz)两种不同模式的振动. 大约在1200 s时太阳翼开始失稳, 此时梁的压缩轴力达到了压杆失稳的临界值. 图7(c)显示的失稳后的弯扭耦合变形图和实际照片非常相似. 由此证明了本文方法的有效性. 图7(a)还说明, 虽然采用线性分析可以模拟出弯扭耦合热诱发振动, 但是无法得到热颤振和热屈曲的正确结果. 可见对于柔性大型空间结构, 必须采用非线性分析.
3. 产生热诱发振动的必要条件
式(1)只适用于评价突加热流作用下纯弯曲梁的热诱发振动剧烈程度. 而实际工程结构的振动模式是非常复杂的, 所受热流也不是理想的阶跃函数. 故需讨论一般情况下产生热诱发振动的必要条件.
3.1 理论分析[37]
为了能得到类似于式(1)的解析定量结果, 本节和Boley[7-8]原始论文一样, 只采用线性结构模型且不考虑热−结构强耦合效应. 为得到一般的结果, 采用了有限元模型进行分析, 此时方程(13)退化为
$$ {\boldsymbol{M\ddot u}} + {\boldsymbol{D\dot u}} + {\boldsymbol{Ku}} = {{\boldsymbol{F}}_T} $$ (16) 式中的刚度矩阵K可以包含预应力. 热载荷
$ {{\boldsymbol{F}}_T} $ 包含了热轴力、热弯矩和热双力矩. 其中热轴力和杆件横截面平均温度有关, 是慢变量; 其他和杆件横截面的摄动温度有关, 是快变量, 可能导致热诱发振动. 为评估产生热诱发振动的必要条件及剧烈程度, 用模态叠加法求解由摄动温度引起的位移[37].根据1.1节, 用模态叠加法得到的摄动温度为
$ \tilde T = {{\boldsymbol{\varPhi }}_T}{\boldsymbol{x}}\left( t \right) $ , 其中$ {{\boldsymbol{\varPhi }}_T} $ 是包含前J阶热特征主模态$ {{\boldsymbol{\varPhi }}_{\;Tj}} $ 的矩阵; 在突加恒定热流作用下各阶热模态矢量x的第j个元素为${x_j} = x_j^*\;\left( {1 - {{\rm{e}}^{ - t/{\tau _j}}}} \right)\;$ ,$ x_j^* $ 是和热流相关的稳态解. 于是可将摄动温度所对应的热载荷表示为热模态空间的矩阵${\tilde {\boldsymbol{F}}_T} = {\boldsymbol{Wx}}\left( t \right)$ , 其中W和杆件横截面形状有关.记
$ {\boldsymbol{\tilde u}} $ 为$ {\tilde {\boldsymbol{F}}_T} $ 产生的热诱发振动位移, 它在振动模态空间表达为${\boldsymbol{\tilde u}} = \displaystyle\sum\limits_{i = 1}^I {{\boldsymbol{\phi}} _{\;i}^{\rm{T}} {{{U}}_i}}$ ,$ {{\boldsymbol{\phi}} _{\;i}} $ 和$ {U_i} $ 分别是系统的第i阶振动模态和模态坐标. 将热载荷$ {\tilde {\boldsymbol{F}}_T} $ 转换至振动模态空间, 于是方程(16)可以解耦为I个独立的方程$$ {\ddot U_i} + 2{\omega _i}{\varsigma _i}{\dot U_i} + \omega _i^2{U_i} = {\boldsymbol{\phi}} _{\;i}^{\rm{T}} {\boldsymbol{Wx}}\left( t \right) = {{\boldsymbol{\varTheta }}_i}{\boldsymbol{x}}\left( t \right) $$ (17) 式中,
$ {\omega _i} $ 和$ {\varsigma _i} $ 分别为系统的第i阶(i = 1, 2, ···, I)振动圆频率和阻尼比, 而$ {{\boldsymbol{\varTheta }}_i}{\boldsymbol{x}}\left( t \right) $ 可以表示为$$ {{\boldsymbol{\varTheta }}_i}{\boldsymbol{x}}\left( t \right) = \sum\limits_{j = 1}^J {{\varTheta _{ij}}x_j^*\left( {1 - {{\rm{e}}^{ - t/{\tau _j}}}} \right)} $$ (18) 将式(18)代入式(17), 可得
$ {U_i} = U_i^* + {\tilde U_i} $ , 其中的稳态解$ U_i^* $ 和忽略了阻尼的时变解$ {\tilde U_i} $ 为$$ U_i^* = \sum\limits_{j = 1}^J {{\varTheta _{ij}}x_j^*} /\omega _i^2 $$ (19) $$ {\tilde U_i} = \frac{{{\varTheta _{ij}}x_j^*}}{{\omega _i^2}}\left[ {1 - \frac{{\cos \left( {{\omega _i}t} \right) + {\omega _i}{\tau _j}\sin \left( {{\omega _i}t} \right)}}{{1 + {{\left( {{\omega _i}{\tau _j}} \right)}^2}}}} \right] $$ (20) 于是可知在突加恒定热流并忽略阻尼的条件下, 第i阶模态热诱发振动的剧烈程度为
$$ \dfrac{{{{\left( {{U_i}} \right)}_{\max }}}}{{U_i^*}} = \dfrac{{\displaystyle\sum\limits_{j = 1}^J {{\varTheta _{ij}}x_j^*\left[ {1 + \dfrac{1}{{\sqrt {1 + {{\left( {{\omega _i}{\tau _j}} \right)}^2}} }}} \right]} }}{{\displaystyle\sum\limits_{k = 1}^J {{\varTheta _{ik}}x_k^*} }} $$ (21) 上式说明第i阶模态热诱发振动的剧烈程度不但和
${\omega _i}$ 有关, 而且和各阶热特征时间${\tau _j}$ 都有关. 若只考虑第一阶热特征模态(J = 1), 则式(21)退化为式(1). 式(21)对具有复杂模态的工程结构仍成立, 而Boley[6-7]所给出的热诱发梁弯曲振动只是一特例.式(21)成立的条件也是非常理想的. 实际结构中会存在阻尼, 可能是非线性的, 所受热流也不是理想的阶跃函数. 此时式(21)的估计往往是保守的.
3.2 算例说明
算例1 文献[4]将哈勃太空望远镜太阳翼的Bi-Stem梁(见图6)简化为闭口薄壁圆管, 给出了其热诱发振动的理论估计. 文献[43]改变帆板伸展梁的厚度和传热参数, 用本文数值解得到不同的Boley参数与帆板端部最大位移占准静态位移之百分比的关系(见图8), 和理论预测非常吻合.
算例2 图9是一个典型的半刚性太阳翼力学模型, 主要由受18.9 N预张力的张紧网、受预压力的支承梁、顶板以及边框等组成. 整个结构由10块太阳能帆板展开而成, 每块帆板分为4个700 mm × 750 mm框. 帆板框由薄壁正方形管构成, 支承梁为薄壁圆管. 帆板框和支承梁固结于卫星舱体上. 有限元分析时将电池片质量计入由凯夫拉纤维构成的张紧网纤维上, 略去舱体运动对帆板的影响并假设二者之间绝热. 从0时刻开始, 结构受到垂直于帆板平面的太阳热流照射, 并向太空辐射散热. 太阳翼各构件的尺寸和材料见表1.
![]()
图 9 突加热流作用下的半刚性太阳翼Figure 9. A semi-stiff solar panel subjected to suddenly applied heat flux表 1 太阳能帆板各构件的尺寸及其材料Table 1. Dimension and material of the solar panelElement Dimension/mm Material frame square tube structure 1: 16 × 0.8 FRP structure 2: 8 × 0.8 boom circular tube R10 × 0.5 aluminum string rod R0.5 kevlar 文献[43]采用杆件Fourier有限元对该太阳翼模型进行离散, 总计682个节点, 1052个单元. 其中支承梁采用薄壁圆管单元, 张紧网绷弦采用实心圆杆单元(其材料参数已经等效了电池片的质量和热容量), 帆板框和顶板为薄壁矩形管单元. 分别用直接积分法和基于Lanczos方法的模态叠加法计算温度场, 所得摄动温度吻合得很好(见图10). 用Lanczos方法得到结构1边框的热特征时间为10.47 s, 而对摄动温度时程曲线以指数函数拟合得到的热特征时间约为10.5 s, 二者一致.
![]()
图 10 帆板框(结构1)和支承杆的摄动温度时间历程Figure 10. Change of perturbation temperatures of the frame and boom图11是采用弱耦合方法得到的A点(见图9)热诱发振动响应. 其中的结构1和结构2可见表1. 表2显示, 对于这个比较复杂的结构实际的热诱发振动剧烈程度比式(1)的估计低. 这是因为该结构的第一阶模态不是理想的纯弯曲梁形式(帆板是板振动模式, 和支承梁的弯曲模式不同), 热特征时间也不能完全由边框代表.
![]()
图 11 突加热流作用下太阳翼端部的法向位移Figure 11. Tip displacement of the solar panel subjected to suddenly applied heat flux表 2 A点的主要特征量 (t = 2000 s) Table 2. Main features of point A at t = 2000 sStructure 1 Structure 2 ${\tau _1}$/s 10.50 2.68 $ {\omega _1} $/(rad·s−1) 0.085 × 2π 0.052 × 2π B 5.620 0.876 $ 1 + 1/\sqrt {1 + {B^2}} $ 1.175 1.752 $ {d_{\max }}/{d_{{\text{st}}}} $ 1.089 1.180 算例3实际的入射热流总是经过一段时间后才逐渐稳定的, 不可能是阶跃函数. 因此经典的热特征时间总是小于实际温差趋于稳定的时间, 热诱发振动的剧烈程度也会小于式(1)或式(21)的估计值. 比如, 在轨航天器在进出地球阴影时, 总会经过一段半影区, 轨道越高半影区越宽, 热流变化越慢, 越不容易产生热诱发振动[3]. 另外, 柔性结构由于发生了大转动, 吸收的热流会比只发生小转动的刚性结构少一些, 热诱发振动也会小一些(没发生热颤振时). 文献[28]根据1.5 m长和2 m长之柔性悬臂梁的热诱发振动试验的测量数据, 进一步讨论了这些因素, 并根据温差变化情况拟合出等效的热特征时间τeff和相应的等效Boley参数Beff = ω1τeff. 从图12的比较结果来看, 经典理论预测的热诱发振动最剧烈, 等效理论的预测结果和实测结果接近.
![]()
4. 热诱发振动的稳定性
热诱发振动发生不稳定的物理机制是结构振动变形使结构所吸收的有效热流
$\bar q$ 因结构振动变形而震荡变化(见式(3)), 从而产生热−结构系统的运动稳定性问题. 稳定性是系统略微偏离平衡状态后表现出的性质, 所以应该在结构变形之后的构型上进行讨论. 另外, 结构发生失稳说明产生了不唯一的解, 意味着这一定是非线性问题. 但根据里雅普诺夫关于运动稳定性问题的一次近似判别准则[44], 可用系统方程在变形之后构型的一阶Taylor展开项的性质来判别是否失稳. 本节先借助悬臂杆件热颤振这个简单问题说明研究热颤振的基本方法, 然后再将它推广至复杂的工程结构.4.1 悬臂梁热颤振准则
以图13所示模型来解释引言中介绍的热颤振准则. 以半径为R的闭口薄壁圆管为例, 根据第1节的理论可得到解耦的平均温度和一阶摄动温度的控制方程
![]()
图 13 悬臂梁受突加热流作用Figure 13. A cantilever beam subjected to suddenly applied thermal flux$$ \dot{\bar{T}} + \frac{\varepsilon \sigma}{c \rho t_{s}} \bar{T}^{4}=\frac{1}{{\text{π}}} \frac{\alpha_{s} q}{\rho c t_{s}} \cos \left(\theta_{0}^{(b)} + \frac{\partial w}{\partial x}\right) $$ (22) $$ \dot{\bar{T}} + \left(\frac{k_{s}}{c \rho R^{2}} + \frac{4 \varepsilon \sigma}{c \rho t_{s}} \bar{T}^{3}\right) \tilde{T}=\frac{1}{2} \frac{\alpha_{s} q}{\rho c t_{s}} \cos \left(\theta_{0}^{(b)} + \frac{\partial w}{\partial x}\right) $$ (23) 对无阻尼的悬臂杆, 令z方向位移w的插值形式为
$ w\left( {x,t} \right) = N\left( x \right)W\left( t \right) $ , 方程(16)变为$$ M\ddot W + KW = {F_T} $$ (24) 式中,
$ {F_T} = N\left( l \right){M_T} $ ,$ {M_T} = - EI{\alpha _T}\tilde T\left( {l,t} \right)/R $ 是热弯矩, EI为梁的弯曲刚度,$ {\alpha _T} $ 是热膨胀系数.在稳态(t→∞)下, 根据式(22)和式(24)不难得到准静态挠度
$ {w_\infty } = - {\alpha _T}{\tilde T_\infty }{x^2}/\left( {2 R} \right) $ 和杆端部的转角$ {\bar \theta ^{(b)}} = - {\left. {\partial {w_\infty }/\partial x} \right|_{x = l}} = {\alpha _T}{\tilde T_\infty }l/R $ (见图13).定义杆端部的状态变量
${\boldsymbol{X}} = {\left( {W,\dot W,\tilde T} \right)^{\text{T}}}$ 和${{\boldsymbol{X}}_\infty }$ =${\left( {{W_\infty },{{\dot W}_\infty },{{\tilde T}_\infty }} \right)^{\text{T}}}$ , 则根据式(23)和式(24), 相对于稳态(热流和杆法线的夹角为$ \theta _0^{(b)} - {\bar \theta ^{(b)}} $ )的扰动量${\boldsymbol{Y}} = {\boldsymbol{X}} - {{\boldsymbol{X}}_\infty }$ 满足一阶近似关系$ {{\dot {\boldsymbol{Y}}}} \approx {{\bar {\boldsymbol{A}}}}{\boldsymbol{Y }}$ , 其中$$ {{\bar {\boldsymbol{A}}}} = \left[ {\begin{array}{*{20}{c}} 0&1&0 \\ { - \dfrac{K}{M}}&0&{ - \dfrac{{EI{\alpha _T}}}{{RM}}N\left( l \right)} \\ { - \dfrac{{{\alpha _s}q}}{{2\rho c{t_s}}}\sin \left( {{\theta _0} - \bar \theta } \right)N'\left( l \right)}&0&{ - \dfrac{{{k_s}}}{{\rho c{R^2}}} - \dfrac{{4\varepsilon \sigma }}{{\rho c{t_s}}}\bar T_\infty ^3} \end{array}} \right] $$ (25) 记
$ {\bar A_{ij}} $ 为矩阵$ {{\bar {\boldsymbol{A}}}} $ 第i行j列元素, 则$ {\bar{\boldsymbol A}} $ 的特征值$ \lambda $ 满足$ {\lambda ^3} - {\bar A_{33}}{\lambda ^2} - {\bar A_{21}}\lambda + {\bar A_{21}}{\bar A_{33}} - {\bar A_{23}}{\bar A_{31}} = 0 $ . 根据Routh-Hurwitz判据[15]可知$ {\bar A_{23}}{\bar A_{31}} > 0 $ (即$ \theta _0^{(b)} > {\bar \theta ^{(b)}} $ )时, 振动才是稳定的. 这意味着热流垂直于悬臂梁入射时是会发生热颤振的, 和实验结果[3, 14]一致.采用同样的思路, 也可以得到太阳帆的热颤振条件[45], 甚至在考虑开口方向的影响后还可以进一步得到开口杆件的热颤振准则[34, 46].
4.2 热颤振的能量解释
发生热诱发振动时, 入射热流的方向ν并未改变, 所以这是一种典型的自激振动. 当系统从外界吸收的热量超过发散的热量时, 就会发生热颤振. 因此, 也可以从能量的角度推导出悬臂梁的热颤振准则[29, 47].
围绕着稳态平衡位置, 梁的转角是振荡的
$ {\theta ^{(b)}}\left( t \right) \approx {\bar \theta ^{(b)}} + {\tilde \theta ^{(b)}}\sin \left( {\omega t} \right) $ , 相应的法向有效热流为$$ \begin{split} \bar q \approx & {q_\infty } + \Delta q = q\cos (\theta _0^{(b)} - {\bar \theta ^{(b)}}) + \\ &q\sin (\theta _0^{(b)} - {\bar \theta ^{(b)}}){\tilde \theta ^{(b)}}\sin (\omega t) \end{split} $$ (26) 将式(26)代入式(23)并注意到
$ {k_s}/{R^2} \gg 4\varepsilon \sigma {\bar T^3}/{t_s} $ , 可得$$ \tilde T \approx \tilde A{q_\infty }\tilde \tau \left( {1 - {{\rm{e}}^{ - t/\tilde \tau }}} \right) + \tilde Aq\sin \left( {\theta _0^{(b)} - {{\bar \theta }^{(b)}}} \right)\frac{{\tilde \tau {{\tilde \theta }^{(b)}}\sin \left( {\omega t - \varphi } \right)}}{{\sqrt {1 + {B^2}} }} $$ (27) 其中
$ \tilde A = {\alpha _s}/(2\rho c{t_s}) $ ,$ \tan \varphi = B = \omega \tilde \tau $ .一个周期
$T = 2{\text{π}} /\omega $ 内, 热弯矩做的功为$$ {W_T} = \int\limits_{t = T/4}^{5T/4} {{M_T}{{\dot \theta }^{(b)}}{\text{d}}t} = \frac{{B{\text{π}} {M_{T\infty }}}}{{1 + {B^2}}}{({\tilde \theta ^{(b)}})^2}\tan \left( {\theta _0^{(b)} - {{\bar \theta }^{(b)}}} \right) $$ (28) 式中
$ {M_{T\infty }} = - E{\alpha _T}{R^2}{t_s}{\tilde T_\infty } $ . 可见当$ \theta _0^{(b)} > {\bar \theta ^{(b)}} $ 时$ {W_T} < 0 $ , 不会发生热颤振.4.3 复杂结构的热颤振分析[38]
对复杂工程结构可先计算其低阶固有频率和对应的热特征时间, 按照第3节“发生热诱发振动的必要条件”, 初步判断是否会发生热诱发振动.
对于柔性结构, 需按第2节, 计算结构振动随时间的变化, 判断会否发生热诱发振动与热颤振.
若结构刚度较大, 可按有限元方程直接分析热诱发振动的稳定性区间[34]. 将式(17)写成矩阵形式
$$ {\ddot{\boldsymbol U}} + {{\boldsymbol{D}}_s}{\dot{\boldsymbol U}} + {{\boldsymbol{\varOmega}} _s}{\boldsymbol{U}} = {{\boldsymbol{\varPhi}} ^{\text{T}}}{\boldsymbol{Wx}}(t) $$ (29) 并在式(8)中考虑振动耦合项BIU的影响
$$ \boldsymbol{C} \dot{\tilde{\boldsymbol{T}}} + [\tilde{\boldsymbol{K}} + \tilde{\boldsymbol{R}}(\overline{\boldsymbol{T}})] \tilde{\boldsymbol{T}}=\tilde{\boldsymbol{Q}} + \boldsymbol{B}_{I} \boldsymbol{U} $$ (30) 定义状态变量
${\boldsymbol{ X}} = {\left( {{{\boldsymbol{U}}^{\rm{T}}},{{\dot {\boldsymbol{U}}}^{\rm{T}}},{{\tilde {\boldsymbol{T}}}^{\rm{T}}}} \right)^{\rm{T}}}$ , 并记式(29)和式(30)的稳态解(特解)为$ {{\boldsymbol{X}}_\infty } $ , 则扰动量${\boldsymbol{Y}} = {\boldsymbol{X}} - {{\boldsymbol{X}}_\infty }$ 满足一阶近似关系${\dot{\boldsymbol Y}} \approx {\boldsymbol{L}}\left( {{\boldsymbol{\nu}} ,{\boldsymbol{\varsigma}} } \right){\boldsymbol{Y}}$ , 其中$$\boldsymbol L\left( {{\boldsymbol{\nu}} ,{\boldsymbol{\varsigma}} } \right){\text{ = }}{\left[ {\begin{array}{*{20}{c}} {\boldsymbol{I}}&{\boldsymbol{0}}&{\boldsymbol{0}} \\ {\boldsymbol{0}}&{\boldsymbol{I}}&{\boldsymbol{0}} \\ {\boldsymbol{0}}&{\boldsymbol{0}}&{\boldsymbol{C}} \end{array}} \right]^{ - 1}}\left[ {\begin{array}{*{20}{c}} {{{\boldsymbol{B}}_I}}&{\boldsymbol{0}}&{ - {\boldsymbol{J}}} \\ {\boldsymbol{I}}&{\boldsymbol{0}}&{\boldsymbol{0}} \\ { - {{\boldsymbol{\varOmega}} _s}}&{ - {{\boldsymbol{D}}_s}}&{{{\boldsymbol{\varPhi}} ^{\text{T}}}{\boldsymbol{W}}} \end{array}} \right] $$ (31) 式中, J为
$ {\boldsymbol{\tilde R}}\left( {{\bar{\boldsymbol T}}} \right)\tilde {\boldsymbol{T}} $ 相对于稳态的一阶近似; ν为入射热流的方向矢量;$ \varsigma $ 为阻尼比向量.根据常微分方程理论, 当矩阵L(ν,
$ \varsigma $ )的最大特征值的实部α < 0时, 系统是稳定的; 当α > 0, 则系统不稳定; 当α = 0, 则系统处于临界状态.现以一个典型的环形天线(见图14)在突加太阳热流作用下的稳定性分析[38,48]来说明本文方法的实际应用. 环形天线由主杆、圈杆、竖杆和上下反射网构成, 其截面尺寸以及材料见表3. 天线各构件都采用薄壁圆管单元模拟, 总计748个单元. 假设天线通过主梁固结于卫星舱体, 并忽略与卫星之间的传热, 只研究环形天线本身的热−动力学问题.
表 3 环状天线的单元类型, 尺寸与材料Table 3. Dimension and material of the circular antennaDiameter × Thickness/(mm × mm) Material beam 30.0 × 5.0 FRP C rod 7.5 × 0.5 FRP V rod 15.0 × 1.0 FRP net 5.0 × 0.1 stainless steel 该天线系统的稳定性取决于三个参数: 入射热流矢量的方位角ϕ1和ϕ2和系统阻尼系数ς . 根据式(35)矩阵的最大特征值实部α (ϕ1, ϕ2, ς) = 0的条件, 可解得临界曲面, 它将三维参数空间划分为稳定区和不稳定区(见图15).
利用图15所给分区, 分别在稳定区、临界面和不稳定区计算环形天线的热诱发振动响应, 可进一步验证本节运动稳定性分析的可靠性. 图16是在给定热流方向时阻尼比和α的关系. 可见阻尼小于临界值ςcr后系统就不稳定, 此结论也可通过图17中A点(见图14)的Z向振动情况得到印证.
5. 其他实际问题
前文已经将热诱发振动和热颤振的机理、特性和分析方法做了详细的介绍, 其关键在于理解慢变量和快变量性质, 并采取相应的应对措施. 本节将进一步讨论分析复杂工程结构时需要注意的问题.
工程结构往往由很多不同种类的构件组成. 对于不同截面形式的杆状构件, 需要用不同的Fourier杆件温度单元. 另外, 还有一些构件不是杆状的, 比如太阳翼的电池板. 为了能和Fourier杆单元配合使用, 可以采用等效的梁格模型[34,49], 即用Fourier杆单元网格来等效电池板. 此时不但要求所有力学特性(比如频率)等效, 还要求热学特性(比如热特征时间)等效.
现代航天器的可展开结构不断向大型化和柔性化方向发展. 当它们的转动惯量和航天器舱体相比不再是一个小量时, 就必须考虑和舱体的耦合动力学问题[38,49]. 此时可以采用强耦合方式进行分析. 但稳定性分析的结果表明, 空间结构和航天器舱体的转动惯量比越大, 系统越容易失稳, 此时需要再考虑和姿态控制系统的协调关系[44].
工程中还特别关心怎样避免产生热诱发振动, 或者怎样抑制热诱发振动. 由于热诱发振动的频率很低(通常在0.1 Hz左右甚至更低), 被动隔振方法效果不明显, 主动控制不但难度很大而且还要付出高能耗和可靠性方面的代价[50-51]. 所以采用优化设计方法, 从源头上尽量避免发生热诱发振动, 更便于工程应用. 最简单的措施就是通过增加反射涂层来减小输入热流. 另外, 也可以根据特定的需要做更复杂的优化设计[52-53]. 当然, 最根本的解决办法是采用零膨胀材料来制造空间结构[24]. 但是在轨航天器的表面温度的变化范围在正负一百多度之间, 怎样保证在这么大的温度范围内材料都有稳定的零膨胀性能, 也是极具工程挑战性的任务.
由于热诱发振动产生的条件非常苛刻, 所以验证试验的难度非常大. 理想情况应该进行在轨试验[24], 但要付出巨大的代价. 所以通常还是开展地面试验, 即便如此其难度和花费都是很高的[26-29]. 这种试验需要在热真空罐中开展, 用红外灯阵模拟太阳热流, 并最好用抽取遮阳板的方式来模拟从阴影到光照 区的加热过程. 此外在试验过程中要特别注意以下几点.
(1) 和空间太阳热流相比, 在有限距离内红外灯阵产生的热流方向不是平行的, 其分布密度还随灯与试件的距离而变化. 因此需要考虑这些因素对试验结果的影响[3,27-29].
(2) 直射在热电偶上的热流会影响对试件表面温度测量的准确性, 从而降低了薄壁构件截面温差测量结果的精度. 因此需要精心设计热电偶的安装和温度补偿措施.
(3) 瞬态热流很难直接测量. 通常是用测量温度根据热传导方程反算出热流随时间变化的曲线. 但是由于热电偶测量的延迟效应, 反演的热流变化速度会比真实情况慢. 此时可通过热诱发振动的幅值对输入热流进行修正[28].
(4) 重力对柔性结构地面试验结果的影响很大. 但由于要求试验必须在真空中进行, 因此不能采用气浮装置消除重力. 通常可采用侧向悬挂或竖直悬挂的方式来部分解决此问题. 但是竖直悬挂试件的振动幅值会因为重力产生的恢复力而远小于真空中的幅值. 侧向悬挂受重力的影响会小一些, 但是因为悬挂装置不可能无限长, 所以也会产生额外的约束, 导致试验结果失真[28].
总之, 由于热环境和重力等各种条件的影响, 地面试验结果不能完全和在轨的真实结果吻合. 在无法开展在轨试验的情况下, 现在只能采用数值模拟和地面试验相结合的方法来预测在轨的热诱发振动响应[28]. 即, 用地面试验的结果修正有限元模型, 然后用有限元方法分析在轨真实环境下(比如失重状态)的结构响应. 因此有限元模型和试验测试的精度就至关重要.
材料参数的准确性对数值模拟结果有很大的影响[54]. 但是一些热学参数测量的分散性比较大(特别是复合材料), 而且力热参数通常是随温度变化的. 采用随机有限元方法可以考虑材料参数的不确定性, 但计算工作量很大, 此时可以采用响应面方法来提高分析效率[53]. 若要在模拟过程中考虑参数随温度的变化, 则需要增加额外的计算工作量. 比如在隐式算法中增加额外的迭代步骤, 或者直接采用时间步长较短(保证数值稳定性)的显式算法. 不过, 考虑到材料参数随温度变化的因素也是一个慢变量, 而且热诱发振动往往是在进出地球阴影的较短时间内发生的(此时平均温度还没有发生太大的变化), 那么也可以取所关心的温度范围内的参数平均值进行分析. 但是材料参数的不确定性所造成的误差是否小于空间望远镜等高精度设备的要求, 还是需要进一步讨论的.
6. 结论
经过几十年的努力, 人们虽然已经理解了热诱发振动现象的机理, 但是在工程应用时还是面临着诸多挑战[6,23]. 特别是在研制空间望远镜等对指向性要求特别苛刻的航天器时, 应该综合考虑光−热−结构−控制系统的一体化耦合分析[54]. 此时, 将多个现有软件系统集成的方式已经不能满足这种高精度的要求, 而应该发展专门的一体化数值模拟方法和与之配套的地面试验和验证技术. 这正是未来需要努力的方向.
致谢
20多年来, 清华大学先后有多位老师与博士、硕士研究生参加了本课题组的工作或给予了宝贵的建议, 他们是胡宁、任革学、岑章志、曹阳、丁勇、程乐锦、姚海民、黄彦文、李军、张晓东、张逸凡、唐羽烨、段进、李伟、范立佳、张军徽、樊孝清、袁小德和金邓格. 他们的创造性研究和辛勤劳动, 为课题的顺利完成做出了重要贡献, 作者谨在此致以衷心的感谢.
-
![]()
图 1 薄壁杆件温度单元及其单元局部坐标系
Figure 1. thin-walled bar temperature element and its local coordinate system
![]()
图 9 突加热流作用下的半刚性太阳翼
Figure 9. A semi-stiff solar panel subjected to suddenly applied heat flux
![]()
图 10 帆板框(结构1)和支承杆的摄动温度时间历程
Figure 10. Change of perturbation temperatures of the frame and boom
![]()
图 11 突加热流作用下太阳翼端部的法向位移
Figure 11. Tip displacement of the solar panel subjected to suddenly applied heat flux
![]()
![]()
图 13 悬臂梁受突加热流作用
Figure 13. A cantilever beam subjected to suddenly applied thermal flux
表 1 太阳能帆板各构件的尺寸及其材料
Table 1 Dimension and material of the solar panel
Element Dimension/mm Material frame square tube structure 1: 16 × 0.8 FRP structure 2: 8 × 0.8 boom circular tube R10 × 0.5 aluminum string rod R0.5 kevlar 
下载: 导出CSV
表 2 A点的主要特征量 (
t = 2000 s) Table 2 Main features of point A at t = 2000 s
Structure 1 Structure 2 ${\tau _1}$/s 10.50 2.68 $ {\omega _1} $/(rad·s−1) 0.085 × 2π 0.052 × 2π B 5.620 0.876 $ 1 + 1/\sqrt {1 + {B^2}} $ 1.175 1.752 $ {d_{\max }}/{d_{{\text{st}}}} $ 1.089 1.180
下载: 导出CSV
表 3 环状天线的单元类型, 尺寸与材料
Table 3 Dimension and material of the circular antenna
Diameter × Thickness/(mm × mm) Material beam 30.0 × 5.0 FRP C rod 7.5 × 0.5 FRP V rod 15.0 × 1.0 FRP net 5.0 × 0.1 stainless steel
下载: 导出CSV
-
[1] Thomson WT, Reiter GS. Attitude drift of space vehicles. The Journal of the Astronautical Sciences, 1960, 7: 29-36
[2] Etkin B, Hughes PC. Explanation of the anomalous spin behavior of satellites with long flexible antennae. Journal of Spacecraft and Rockets, 1967, 9(4): 1139-1145
[3] Foster RS. Thermally induced vibrations of spacecraft booms. [PhD Thesis]. Virginia, US: University of Virginia, 1998
[4] Thornton EA, Kim YA. Thermally induced bending vibrations of a flexible rolled-up solar array. Journal of Spacecraft and Rockets, 1993, 30(4): 438-448 doi: 10.2514/3.25550
[5] Hawley SA. Hubble space telescope solar array concerns and consequence for service mission 2. Journal of Spacecraft and Rockets, 2016, 53(1): 15-24
[6] Corbacho VV, Kuiper H, Gill E. Review on thermal and mechanical challenges in the development of deployable space optics. Journal of Astronomical Telescopes, Instruments, and Systems, 2020, 6(1): 010902
[7] Boley BA. Thermally induced vibrations of bars. Journal of the Aeronautical Sciences, 1956, 23(2): 179-181
[8] Boley BA. Approximate analyses of thermally induced vibrations of bars and plates. Journal of Applied Mechanics, 1972, 39(3): 212-216
[9] Boley BA. A communication quasi una fantasìa. Mechanics Research Communications, 2015, 68(9): 2-4
[10] Augusti G. Instability of struts subject to radiant heat. Meccanica, 1968, 3(9): 167-176
[11] Donohue JH, Frisch HP. Thermoelastic instability of open-section booms. NASA Technical Note, D-5310, 1969
[12] Yu YY. Thermally induced vibration and flutter of a flexible boom. Journal of Spacecraft and Rockets, 1969, 6(8): 902-910 doi: 10.2514/3.29725
[13] Graham JD. Solar induced bending vibrations of a flexible member. AIAA Journal, 1970, 8(11): 2031-2036 doi: 10.2514/3.6042
[14] Rimrott F, Abdel-Sayed R. Flexural thermal flutter under laboratory conditions. Transactions of the Canadian Society for Mechanical Engineering, 1976, 4(4): 189-196 doi: 10.1139/tcsme-1976-0027
[15] Zhang JH, Xiang ZH, Liu YH, et al. Stability of thermally induced vibration of a bar subjected to solar heating. AIAA Journal, 2014, 52(3): 660-665 doi: 10.2514/1.J052574
[16] Mason JB. Analysis of thermally induced structural vibrations by finite element techniques. NASA Technical Memorandum, X-321-68-333, 1968
[17] Namburu RR, Tamma KK. Thermally-induced structural dynamic response of flexural configurations influenced by linear/non-linear thermal effects. AIAA Paper, 91-1175, 1991
[18] Givoli D, Rand O. Harmonic finite element thermo-elastic analysis of space frames and trusses. Journal of Thermal Stresses, 1993, 16(3): 233-248 doi: 10.1080/01495739308946228
[19] Xue MD, Ding Y. Two kinds of tube elements for transient thermal–structural analysis of large space structures. International Journal for Numerical Methods in Engineering, 2004, 59(10): 1335-1353 doi: 10.1002/nme.918
[20] Xue MD, Duan J, Xiang ZH. Thermally-induced bending-torsion coupling vibration of large scale space structures. Computational Mechanics, 2007, 40(4): 707-723 doi: 10.1007/s00466-006-0134-x
[21] Duan J, Xiang ZH, Xue MD. Thermal-dynamic coupling analysis of large space structures considering geometric nonlinearity. International Journal of Structural Stability and Dynamics, 2008, 8(4): 569-596 doi: 10.1142/S0219455408002806
[22] Thornton EA. Thermal Structures for Aerospace Applications. Virginia, US: AIAA Education Series, 1996
[23] 胡斌, 李创, 相萌等. 可展开空间光学望远镜技术发展及展望. 红外与激光工程, 2021, 50(11): 20210199 Hu Bin, Li Chuang, Xiang Meng, et al. Development and prospects of deployable space optical telescope technology. Infrared and Laser Engineering, 2021, 50(11): 20210199(in Chinese))
[24] Chamberlain MK, Kiefer SH, LaPointe M, et al. On-orbit flight testing of the Roll-Out Solar Array. Acta Astronautica, 2021, 179(2): 407-414
[25] [26] Su XM, Zhang JH, Wang J, et al. Experimental investigation of the thermally-induced vibration of a space boom section. Science China Physics, Mechanics & Astronomy, 2015, 58(4): 044601
[27] Fan C, Bi YQ, Wang J, et al. Experimental investigation of heat flux characteristics on the thermally induced vibration of a slender thin-walled bar. International Journal of Applied Mechanics, 2020, 12(5): 2050053 doi: 10.1142/S1758825120500532
[28] Wang J, Jin DG, Fan C, et al. Predicting the on-orbit thermally induced vibration through the integrated numerical and experimental approach. Acta Astronautica, 2022, 192(3): 341-350
[29] Jin DG, Fan C, Wang J, et al. Experimental verification of the thermal flutter criterion for a slender cantilever boom. AIAA Journal, 2022, in press
[30] Liu L, Sun SP, Cao DQ, et al. Thermal-structural analysis for flexible spacecraft with single or double solar panels: A comparison study. Acta Astronautica, 2019, 154(1): 33-43
[31] Cao YT, Cao DQ, He GQ, et al. Thermal alternation induced vibration analysis of spacecraft with lateral solar arrays in orbit. Applied Mathematical Modelling, 2020, 86(10): 166-184
[32] 丁勇. 大型空间结构的热-结构有限元分析. [博士论文]. 北京: 清华大学, 2002 Ding Yong. Thermal-structural analysis of large space structure with finite element method. [PhD Thesis]. Beijing: Tsinghua University, 2002 (in Chinese))
[33] 段进. 大型柔性空间结构的热-动力学耦合有限元分析. [博士论文]. 北京: 清华大学, 2007 Duan Jin. The thermal-dynamic coupling analysis of large flexible space structures by finite element method. [PhD Thesis]. Beijing: Tsinghua University, 2007 (in Chinese))
[34] 袁小德. 两种热致响应分析新单元及弯扭耦合热颤振准则. [硕士论文]. 北京: 清华大学, 2019 Yuan Xiaode. Two new elements for thermally-induced responses and the bending and torsion coupled thermal-flutter criterion. [Master Thesis]. Beijing: Tsinghua University, 2019 (in Chinese))
[35] Huebner KH, Donald LD, Douglas ES, et al. Finite Element Method for Engineers, 4th ed. New York: John Wiley & Sons Inc, 2001
[36] Nour-Omid B. Lanczos method for heat conduction analysis. International Journal for Numerical Methods in Engineering, 1987, 24(1): 251-262 doi: 10.1002/nme.1620240117
[37] 程乐锦. 大型空间结构的热诱发振动有限元分析. [博士论文]. 北京: 清华大学, 2003 Cheng Lejin. Finite element analysis for thermally induced vibrations of large space structures. [PhD Thesis]. Beijing: Tsinghua University, 2003 (in Chinese))
[38] 李伟. 卫星刚体-结构附件耦合系统热-动力学有限元分析. [博士论文]. 北京: 清华大学, 2007 Li Wei. The thermal-dynamic analysis of satellites with flexible appendages by finite element method. [PhD Thesis]. Beijing: Tsinghua University, 2007 (in Chinese))
[39] Yang YB, Chiou HT. Rigid body motion test for nonlinear analysis with bar elements. Journal of Engineering Mechanics, 1987, 113(9): 1404-1419 doi: 10.1061/(ASCE)0733-9399(1987)113:9(1404)
[40] 黄彦文. 含开口薄壁杆件的大型空间结构热诱发弯扭耦合振动有限元分析. [硕士论文]. 北京: 清华大学, 2004 Huang Yanwen. Thermally induced vibration analysis of large space structures including thin-walled open section beam by FEM. [Master Thesis]. Beijing: Tsinghua University, 2004 (in Chinese))
[41] 黄彦文, 薛明德, 程乐锦等. 含开口薄壁杆的大型空间结构热诱发弯扭振动. 清华大学学报(自然科学版), 2005, 45(2): 262-266 doi: 10.3321/j.issn:1000-0054.2005.02.032 Huang Yanwen, Xue Mingde, Cheng Lejin, et al. Thermally induced vibrations of large space structures including thin-walled open bar sections. Journal of Tsinghua University (Sci and Tech), 2005, 45(2): 262-266(in Chinese)) doi: 10.3321/j.issn:1000-0054.2005.02.032
[42] 黄克智, 薛明德, 陆明万. 张量分析, 第3版. 北京: 清华大学出版社, 2020: 71 Huang Kezhi, Xue Mingde, Lu Mingwan. Tensor Analysis, 3rd ed. Beijing: Tsinghua University Press, 2020: 71 (in Chinese)
[43] 程乐锦, 薛明德, 唐羽烨等. 大型空间结构的热-结构动力学分析. 应用力学学报, 2004, 21(2): 1-9 doi: 10.3969/j.issn.1000-4939.2004.02.001 Cheng Lejin, Xue Mingde, Tang Yuye, et al. Thermal-dynamic analysis of large scale space structures by FEM. Chinese Journal of Applied Mechanics, 2004, 21(2): 1-9(in Chinese)) doi: 10.3969/j.issn.1000-4939.2004.02.001
[44] 中国大百科全书, 力学卷, 运动稳定性. 北京: 中国大百科全书出版社, 1985: 570-573 Encyclopedia of China, Mechanics, Stability of Motion. Beijing: Encyclopedia of China Publishing House, 1985: 570-573 (in Chinese))
[45] Zhang JH, Wang PH, Liu YF, et al. Can boom-supported solar sails flutter? AIAA Journal, 2020, 58(10): 4600-4603
[46] Yuan XD, Xiang ZH. A thermal-flutter criterion for an open thin-walled circular cantilever bar subject to solar heating. Chinese Journal of Aeronautics, 2018, 31(9): 1902-1909 doi: 10.1016/j.cja.2018.07.002
[47] Shen Z, Hu G. Thermoelastic−structural analysis of space Thin-Walled beam under solar flux. AIAA Journal, 2019, 57(4): 1781-1785 doi: 10.2514/1.J057793
[48] Li W, Xiang ZH, Chen L, et al. Thermal flutter analysis of large-scale space structures based on finite element method. International Journal for Numerical Methods in Engineering, 2007, 69(5): 887-907 doi: 10.1002/nme.1793
[49] 樊孝清. 舱体-挠性附件系统的热诱发振动分析与控制. [硕士论文]. 北京: 清华大学, 2016 Fan Xiaoqing. Thermally induced vibration analysis and control of the rigid hub-flexible attachment system. [Master Thesis]. Beijing: Tsinghua University, 2016 (in Chinese)
[50] Zhang JH, Xiang ZH, Liu YH. Control of the thermally induced vibration of space structures by using heaters. Journal of Spacecraft and Rockets, 2014, 51(5): 1454-1463 doi: 10.2514/1.A32601
[51] Liu CC, Jing XJ, Daley S, et al. Recent advances in micro-vibration isolation. Mechanical Systems and Signal Processing, 2015, 56-57: 55-80 doi: 10.1016/j.ymssp.2014.10.007
[52] Fan LJ, Xiang ZH, Xue MD, et al. Robust optimization of thermal-dynamic coupling systems using a kriging model. Journal of Spacecraft and Rockets, 2010, 47(6): 1029-1037 doi: 10.2514/1.49307
[53] Fan LJ, Xiang ZH. Suppressing the thermally induced vibration of large-scale space structures via structural optimization. Journal of Thermal Stresses, 2015, 38(1): 1-21 doi: 10.1080/01495739.2014.950529
[54] Levine M, Fanson J. Advanced thermo-structural technologies for the NASA terrestrial planet finder mission. Structural Control and Health Monitoring, 2006, 13(1): 190-209 doi: 10.1002/stc.136
-
期刊类型引用(2)
1. 刘华雩，郑永彤，彭海峰，高效伟，杨恺. 基于径向积分边界元法的周期性复合材料面内等效热弹性参数分析方法. 力学学报. 2024(08): 2294-2303 . 
本站查看
2. 吴磊，韩飞，邓子辰. 计及轴向力作用的柔性梁热——形耦合振动研究. 航空工程进展. 2024(06): 216-223+254 . 百度学术
其他类型引用(2)
计量
- 文章访问数: 14
- HTML全文浏览量: 0
- PDF下载量: 2
- 被引次数: 4
