（工程力学专业论文）基于振动响应信号的复合材料层合壳体结构遗传算法损伤检测.pdf

摘要 摘要 本文研究了使用遗传算法对含裂纹悬臂复合材料层合圆柱壳体 进行损伤振动检测的问题。 首先，建立表面粘贴有压电片材料的悬臂完好复合材料层合壳体 及损伤复合材料层合壳体的有限元动力学模型，并求出其动力响应。 其次使用相关系数来判定各种损伤情况响应之间的相似程度。 并通过使用遗传算法来搜索各种损伤情况，找到其响应与待识别结构 响应一致的损伤情况；从而实现损伤识别。 最后，对使用遗传算法进行损伤识别时的所需计算量进行了分 析。分析结果表明遗传算法能有效进行损伤识别。而且遗传算法方法 所具备的一些特点使其在对大型复杂结构和多损伤情况进行损伤识 别时将体现出显著的优势。并提出，通过分析结构损伤识别所用遗传 算法迭代次数等数值的分布规律，可以在对新结构进行损伤识别时 较准确迅速地预计所需的计算量，并依此改进结构损伤识别计算的效 率。 本文的工作对于进一步研究使用遗传算法进行损伤识别，具有参 考意义。 关键词：复合材料层台壳体，损伤检测，遗传算法 a b s t r a c t a b s t r a c t i nt h i sp a p e r , o n l i n ed e t e c t i o no f c r a c kd a m a g ef o rac a n t i l e v e rl a m i n a t e dc o m p o s i t e c y l i n d r i c a ls h e l li ss t u d i e db yu s i n g g e n e t i c a l g o r i t h m s ( g a s ) f i r s t ，t h e f i n i t ee l e m e n td y n a m i cm o d e l so ft h ec a n t i l e v e ri n t a c ta n dd a m a g e d l a m i n a t e d c o m p o s i t ec y l i n d r i c a l s h e l l sa r ee s t a b l i s h e d t h e p i e z o e l e c t r i c s m a r t m a t e r i a l sb o n d e do nt h el a m i n a t e dc o m p o s i t es h e l la r cu s e da sb o t ha c t u a t o ra n d s e n s o r i n t r o d u c i n g c r a c kd a m a g ei n t os h e l lf e mm o d e l ，v a r i a t i o n so fp h y s i c a l p a r a m e t e r s d u et oc r a c ka r ed i s c u s s e d t h ed y n a m i cr e s p o n s es i g n a l sa r ec a l c u l a t e d s e c o n d ，b yt a k i n g c o r r e l a t i o nc o e f f i c i e n t sa st h ec r i t e r i o nt h a tm e a s u r e st h e c o r r e l a t i o nb e t w e e nd y n a m i cr e s p o n s es i g n a lo fi n t a c ts h e l la n dd a m a g e ds h e l l s ，g a s i su s e dt of i n dt h ed a m a g ec a s et h a th a st h es a m ed y n a m i cr e s p o n s es i g n a lw i 也t h e n e e d i n gi d e n t i f i c a t i o nd a m a g ec a s e w h e ns u c hac a s ei sf o u n d ，t h ed a m a g e d e t e c t i o n i sa c h i e v e db yt a k i n gt h i sc a s ea st h ed e t e c t i o nr e s u l t a tl a s t ，t h ec o m p u t ew o r k l o a do fc r a c kd a m a g ed e t e c t i o nf o rac a n t i l e v e rl a m i n a t e d c o m p o s i t ec y l i n d r i c a ls h e l lb yu s i n gg a s i sa n a l y z e d t h er e s u l ts h o w st h a tg a sc a l l g e tt h ed a m a g e i n f o r m a t i o n e f f i c i e n t l y t h ea n a l y s i sa l s os h o w s t h a ts o m ec h a r a c t e ro f g a sw i l lm a k ei t s u p e r bi n t h ep r o b l e mo fd a m a g ed e t e c t i o no fl a r g es t r u c t u r eo r m u l t i - d a m a g ed e t e c t i o n m o r e o v e r , t h ed i s t r i b u t i o nf u n c t i o n so f t h eg a sc h a r a c t e r s s u c ha si t e r a t e dt i m ei sf o u n d f u r t h e r m o r e ，t h ee x p e c t e dc o m p u t ew o r k l o a do f d a m a g e d e t e c t i o ni nan e ws t r u c t u r ec a nb eq u i c k l ya c q u i r e du s i n gt h i sf u n c t i o n t h er e s u l t sp r e s e n t e di nt h i sp a p e ra r eu s e f u lf o rt h ef u t u r es t u d yo fu s i n gg a si n d a m a g e d e t e c t i o n k e y w o r k s ：l a m i n a t e dc o m p o s i t es h e l l ，c r a c kd a m a g ed e t e c t i o n ，g e n e t i ca l g o r i t h m s i i 第一章绪论 1 1 结构损伤检测概述 2 0 世纪9 0 年代以来，随着大跨度铁路公路两用桥、大型水利发 电机组、大型核能发电机组、国际空间站计划、火星无人探测器计 划等大型工程项耳的实施，大型结构体系的健康自动诊断问题迅速 成为了国际学术界和工程界关注的热点问题。以国际空间站的建造 为契机，美、英、法、德、日等发达国家相继投入了大量的研究与 发展资金用于支持该领域的探索性研究。我国自9 0 年代中期以来， 也先后在一些航天高科技计划和重大项目中投入了一定的资金用于 支持大型结构体系健康自动诊断问题的探索性研究。 调查研究发现，导致复杂结构体系发生故障和失效的原因是多 方面的，其中，循环载荷作用下的疲劳损伤累积，及有损伤结构在 动力荷载作用下的裂纹扩展是造成大型结构体系发生灾难性事故的 主要原因之一。初步研究表明损伤检测、预警及适时维修制度的 建立有助于在一定程度上消除隐患及避免灾难性事故的发生。建 立损伤检测、预警及适时维修制度涉及到两个核心技术问题：是 如何对结构体系中已经出现的损伤进行有效的识别、定位和量化标 定，二是如何对有损伤结构系统的剩余寿命进行定量评估。作为第 一个方向螅研究，自7 0 年代“。1 就已经开始进行研究的结构损伤识别 技术近年里得到了巨大发展。 结构损伤检测可分为两类：一类是利用染色渗透、x 射线、y 射线、 光干涉、超声波和电磁学监测等技术对结构的某些局部进行定期检 查。染色渗透技术是对结构的局部表面进行涂层，涂料则渗透到裂缝 里，观察表面就可发现表面的裂纹；x 射线、7 射线探损技术是利用构 件的x 射线、y 射线照片进行损伤识别；超声波技术是向构件发射高 频声波并测量折射情况，从而识别出复合材料的分层或其他损伤；电 磁学检测技术是利用涡流和磁场的原理进行结构损伤的识别“1 。这类 技术在建筑、航天和船舶等领域有着广泛的应用“1 。但这类技术在 应用上有很多缺点：一是对一些不可见、不开敝的部件难以监测；二 是对于一些大型结构特别是比较复杂的大型结构检测其损伤是不可 能的；三是这类技术要求监测人员必须到现场才能检测。可见这类 技术仅适用于小型结构的检测。对于一些不可见、不开敞的部位，该 类技术不仅无法实施，甚至要求结构的一些功能停止使用或停止工 作。例如，要求以民用客机停飞、发电机组停止运转的情况下进行测 试，这是非常不经济的，甚至无法做到。于是出现另外一种局部损伤 识别技术，这种技术把传感器( 如光纤传感器或压电传感器) 固定在 一些重要部件中，从而对这些部件进行远距离检测。这类技术因为 具有可以直接确定构件的裂纹及其位置的优点，所以在公路、桥梁 和建筑等许多方面得到应用。在这一类技术中基于动力特性的损伤 识别技术受到了广泛的重视。 结构动力学损伤振动检测的基本问题是如何用给定的结构动力 特性确定损伤的出现、位置和程度。通常，结构损伤位置的确定等 价于在结构中用一个可测量的量来确定结构的刚度和承载能力有所 下降的区域。在线损伤探测法可以在结构服役期间通过周期性的参 数识别来探测结构损伤的位置和程度。 1 2 动力学结构损伤检测的研究现状 从7 0 年代初起，在基于动力学特征的结构早期损伤识别问题上， 已经完成了众多的研究工作。这些研究工作可以分为两个步骤：研 究可以用作损伤识别指标的结构动力学特征变化和研究依据相应的 动力学特征变化得到损伤信息的方法。在很多工作中，这两个步骤 往往是同时或相继完成的，在理论分析证明某动力学特征的变化可 以用作损伤指标的过程中，实际上给出来了依据相应变化得到损伤 信息的方法。而近期的一些研究工作则集中在其中的某方面进行。 由于结构动力学损伤识别技术是一门新兴的多门学科前沿知识 交叉的学科，而且目前正处于蓬勃发展之中，所以结构动力学损伤识 别技术有很多种分类方法，主要有按基于损伤检测所用的信息进行 分类、按用不用有限元( f e m ) 建模进行分类、按损伤检测的功能进行 分类、按用不用先验知识进行分类和按损伤识别技术的使用范围进 行分类等等。现依研究所用信息对动力学损伤识别进行分类，并分别 叙述这些方法的特点。 1 基于固有频率变化的损伤识别技术 固有频率是模态参数中最容易获得的种参数。基于固有频率 的变化迸行损伤识别有很多方法 。8 1 。文献”3 对这种技术给出了很好 的总结。这类方法的共同特点是：认为结构发生损伤时，仅结构的刚 度降低，而忽略结构质量的变化，并在结构早期损伤之前建立一个修 正的理论模型。其损伤识别的步骤为：第一，根据理论模型，先假设结 构可能有一组损伤位置的方案，并计算每个理论方案所对应的任意 二个频率改变量之比i 第二，计算实验测量的任意二个频率改变量之 比：第三，将上述理论比值与实验比值进行比较，找出与实测最为接 近的理论值，则该值对应的损伤方案即为实际结构的损伤状态。利用 固有频率的变化进行损伤识别的优点是：测试时，固有频率容易获得 且测试的精度比较高。但是，很多实践表明该类技术在应用上有一些 不足：( 1 ) 固有频率对结构早期损伤有时并不十分敏感，往往只能发 现破损，而无法确定破损的位置。这是因为不同位置的损伤可能引起 相同量的频率变化”：( 2 ) 虽然当损伤的位置在结构的高应力区域时， 利用固有频率的变化进行损伤识别比较可靠，但是当损伤位置在结 构的低应力区域时，利用固有频率的变化将无法进行损伤识别”3 ：( 3 ) 结构早期损伤往往是局部的和微小的，受影响的一般是高阶固有频 率，雨高阶固有频率的变化其计算或测量精度一般较差。所以，利用 固有频率的变化无法识别结构的小损伤。 2 基于振型变化的损伤识别技术 虽然振型的测试精度低于固有频率，但振型包含更多的损伤信 息。利用振型变化识别结构早期损伤的方法很多。“。常用的方法 有： ( 1 ) 模态置信度判据法模态置信度判据法是利用模态置信判 据进行损伤识别( 如m a c 、c o m a c ) 。其原理是：当损伤未发生时，模态 置信度判据为一。可一旦破坏发生，由于振型的变化，模态置信度判 据不等于l 。 ( 2 ) 模态正交法顾名思义，模态正交法是利用模态的正交条 件进行损伤识别。当结构无损伤时，模态满足正交条件。当结构发生 损伤时，则模态不满足正交条件。当然，该方法要用到模型矩阵( 如刚 度矩阵、质量矩阵) ，这就涉及到测量模态的插值扩阶或模型减缩问 题。 ( 3 ) 振型曲率法如果结构出现破损，则破损处的刚度会降低， 而曲率便会增大。振型曲率的变化随着曲率的增大而增大。因此，可 以根据振型曲率的变化确定损伤发生的位置。这种方法以振型曲率 作为定位参数。该方法的不足之处是需要非常邻近的测点，以便利用 中心差分法求取曲率模态。这样就要求足够密的测点，或者要求精度 非常好的插值扩阶模态，否则将增大曲率模态振型的误差“。 ( 4 ) 振型变化图形法该方法是以振型相对变化量作为定位参 数，即损伤前后振型的差值与损伤前振型的比值。当发生破损时，受 到影响的自由度上的振型相对变化量在损伤区域内就会出现比较大 的值。所以，利用振型相对变化图可以识别损伤的位置。 实际上，上述基于振型变化的损伤识别技术在应用中面临着测 量振型不完整( 不仅测量的振型个数少于分析模型的个数，而且测量 的自由度个数也少于分析模型的自由度数) 和噪声的影响问题，当缺 少破损影响较大的测量模态时，该类技术将不能识别结构的损伤。 3 基于柔度变化的损伤识别技术 许多研究者在利用柔度变化进行损伤识别方面做了大量有益的 研究“。主要原理是，在模态满足归一化的条件下，柔度矩阵是频率 的倒数和振型的函数。随着频率的增大，柔度矩阵中高频率的倒数影 响可以忽略不计。这样只要测量前几个低阶模态参数和频率就可获 得精度较好的矩阵。根据获得损伤前后的二个柔度矩阵的差值矩阵， 求出差值矩阵中各列中的最大元素，通过检查每列中的最大元素就 可找出损伤的位置。 4 基于刚度变化的损伤识别技术 当一个结构发生损伤时，剐度矩阵一般提供的信息比质量矩阵 多。利用刚度矩阵的变化进行损伤识别有很多人在研究，因为结构发 生较大的损伤时，其刚度将发生显著的变化。但是，结构发生微小的 损伤( 比如小于百分之一) 时，这类方法将无法进行损伤识别。 5 基于能量变化的损伤识别技术 在利用能量变化识别损伤技术中，由于表达能量所用的参数不 同，这就产生了多种方法。有些方法利用模态参数表达能量，有些方 法不仅用到模态参数同时还引入了有限元模型信息。 文献 1 8 一1 9 介绍了能量传递比( e t r ) 法，把e t r 作为损伤识别的 定位参数。因为有些结构如桥梁，用传统的模态参数作为定位参数不 能有效地进行损伤识别。从文献 1 8 可以得出e t r 法有很多优点：( 1 ) 在损伤区域或靠近损伤区域时，e t r 的值比较大，而远离区域时，e t r 变化不大：( 2 ) e t r 对损伤的敏感程度远远超过固有频率、阻尼比对 损伤的敏感程度，e t r 法不仅能判断损伤的存在，而且能识别损伤的 位置：( 3 ) e t r 法不需要有限元模型，可以用于在线损伤检测。但e t r 法的不足之处有：一是没有考虑噪声的影响，二是假定系统仅有少量 几阶模态，在这个假定的基础上，仅测量前三阶模态，其余模态被截 断。 文献 2 0 描述了另一种能量法，即应变能方法。与文献 1 8 1 9 不同的是该方法同时利用模态参数和有限元信息，其基本原理是利 用结构早期损伤前后应变能发生交化的差值进行损伤识别。应变能 方法成功地应用于桥、板和悬臂结构的损伤识别。文献 2 1 给出结 构单元模态应变能的概念，提出了基于单元模态应变能变化率的结 构早期损伤位置的识别方法。 6 基于传递函数( 频响函数) 变化的损伤识别技术 j i a n ns h i u n l e w ”根据传递函数的变化，提出了一种损伤识别的 方法。因为由于损伤引起的传递函数的变化唯一地由损伤的类型和 位置确定。虽然传递函数或频响函数的信息量大，但是损伤识别仅利 用频响函数的一列数据。d a y i dc 等学者“”利用频响函数数据和有 限元模型，结合最小秩摄动理论，并假设所有产生损伤的方案情况， 成功地对一桁架结构进行了损伤识别。 m a i anmm 等人”4 1 也提出了一种传递函数识别损伤的方法，即频 响函数曲率法。其原理类似于振型曲率法，但不需测试振型，比振型 曲率法识别效果好。该文献中的实例表明，频响函数曲率法可以很好 地识别梁的损伤至少能识别杨氏模量降低2 5 这样的损伤量，并考虑 了5 噪声的影响。不足之处是识别的位置还不够精确。 1 9 9 7 年，m a r kjs 等人”钉提出了另一种传递函数识别损伤的方 法，即传动函数法( t r a n s i n i t t a n c ef u n c t i o n ) 。t f 的定义是结构上任 意二点加速度的互谱与二点中任意一点的自谱的比值。对于相同均 方根幅值的随机激励，t f 是频响函数矩阵列的函数，与输入( 力、位 移速度、加速度无关) 。那么，复传动函数就描述了结构性质的变化。 这样，结构的复传动函数的最大变化反映了结构的损伤情况。因此， 根据t f 的变化就可识别结构早期损伤。t f 法的优点为：一是不需要 数学模型( 如有限元模型等) ：二是可以不需要先验的知识( 如以往的 数据等) ：三是可用于在线损伤监测。t f 法的不足之处是测量点的数 量和位置影响损伤识别的精度。 7 基于统计信息的结构早期损伤识别技术 在系统分析中，损伤可以看作是系统的一种附加激励，它引起系 统输出信号的改变。我们测量的是结构的输出信号，根据输出信号求 附加激励，必须考虑输出信号中的噪声问题，特别是结构的小损伤引 起输出信号的变化可能被噪声信号所淹没。于是，出现了基于统计信 息方法的损伤识别技术托“。 8 基于响应变化的结构损伤识别技术 损伤引起的响应变化一般很小，但是由于测量结构具体位置上 的响应相对容易，因此在这方面有很多研究成果。如在分析响应的 微小变化方面采用小波变换法”7 。“，在识别响应变化方法提出的使 用神经网络法”“1 和遗传算法n ”1 等方法。 除了上述损伤识别技术之外，还有灵敏度法”7 。等方法。 1 3 动力学结构损伤检测技术的发展趋势 动力学结构损伤检测技术总的发展趋势围绕下面四个问题：判 断损伤是否存在、确定损伤的位置、计算损伤量的大小和解决损伤 对结构剩余寿命的影响问题。具体地讲，损伤识别技术的发展趋势有 以下几个方面： 1 结构早期损伤机理的研究。可以讲，到目前为止，判断损伤是 否存在、如何确定损伤的位置和计算损伤量的大小这三个问题仍然 是主要问题，而且目前还无法解决损伤对结构寿命的影响问题： 2 识别微小损伤的研究。对于大型结构，大损伤对系统参数有 较大的影响，而小损伤对系统参数影响较小，再加上噪声的影响，上 述许许多多的方法就无法识别小损伤。然而，如本文上面所述有些大 型结构如航天器、大型客机大型发电机组等，小损伤能在极短的时间 内发展为大损伤，在人们还未察觉的情况下，整个结构已经毁坏。虽 然小损伤对系统参数影响较小，但对局部参数影响较大。所以，在实 用和理论上，这是一个很值得迸步研究的课题； 3 传感器的布置及其数量优化的研究。传感器的数量和布置与 损伤识别技术的成本的精度有关，即使损伤方法的精度再高也无法 改变响应数据匮乏所带来结果的粗糙性。无论在经济上还是在理论 上，这是一个很有意义的课题； 4 研究利用响应数据直接进行损伤识别。利用响应数据直接进 行损伤识别的研究具有经济和实用意义。比如，如何利用地震响应数 据识别地震后的大型桥梁、高耸建筑、核电站和大型水利工程等结 构的损伤问题。 5 在线损伤识别方法的研究。在线损伤识别方法具有实时性、 连续性和预报性，有巨大的实用价值和理论价值且应用前景广阔。目 前，对该类方法的研究是一个发展趋势，但需要解决的理论问题很多 比如，如何处理工作环境激励中存在非自噪声的信号、如何处理环境 激励以面力或体力为主的情况等问题。 6 基于现代信号处理技术的损伤识别方法的研究。小波变换法 等现代信号处理技术在损伤识别方法中的作用有巨大潜力，但将先 进的信号处理技术应用到实际工程结构的损伤识别中，还有很多问 题需要解决。 7 对大型复杂结构的损伤识别。现有的损伤识别方法在应用于 大型复杂结构损伤识别时，由于计算量等原因，不能进行有效识别， 使得这方面的研究和简单结构的损伤识别相比，成为一个具有不同 特点的问题。而实际工程中所需进行损伤识别的对象往往是大型复 杂结构，因此很有必要针对大型复杂结构提出新的方法来解决这一 问题。 8 对多损伤情况的研究。现有的一些损伤识别方法往往是在基 于结构中只有个损伤的假设下进行的，而在工程实际中，由于各 种原因，结构中同时存在多个不同类型的损伤是完全可能的。此时， 现有的方法是否能够有效使用和发展新的方法都是需要进行研究的 方向。 9 非线性损伤识别技术的研究。目前，不仅对线性损伤识别技 术研究比较多，而且非线性损伤识别技术的研究也是一个发展趋势。 实际工程结构是非线性结构，只是非线性的强弱不同；而且即使是线 性的结构损伤发生后也可能转为非线性。虽然非线性技术研究中出 现的困难比线性中出现的困难多，但非线性识别技术符合实际，应用 范围广。 10 基于多学科交叉的损伤识别技术的研究。由于各类大型复 杂结构都有自己的特点，包含了机械工程、建筑工程、材料结构和振 动理论等多门学科知识。因此，研究结构的损伤识别，把振动理论与 信号处理、模式识别、人工智能、控制理论和材料结构等多学科技 术结合起来是一个发展趋势。 1 4 本文所做的工作 本文基于遗传算法研究了含裂纹的悬臂复合材料层合壳体的在 线损伤检测理论与方法。 复合材料层台壳体的在线损伤检测可以通过粘贴在层合壳体表 面的压电智能材料和基于振动分析基础上的信号处理来完成。 首先，建立表面粘贴有压电片材料的悬臂完好壳体及损伤壳体 的有限元动力学模型，并求出其动力响应。 其次，使用相关系数来判定各种损伤情况响应之间的相似程度。 并以此作为适应度指标来通过使用遗传算法来搜索各种损伤情况， 找到其响应与待识别结构响应一致的损伤情况；从而实现损伤识别。 最后，在遗传算法基本框架上，针对损伤识别时适应度函数的 分布特点，对遗传算法进行了优化工作。通过对使用遗传算法进行 损伤识别时的所需计算量进行分析，证明了优化工作有效提高了搜 索效率，使得遗传算法能有效进行损伤识别。并通过对大量识别过 程效率数据进行分析发现识别时所需迭代次数等数值的分布服从一 定的规律。通过对这些规律的进一步分析，不但可以更深的理解遗 传算法的搜索过程。而且在对新结构进行损伤识别时，运用这些规 律可以迅速得到预期识别计算量。 第二章复合材料壳体结构的有限元建模 2 1 引言 复合材料层和结构与传统材料相比，由于其优良的比刚度，更 好的抗疲劳和抗损伤的特性，已经在航空航天工业、汽车工业、海 洋工程、甚至超导工程中得到了广泛应用。但复合材料容易老化， 使用中经常出现裂纹、分层损伤，致使复合材料结构性能下降，严 重时将无法继续使用。因此，需要知道复合材料结构在使用过程中 是否已有损伤、损伤程度如何、是否还能继续使用，同时也需要知 道损伤发生的位置、结构剩余寿命。所以损伤检测与在线监测是复 合材料结构使用中的重要工程技术。目前比较成熟的检测技术有x 射线法、c 扫描法、超声波法等。 但是这些损伤检测方法都是非在线的和局部的，多数手段成本 较高，设备复杂，不能实现实时在线大面积监测。在结构面积较大， 很难以目视确定损伤位置的情况下，如果采用传统的无损检测方法 ( 如超声、c 扫描等) 不仅费时费工，而且成本也较高：在某些如飞 行器飞行中飞机结构的损伤状况检测、对地下管道的损伤状况检测 或者在有毒有害的环境下进行的检测等情况下，人们无法接近，普 通的结构材料损伤检测方法都不适用，使得很多损伤不能及时发现、 维修。这些情况大大限制了复合材料的使用范围。 近年来发展起来的基于动力特性的损伤识别技术，能在线对损伤 进行分类、定位、评估大小，并报警、预测寿命、进行自适应等功 能，给复合材料结构损伤的在线检测提供了广阔的发展前景。 为了能够实现基于动力学特征变化的损伤检测，一般需要对结构 进行动力学分析。针对各种不同类型工程结构的动力学特性分析方 法因此得到了广泛的研究。 在实际工程结构中，壳体结构是非常重要的一种结构，也是大量 应用复合材料的一种结构。所以在壳体结构的动力学特性理论分析 方面已经发表了为数非常众多的研究成果。理论分析能够揭示各种 动力学现象的本质，而且可以作为数值分析的校核准则，但是，由 于弹性力学方程组的复杂性，理论分析只能用来求解一些几何外形、 材料性质和支撑条件都相对比较简单的壳体结构。所以为了能够分 析在工程实际中遇到的具有比较复杂外形、材料性质和支撑条件的 大型结构，在壳体结构动力学理论分析发展的基础上，数值分析方 法和半数值分析方法领域也涌现出了非常众多的研究成果。 有限单元法作为目前最重要的一种数值分析方法，在这方面先后 基于经典壳体理论和剪切变形理论提出了多种单元，包括基于薄壳 理论的轴对称壳体单元、位移和转动各自独立插值的轴对称壳体单 元、轴对称超参数壳体单元、平板壳体单元和相对自由度壳体单元 等等。”1 。在本文所研究复合材料层合壳体振动方面，文献 4 0 - 4 2 中 系统的总结了这方面的理论分析研究工作，文献 4 3 4 6 中总结了 这方面的数值分析研究工作。 在本文的工作中，通过有限元方法获得复合材料层合壳体的强迫 振动响应分为三步进行，首先在2 2 节中确定了单元类型，并讨论 了这种单元的构造方法：其次在2 3 节中讨论了如何使用这种单元 对层舍结构进行计算；最后在2 4 节中叙述了通过m a t l a b 得到强迫 响应的计算过程。 2 2 非轴对称超参数壳体单元 由于考虑到圆柱壳体在产生损伤后，将不能再保持轴对称特性 等因素，所以在本文的有限元模型中选用了sa h m a d ，bmi l o l l s 和 0cz i o i q k i e w ic z 提出的非轴对称的超参数壳体单元“。 由于对于厚板和厚壳，中面法线在变形后仍基本保持为直线， 所以在这种单元中假定中面法线在变形后仍保持为直线，并忽略垂 直于中面的正应力所引起的应变能，每一节点有五个自由度。通过 利用形函数作坐标变换，壳体中面可以是任意曲面。在这种单元计 算中还考虑了横向剪应力影响，故比般基于薄壳理论的单元准确。 因此这种单元不但可以分析厚壳，而且可以分析薄壳。大量的研究 和计算实例表明，这种单元具有计算精度高、适用范围大等优点， 是目前为止最好的壳体单元之一。 2 2 1 曲面单元与映射 如图2 1 所示，单元的局部坐标系为售，刁，f ) ，f = 0 为单元的中面， f = 1 为单元的上、下表面，均为曲面。善= l 和叩= l 是由直线产生 的四个截面。单元节点取在中面上，共取八个节点，即四个角点和四 个边中点。 5 3 图2 1 设结点i 的直角坐标为g 。，y ，) ，则中面0 = o ) 上任一点掌，叩的整体 坐标可以表示如下： x = n ，一，j ，= yz = ，乞， ( 2 1 ) 式中n ，g 川) 为单元的形函数，具体表达式为： 对于角点：n 。= 三( 1 + 毒孝x l + 矾吁x 参古+ 仇叩一1 ) 对于边中点：鲁：o ， ，：妻( 1 一毒z x l + r l , r 1 ) 野。= o ，n i = 三( 1 + g ，毒) ( 1 一，7 2 ) ( 22 ) 设f = + l 代表单元的下表面，f = 1 代表单元的上表面，单元内任 一点的整体坐标( z ，y ，z ) 可以用局部坐标皓，r ，f ) 表示如下： 耻盹呷，掣阱北叩，掣烈 ， 令血，缈，a z 。，为结点f 在厚度方向的坐标差值： 斟懈 代入( 2 3 ) ，得到单元内任一点的整体坐标如下 ( 2 4 ) 髀姒唧蚓 s ， 在( 2 3 ) 一( 2 5 ) 式中，结点坐标向量中的下标t 表示下表面， 表明此时的结点坐标是指与该结点相对应单元下表面处的坐标，类似 的，b 表示上表面、m 表示中面。这样定义出近似的垂直于中面的局 2 2 2 位移函数 在每一结点i 有三个线位移和两个角位移，为了定义角位移，如 图2 2 所示，作三个正交向量日，h ：， 首先，作垂直于中面的正交向量日。： h “= 缸。妙，& ，y( 2 7 ) 1 2 h 。的方向余弦为 k = 等，= 等，b = 等， s ， 圈2 2 其中f 。= ( 出? + a y ? + a z ? y “为结点f 的壳体厚度。令在结点i 再作两 个向量甘。和日：。，它们均切于中面，并正交于日 为了使h ，日2 ，h 。 与方向x ，y ，z 大体一致，采用如下做法，令 h 2 f = h 3 ，f ，h i 。= h 2 l h 3 ( 29 ) 即日：，正交马。于和x 轴，。正交于打：和也，根据向量运算法则 不难由上式求出h 。的方向余弦，。m ”聆。：和：。的方向余弦，2 ；，m ：。， ：，。 如图2 2 所示，结点i 在x ，y ，z 方向的线位移分量分别为珥，r 。w ，。 角位移如图2 3 所示，矽，为壳体中面法线向量- 。绕向量h ：，的转角， 为向量也，绕向量q ，的转角。 图2 3 由于转角痧，在q 。方向产生的线位移为塑2 ，它在儿；方向的 位移分量为孚“，竽聊，孚，由于转角，在也，方向产生的线 位移为迎2 ，它在t ，z 方向的位移分量为墅誓，。华，毕。 利用形函数m 船，卵) ，单元内任意一点的位移可以用中丽螽点位移 矽j 7 = “，v 。，w f ，谚，1 表示如下： 弘f i 比较( 2 5 ) 和( 2 1 0 ) 两式可见，定义单元几何尺寸的( 2 5 ) 一一 1 4 比定义单元位移的( 2 1 0 ) 具有较多的自由度，因此这种单元属于超 参数单元，不会自动满足常应变条件。但从后面对应变的定义中可以 看出，关于刚体位移和常应变的条件是可以满足的。 2 23 整体坐标中的应变 整体坐标g ，y ，：) 中的应变取决于下列矩阵 孔加o w o xo x缸 o uo vo w 矽砂砂 o uo vo w 出出瑟 将位移函数( 2 1 0 ) 代入分别对x ，y ，z 求微商，可得 o u 却o w o xo x 苏 o uo vo w 砂砂砂 锄o v 却 瑟庇& 其中j 为雅可比矩阵 ，= = j 一 2 2 4 局部坐标中的应变与应力 o u 加跏 8 8 ta 考 o u 加细 o r lo r lo r l 抛o v 却 a 8 a 考 ( 2 1 1 ) ( 2 1 2 ) ( 213 ) 在厚壳的分析中，由于壳体中面的法线方向与整体坐标系的z 轴 鱼苗鱼却丝暂堂凿砂一却砂暂出一西西一劬锄西 并不一致，而且法线的方向随着点的位置不同而变化，既然假定垂直 于壳体中面的正应力等于零，就必须求出局部坐标系中的应力与应 变。现假定局部坐标系为( x ，y ，z ) ，其中z 轴正交于壳体中面。 前面已经求成了结点i 的局部正交坐标日峨，打 现在建立 单元内任意一点的局部坐标系b 7 y ，z ) 。 如图2 4 所示，在皖，吼，氕) 点作曲面f = 4 0 ，在此曲面上有两条 空间曲线，一条是7 7 = ，另一条是舌= 磊。再作曲面f = 氕的两个切 向量和d 亏，其中d 手切于曲线7 7 = 7 7 。，d 厅切于曲线孝= 彘，由此可知， d 手2 l i 赛+ 7 毒“奏户孝，斯5 ( i 岛+ 歹考+ 嚆p 叩。式中i ，歹，云为 x ，y ，z 方向的单位向量。 图2 4 作z 垂直于西和d 彳，即正交于f = 厶曲面 ll ：= a ( x a 露2 陲 j 苏 i o r 由此可以求出z 的方向余弦厶，m ，鸭。为使局部坐标g ，y ，一) 与 整体坐标b ，y ，= ) 大体一致，以便于结果的整理和边界条件的处理，对 于其他两个坐标x 和y 用如下方法选取。令 由此求出y 的方向余弦，2 ，m 2 ， y = z “y 。最后，作工正交于y 和z ， 方向余弦1 l 。m i ，月l 。 心。如果z 和x 轴平行，可改用 即x = y7 z ，由此可以求出x 的 现在得到局部坐标系( x ，y ，z ) 的方向余弦矩阵如下： ，。 ，： 厶 0 1 = j 聊1 m 2 肌3 l ( 2 1 6 ) ，l “i n 2 悔j 利用矩阵吲可将整体坐标系中的矩阵( 2 1 1 ) 转换到局部坐标系 中去。 3 u 缸 o u 砂 o u 出 却 a o v 7 砂 a v 庇 o w a a w 砂 o w a z = o f 由此得到局部坐标系中的应变如下： ( 217 ) *2 一。 如暂玉劝 ，、砂一磐砂却 )52( 他0 o t，l 1 j i x= = y 口 1，l，l，lj 伽一良挑一砂加一如加一融加一砂却一功抛一缸如一砂抛一出 s j y 二y y 二 y 口 将矩阵陋】分块，得到 o u 缸 加 乱 c 3 u 西 o y 赢 型型 妙7 。 a a “ 苏岔 = 陋怡 = b ) 式中 点广= 函。，v ，谚，y ，) 。陋；】为如下5 5 矩阵： 陋】= 式中 i i 口【 1 2 口2 t i 口2 + ? 2 口l l z a 3 + 1 3 c t 2 1 3 a l + ，l 口3 其中m 。= n ，f 。 m i r t t o ：i届n 崩 2 6 t 2n 2 口2岛y 2 屈五2 坍1 1 2 2 + m 2 口i n l c t 2 + n 2 t t l l 层，2 + 屈，l 届也+ 屈 m 2 0 f 3 + m 3 1 2 2 ”2 a 3 + 押3 口2 2 y 3 + 卢3 ，z 芦：丑+ 芦j 2 n 3 a 1 - + m l 钙1 3 2 , l + l i t z 3 及y i + 屈y 3 热五l + 屈南 铲t 掌4 - - l l s _ a n , 饥_ a n , 蹦钟宓 屈2 一警m 警饥掣o z , 妊 靠却j2 y ，= l , l “+ m 。m i ，+ g l s g l t = ，s ，“+ 肼，m 2 j + 挖j 玎2 。 ( 2 1 8 ) ( 2 19 ) ( 22 0 ) ( 2 。2 1 ) - 1 8 4 ；4 ：磊 r、unl，l 酶日 = - - 皓 在应力一应变关系式中令盯：= 0 ，可得各向同性弹性固体的5 5 弹性矩阵如下： 【。】= 去 1 卢 1 00 00 00 o o l 一“ o 0 0 o o 1 一t 2 七 o o o 0 0 i 一“ 2 k ( 2 2 2 ) 其中占和u 分别是弹性模量和泊松比，最后两个剪应力项中包含 的系数k 可取为1 2 0 。引入系数k 是为了考虑剪应力分布不均匀的影 响。根据前面所述位移函数，剪应力沿厚度方向接近均匀分布，实际 上是抛物线分布，k 的数值就是两种应变能的比值。 由广义虎克定理，局部坐标系中的单元应力为： 扫 = 吒 盯y f r ， f ，一 r z - = p 始 一 ) + 娥 ( 2 2 3 ) 式中弛) 为初应变，双 为初应力，【d 】为弹性矩阵。在刚度系数计算 中，直接利用局部坐标系中的应力p ，在计算输入结果是，除f 盯 外， 还应输出整体坐标系中的应力。 臣 2 2 5 单元刚度矩阵与结点荷载 r ， 盯。， 彳， 单元刚度矩阵直接在局部坐标系中使用下式计算 一1 9 - “2 r 纠 1，j ： ： o o p y r 盯o o r0j-1l p j j 1，j k k t 锄“和 【七。j ：阻jl d l 反j ，l 喇，7 彤 ( 2 2 5 ) 可以通过数值积分来求出刚度矩阵。需要指出的是，正如轴对称 壳体情况一样，超参数单元在引入一定的几何假设以后和位移及转角 各自独立插值的壳单元是等价的。迸一步的研究工作指出“，对于 作为三维蜕化实体元的超参单元，( 2 1 8 ) 式中所表达的应变分量s 。 。和，。，实际上包含着沿厚度( 即局部坐标系。) 均匀分布和线形分 布的两部分。后者是壳体的弯曲应变，它相当于m in d n 板中的曲率 变化茁引起的应变；前者是壳体的薄膜应变，是壳体中面内的变形( 在 板弯曲问题中，它被忽略) 所引起的。通过量纲分析，可以认为，薄 膜应变能项和剪切应变能项相同，在泛函中也具有罚函数的性质。因 此它相关的刚度矩阵也应是奇异的。否则，在壳体越来越薄时，它也 会造成“锁死”现象。因为这种锁死是由过分的虚假薄膜应变能引起 的，所以称之为薄膜“锁死”。 为了保证非奇异性，且避免剪切“锁死”和薄膜“锁死”，从原 则上说，在m i n d l i n 板单元中讨论的各种方法都可以用于此超参数单 元。文献1 4 引中提出了一种加上应变的方案。对占。，s ，和h ，的弯曲应 变部分和薄膜应变部分以及横剪应变，和凡：，采用不同取样点的插 值表示。但是对于一般形状的超参数单元，整个推导过程和表达格式 相当复杂。因此，从实用角度考虑，一般采用的方法是在对( 2 2 5 ) 式进行数值积分时使用缩减积分方案。 数值积分是在积分区间内选取一定的积分点，累加这些积分点上 的函数值与相应系数的乘积得到积分的近似值。根据选择积分点的位 置和系数计算的方法不同，数值积分也分为多种方法，包括牛顿一科 斯特积分和高斯积分等等。高斯积分法能在同样的积分点数目下实现 更高的精度，是一般采用的数值积分方法。在使用数值积分时，首先 要选择积分点的数目，较多的点数可以获得更高的精度。但是分析表 明，在对( 2 2 5 ) 式进行数值积分时，如果毒，7 方向采用3 3 点的积 分，在对厚板和厚壳进行计算时，能得到较好的结果。而用于薄板和 薄壳时，结果反而不好。其原因是在厚度变小时，发生了“锁死”现 象。而改用2 2 点的积分方案，则可以避免“锁死”现象，精度可以 得到明显的改进，且节省了计算时间。所以对本文中的8 结点单元， 在善，珂和f 方向都采用这样的2 点缩减积分方案。 结点力与结点位移之间存在如下关系： 织 = u ， 一 矿 m ! n ： 8 = s = t 硝 硭砖 七等 潍 ( 22 6 ) 设壳体中面承受了分布载荷，在结点i 分布载荷在z ，y ，z 方 向的分量分别为p j ，p j 。中面上任一点的分布荷载可表示为 p ，= 。p j ，p ，= ，p ? ，p ：= 。p ； ( 2 2 7 ) 在曲面f = 0 上取向量孵和d 审所构成的微小面积撒，由向量运算 法则可知：艘= i d 手斫l = m d 髟叩，其中 川= ( 毒一面o y 一嵩一嚣 2 + 毒一高一考一毒 2 + ( 嵩一毒一毒一骞) 2 ， cz z s ， 由虚功原理，得到分布荷载p 产生的结点荷载如下： 2 2 6 单元质量矩阵 ( 2 2 9 ) 在本文中采用集中质量矩阵，整个单元质量为，每个结点集中 单元八分之一的质量，并略去转动项，得到4 0 4 0 阶的矩阵m 。 咖脚彬忡懈朝 力力办矿 哟，所m 舾眦胁q 肛肛肛仁 =吖il m 鄙弘髟 阻】 m0 0 om - o 00 ，” 在m 对角线上的8 个子矩阵m 是5 5 的对角矩阵。 h 】= 警 0 o 2 3 复合材料层合结构的有限元计算 ( 23 0 ) ( 2 3 1 ) 在上述壳体单元的构造过程中，壳体结构的材料是各项同性单层 材料。下面将叙述通过一定的调整，可以使用这种单元计算复合材料 层合壳体结构。 在复合材料层合壳体中，由于复合材料是各项异性材料，所以其 第j 层主轴弹性矩阵q l 不同于2 2 2 式中各项同性材料的弹性矩阵。q 可由主轴的弹性模量巨o j 、：w 、剪切模量g 和材料的泊松比“，叱 a 2 l l 表示为： q 。= 硝r ) f 刀刃 # , a ( 0 e 0 ) 1 一譬雾 o o oo 0oo 0 g o ) 0 o o 0 g ( r ) k o o o o g ( 0 k ( 23 2 ) 同样引入系数意= 1