振动与冲击 第３Ｏ卷第２期 ＪＯＵＲＮＡＬ ＯＦ ＶＩＢＲＡＴＩＯＮ ＡＮＤ ＳＨ０ＣＫ 一种面向平面多刚柔系统的冲击响应建模和计算方法 贺少华，谢最伟，吴新跃 （海军工程大学机械工程系，武汉４３００３３） 摘 要：结合传递矩阵方法建模灵活和计算效率高的优点，提出了一种基于“传递矩阵”概念的多体系统冲击响 应建模和计算方法。以受基础冲击的平面多刚柔系统为研究对象，采用Ｎｅｗｍａｒｋ 法对元件的方程高阶项进行线性化， 用模态方法处理柔体的变形，建立了一般刚体和典型刚体（刚性均质矩形板、带弹性支撑的刚性均质矩形薄板）、＿＿般柔 体和典型柔体（Ｅｕｌｅｒ—Ｂｅｒｎｏｕｎｉ梁）的冲击扩展传递矩阵，冲击激励包含平动和转动两种成分，给出了基于Ｎｅｗｍａｒｋ１Ｂ法的 系统响应数值迭代求解算法程序。用一个工程实例，通过与有限元方法的对比，验证了方法的准确性，得出了转动冲击激 励成分对总体响应的贡献不能忽略的结论。方法的研究对象虽然只是平面多刚柔系统，但很容易推广到三维的情况。 关键词：冲击；多刚柔体系统；传递矩阵方法；Ｎｅｗｍａｒｋ 法；冲击响应分析 中图分类号：ＴＢ１２２ 文献标识码：Ａ Ｍｏｄｅｌｉｎｇ ａｎｄ ｎｕｍｅｒｉｃａｌ ｓｉｍｕｌａｔｉｏｎ ｏｆ ｓｈｏｃｋ ｒｅｓｐｏｎｓｅ ｆｏｒ ａ ｐｌａｎａｒ ｍｕｌｔｉ－－ｒｉｇｉｄ－ｌｅｘｉｆｂｌｅ・－ｂｏｄｙ ｓｙｓｔｅｍ ＨＥ Ｓｈａｏ—ｈｕａ，ＸＩＥ Ｚｕｉ—ｗｅｉ， ＸＳｎ—ｙｕｅ （Ｄｅｐａｒｔｍｅｎｔ ｏｆ Ｍｅｃｈａｎｉｃａｌ Ｅｎｇｉｎｅｅｒｉｎｇ，Ｎａｖａｌ Ｕｎｉｖｅｒｓｉｔｙ ｏｆ Ｅｎｇｉｎｅｅｒｉｎｇ，Ｗｕｈａｎ ４３００３３，Ｃｈｉｎａ） Ａｂｓｔｒａｃｔ： Ａ ｎｅｗ ｍｅｔｈｏｄ ｗａｓ ｐｒｏｐｏｓｅｄ ｆｏｒ ｓｈｏｃｋ ｒｅｓｐｏｎｓｅ ｍｏｄｅｌｉｎｇ ａｎｄ ｃｏｍｐｕｔｉｎｇ ｆｏｒ ａ ｐｌａｎａｒ ｍｕｌｔｉ—ｉｇｉｒｄ—ｌｆｅｘｉｂｌｅ ｂｏｄｙ ｓｙｓｔｅｍ ｂｙ ｉｎｔｅｇｒａｔｉｎｇ ｔｈｅ ｔｒａｎｓｆｅｒ ｍａｔｒｉｘ ｔｅｃｈｎｉｑｕｅ ｗｉｔｈ ｉｔｓ ｆｌｅｘｉｂｉｌｉｔｙ ａｎｄ ｈｉｇｈｅｒ ｃｏｍｐｕｔａｔｉｏｎａｌ ｅｆｉｃｉｆｅｎｃｙ．Ａｄｏｐｔｉｎｇ Ｎｅｗｍａｒｋ—Ｂ ａｌｇｏｒｉｔｈｍ ｔｏ ｌｉｎｅａｒｉｚｅ ｈｉｇｈｅｒ ｏｒｄｅｒ ｔｅｒｍｓ ｉｎ ｄｙｎａｍｉｃ ｅｑｕａｔｉｏｎｓ ａｎｄ ｍｏｄａｌ ｓｈａｐｅｓ ｔｏ ｅｘｐｒｅｓｓ ｆｌｅｘｉｂｌｅ ｂｏｄｙ ｄｅｆｏｒｍａｔｉｏｎ，ｅｘｔｅｎｄｅｄ ｔｒａｎｓｆｅｒ ｍａｔｉｒｘｅｓ ｏｆ ｇｅｎｅｒａｌ ｒｉｇｉｄ－ｂｏｄｉｅｓ ａｎｄ ｓｐｅｃｉａｌ ｒｉｇｉｄ ｂｏｄｉｅｓ（ｒｉｇｉｄ ｈｏｍｏｇｅｎｅｏｕｓ ｒｅｃｔａｎｇｕｌａｒ ｐｌａｔｅ，ｒｅｃｔａｎｇｕｌａｒ ｓｈｅｅｔ ｐｌａｔｅ ｗｉｔｈ ｅｌａｓｔｉｃ ｓｕｐｐｏｒｔｓ），ｇｅｎｅｒａｌ ｆｌｅｘｉｂｌｅ－ｂｏｄｉｅｓ ａｎｄ ｓｐｅｃｉｌａ ｆｌｅｘｉｂｌｅ ｂｏｄｉｅｓ（Ｅｕｌｅｒ—Ｂｅｒｎｏｕｕｉ ｂｅａｍ）ｉｎ ｓｙｓｔｅｍｓ ｕｎｄｅｒ ｂａｓｅ－ｔｒａｎｓｆｅｒｒｅｄ ｓｈｏｃｋ ｗｅｒｅ ｅｓｔａｂｌｉｓｈｅｄ．Ｔｈｅ ｓｈｏｃｋ ｅｘｃｉｔａｔｉｏｎｓ ｉｎｃｌｕｄｅｄ ｂｏｔｈ ｔｒａｎｓｌａｔｉｏｎａｌ ａｎｄ ｒｏｔａｔｉｏｎａｌ ｃｏｍｐｏｎｅｎｔｓ．Ｎｕｍｅｒｉｃａｌ ｉｔｅｒａｔｉｖｅ ａｌｇｏｒｉｔｈｍ ｐｒｏｇｒａｍ ｂａｓｅｄ ｏｎ Ｎｅｗｍａｒｋ—Ｂ ｍｅｔｈｏｄ ｗａｓ ａｌｓｏ ｇｉｖｅｎ．Ｗｉｔｈ ａｎ ｅｎｇｉｎｅｅｒｉｎｇ ｅｘａｍｐｌｅ，ｔｈｅ ａｃｃｕｒａｃｙ ｏｆ ｔｈｉｓ ｐｒｏｐｏｓｅｄ ｍｅｔｈｏｄ ｗａｓ ｖｅｒｉｆｉｅｄ ｂｙ ｆｉｎｉｔｅ ｅｌｅｍｅｎｔ ｍｅｔｈｏｄ．Ｔｈｅ ｒｅｓｕｌｔｓ ｓｈｏｗ ｔｈａｔ ｔｈｅ ｃｏｎｔｉｒｂｕｔｉｏｎ ｏｆ ｔｈｅ ｒｏｔａｔｉｏｎａｌ ｃｏｍｐｏｎｅｎｔｓ ｏｆ ｓｈｏｃｋ ｔｏ ｔｈｅ ｏｖｅｒａｌｌ ｒｅｓｐｏｎｓｅ ｃａｎ ｎｏｔ ｂｅ ｉｇｎｏｒｅｄ．Ａｌｔｈｏｕｇｈ ｔｈｅ ｏｂｊｅｃｔｉｖｅ ｏｆ ｔｈｉｓ ｓｔｕｄｙ ｗａｓ ａ ｐｌａｎａｒ ｍｕｌｔｉ－ｉｇｉｒｄ－ｌｆｅｘｉｂｌｅ ｓｙｓｔｅｍ，ｔｈｅ ｐｒｏｐｏｓｅｄ ｍｅｔｈｏｄ ｃｏｕｌｄ ｅａｓｉｌｙ ｂｅ ｅｘｔｅｎｄｅｄ ｔｏ ｔｈｒｅｅ—ｄｉｍｅｎｓｉｏｎａｌ ｃａｓｅｓ． Ｋｅｙ ｗｏｒｄｓ：ｓｈｏｃｋ；ｍｕｌｔｉ－ｒｉｇｉｄ—ｌｆｅｘｉｂｌｅ—ｂｏｄｙ ｓｙｓｔｅｍ；ｔｒａｎｓｆｅｒ—ｍａｔｉｒｘ ｔｅｃｈｎｉｑｕｅ；Ｎｅｗｍａｒｋ—Ｂ ｉｔｅｒａｔｉｖｅ ａｌｇｏｒｉｔｈｍ； ｓｈｏｃｋ ｒｅｓｐｏｎｓｅ ａｎａｌｙｓｉｓ 对于复杂系统的冲击瞬态响应研究，目前采用的 建模和计算方法一般包括有限元方法（包括有限元子 方程，这往往是一般技术人员感到费神的地方 。Ｊ。 传递矩阵法由于其建模灵活、计算效率高、无需建 立系统的总体动力学方程等优点，已在包括光学、声 学、电子学、表面科学、机器人学、兵器、航空、航天、机 械等诸多现代工程、科学和技术领域得到了广泛应用。 传递矩阵法在结构力学特别是像转子这样的链式系统 动力学中得到了长期普遍应用，原因在于这种方法频 结构法）和多体系统方法，结合各种数值迭代算法进行 求解。有限元法是结构分析中应用最广泛和强有力的 工具，但对复杂结构使用大量结点导致必须处理和管 理非常大的矩阵，计算工作量大，大刚度梯度系统往往 伴随计算“病态”，导致计算困难。复杂多体系统总体 动力学方程所涉及的矩阵阶次高，通常不小于系统的 自由度数，计算量大，而且必须建立系统的总体动力学 率行列式的阶次低，非常便于编程和数值计算。 针对传统有限元方法和多体系统动力学方法的应 用瓶颈，能否建立无需多体系统总体动力学方程、计算规 模小、程式化高的复杂系统动力学新方法呢？受经典传 递矩阵法进行状态矢量“传递”等思想的启发，芮筱亭教 授 结合多体系统动力学方法和传递矩阵法，在国内外 基金项目：国防“十一五”预研项目（４０１０１０６０１０３０３） 收稿日期：２０１０—０３—２２修改稿收到日期：２０１０—０７—０１ 第一作者贺少华男，博士生，１９８１年生 振动与冲击 ２０１１年第３Ｏ卷 首次建立了多体系统动力学全新方法一多体系统传递矩 阵法，并将其应用在发射动力学与飞行动力学上。 结合瞬态运动学和动力学理论，将多体系统传递 矩阵法理论加以演变和创新从而运用在复杂系统冲击 响应建模和计算中是本文研究要开展的工作。 我们知道，对于复杂连续系统的冲击响应特性研 究，为节约计算资源，提高计算速度，一般采用有限元 缩聚建模（有限元子结构法）和多体系统方法对系统进 行离散化和简化 Ｊ，而目前采用的大多是前一种方法， 将后者引入复杂系统冲击响应建模和计算的研究则比 较少见，结合传递矩阵法的多体系统传递矩阵法在这 方面的应用研究更是无资料和文献可查。本文以一基 础受冲击平面刚柔混合系统为对象，将建立基于多体 系统传递矩阵法冲击响应建模和计算的一般理论，在 此基础上进一步建立转子（轴系）一类系统的冲击响应 多体系统传递矩阵算法，最后，通过与经典有限元方法 的对比，验证本文理论和算法的准确性和工程可应用 性，为复杂系统冲击（瞬态）响应计算提供一种普遍意 义上的理论建模和计算算法。 １线性化方法 一般多体系统元件的动力学方程中，未知量有速 度（角速度）、加速度（角加速度）、三角函数、铰接点处 的内力，以及相应的交叉项、高次项等。对线性时不变 多体系统自由振动，状态矢量满足关系式： ＝Ｚｅ ，其中ｚ为物理坐标下的量，ｚ为对应模态 坐标下的量。所以有： ＝ｉｗＺｅ ． ＝一 Ｚｅｉ 对非周期运动多体系统，上式不成立。为了建立元 件的传递矩阵，必须对元件的动力学方程线性化。现代 数值积分方法将速度ｊ（角速度）、加速度 （角加速度） 表示为位置坐标（转角）的线性函数，即： （ｔ ）＝Ａ（ｔＨ） （ｔ ）＋ （ｔ ） ｚ（ｔ ）＝Ｃ（ｔ ） （ｔ ）＋Ｄ （ｔ ） 式中：Ａ（ｔ ）、Ｂ （ｔ ）、Ｃ（ ）、Ｄ （ｔ ）在ｔ 时刻是 时问ｔ 的已知函数。在速度和加速度等变量采用时 问分段线性化的基础上，在步长△ 时间间隔内，可将ｔ 时刻的非线性动力学方程用 时刻的运动参量线性表 示。在几种数值积分方法中，Ｎｅｗｍａｒｋ 在瞬态响应计 算中最为常用，因为它稳定性和精确度均较高。本文 也将采用Ｎｅｗｍａｒｋ 法（其中 和 均取等于０．５，这 样保证了算法的无条件稳定）来对动力学方程线性化 和最终的数值求解。Ｎｅｗｍａｒｋ］Ｂ线性化方程公式为： ．． ２ ２ ２． Ｚ ｎ＋ｌ 川一 ｎ一＿ｚ一 ： 一 ． 一△ｚ △ｚ 。２ ’ ２平面运动刚柔体元件输入输出端之间扩展 传递矩阵的建立 ２．１平面运动刚体 ２．１．１一般情况 如图１所示，平面运动刚体，输入端为，，输出端为 ０，质心为Ｃ。全局绝对坐标系为 ｘｏｙ；基础相对全局绝 对坐标系的平动坐标系为 ｏ，Ｙ。，转动坐标系为 Ｘ２０ Ｙ ；刚体输入端相对基础的平动坐标系为Ｘ３０ Ｙ ，转 动坐标系为Ｘ４０ Ｙ ；在原点固定在输入端，的连体坐标 系Ｘ４０ Ｙ４中，输出端０的坐标为（ｂ ，ｂ ），质心Ｃ的坐标 为（ｃ ｃ以）。设扩展状态传递矢量为［ ，），，０，Ｍ，ｑ ，ｑ ， １ ｒ，它们的含义分别为相对基础的位置坐标、转动位移、 内力矩和内力。 图１ 考虑基础冲击运动的平面运动刚体 Ｆｉｇ．１ Ｍｏｔｉｏｎ ｍｏｄｅｌ ｏｆ ａ ｐｌａｎａｒ ｒｉｇｉｄ ｕｎｄｅｒ ｂａｓｅ—ｔｒａｎｓｆｅｒｒｅｄ ｓｈｏｃｋ ｅｘｃｉｔａｔｉｏｎｓ 对于基础冲击下的系统响应，考虑系统相对基础 的位移、速度和加速度响应才有意义，所以，考虑在 Ｘ２０２Ｙ２坐标系中： ｃ＝ ，＋Ｃｃ２ｓｉｎ０２＋Ｃｃ１ｃｏｓＯ２＝ ，＋ ，．ｃ ｙｃ ｙ，＋Ｃｃ２ＣＯｓＯ２一ｃｄ ｓｉｎＯ２ ，＋ｙ，，。（１） ０＝ ，＋ｂ２ｓｉｎＯ２＋ｂ１ｃｏｓＯ２＝ ，＋戈，．ｏ Ｙ０ Ｙ ＋６２ｃｏｓ０２一ｂ１ ｓｉｎ０２ Ｙ，＋Ｙ，０ 考虑到０ 、０ 一般均较小，故由式（１）得出刚体输出端 的位置坐标 。、Ｙ。与输出端位置坐标 ，、Ｙ，的传递方程 为： ０： ，＋ｂ２０２＋ｂ１， Ｙｏ＝Ｙ，＋ｂ２一ｂｌ０２ （２） 将 ｃ、Ｙｃ、 ０、Ｙｏ投影到平动坐标系 １０１Ｙ１： ｛Ｘｃ ＝ ｔ ｘｃ，ｘｃ｝ ｛Ｘｏ，Ｙｏ｝ ＝ ｛ＸＯ，Ｘ０｝ （３） 结合基础的运动 、 ，于是由牛顿定律得到在全局绝 对坐标系中： Ｊｍ（Ｊｘｃ＋ ） ｑ ，，一ｇ ，。 （４） 【ｍ（ｙｃ＋Ｙ８）＝ ，一 。 第２期 贺少华等：一种面向平面多刚柔系统的冲击响应建模和计算方法 ９５ 其中ｑ．，、ｑ分别为刚体输入输出端的内力，进一步地，ｏ ，２．１．２刚性均质矩形板 上一小节中，我们研究得出了一般平面运动刚体 将式（３）代入式（４）可得： 』 佗（ ｃ＋２ ｌｙｃ＋ ｌ ｃ＋ １ Ｙｃ＋ Ｂ） ｇ ，，一ｇ ，。（５） 【，孔（ＹＧ一２０１ ｃ一 ｌ ｃ一 ｌ ｃ＋ ）＝ｇｙ，一ｇｙ．０ ．的传递矩阵形式，根据这一一般原理和形式，在这一小 节中，我们将对一具体的平面运动刚体一一一刚性均 质圆盘的传递矩阵进行推导，为后面的计算实例做铺 垫。我们知道，刚性均质圆盘单元也是一类常见的结 构模型单元，特别是在转子（轴系）系统冲击建模中应 用普遍。在二维平面，刚性均质圆盘即为矩形平板单 将式（１）代人式（５），结合本文第二小节中描述的 线性化方法，得到ｔ 时刻ｑ 、ｑＹｌ与ｑ 、ｑＹＯ之间的传 ，，递方程为： 』ｇ ，。 Ａ ｓ，ｘｘｌ＋Ａ ５，２Ｙ，＋Ａ ５，３０２＋ｇ ，，＋Ａ ５， （６） ＬｇｙＤ＝Ａ６１ ，＋Ａ６２Ｙ，＋Ａ６．．．，３０２＋ｇｙ，＋Ａ６，， 式中：各状态变量的系数Ａ 一， ，Ａ ． 和扩展系数Ａ ． ， Ａ ． 为ｔ 时刻系统各参数值以及ｔ 时刻已知的冲击激 励值的函数，其中冲击激励包含了平动激励 、 、 、 、 、 和旋转激励 、 、 两种成分，相对于传统 的复杂系统冲击响应只能考虑平动激励情况，本文建 立的冲击模型更加全面。 根据刚体相对于动点，的绝对动量矩定理： ｄＧ， ＋ｍＲ，＇ｃ×口，＝ （７） 式中，Ｇ，：Ｊｌ ：为刚体相对于连体坐标系原点，的绝 对动量矩，Ｉ，，为刚体相对于，的转动惯量，尺 为刚体 质心Ｃ相对点，的绝对矢径，ｎ，为点，的绝对加速度： Ｒ１ｃ＝ＸＩ，，ｃ ＝ ｛ＸＩ，Ｃ，Ｙｌ，ｃ｝ ０ ＝｛Ｘ，， ｝ ＝ ＼ｘＩ＋２０ｌＹｌ＋０１Ｙｌ＋８１ｙＩ＋ＸＢ，Ｙｌ一 ２０ｌ ，一 １ ，一０１ ，＋ｙ ｝ 将各表达式代人式（７），根据Ｎｅｗｍａｒｋ一 线性化 方法，得到刚体输出端内力矩ｍ。与输入端内力矩ｍ， 之间的传递方程： Ｍ０ Ａ４１ｘｌ＋Ａ４２Ｙ（＋Ａ４，，３０２＋Ｍｌ＋ ，Ａ４ｓｑ ．，＋Ａ４６ｑｙ．，＋Ａ４７ ，，，（８） 同样，和式（６）类似，上式中各传递状态变量的系数及 扩展系数为ｔ 时刻系统各参数值以及ｔ 时刻已知的 冲击激励值的函数。 综合式（２）、式（６）、式（８），同时容易知道 ，＝ 。 ＝ ，于是有受基础冲击平面运动刚体输入输出端状 态矢量［ ，Ｙ， ，Ｍ，ｇ ，ｇ ，１］ 之间ｔ 时刻的扩展传递矩 阵为： Ｕ ＝ （９） 元，如图２所示。设圆盘单元厚度为 ，同样的，考虑在 Ｘ２０２Ｙ２坐标系中： ＸＣ：Ｘ１＋＋—导：Ｘ １＋Ｘ＋ ，１，。ｃ，Ｙ， ｃ：）），，一导 ：Ｙ，一— ６Ｉ２ ｌ＋Ｙ＋ ，ｌ，。ｃ， Ｄ＝ ，＋６＝ ，＋ ，０， Ｙｏ＝Ｙ，一６ ２＝Ｙ，＋Ｙ，．．０ 后续的推导原理和步骤同上小节中的一般情况， 由于篇幅问题，这里不再赘叙。 Ｄ１（０２） 图２考虑基础冲击运动的平面运动刚性均质矩形板 Ｆｉｇ．２ Ｍｏｔｉｏｎ ｍｏｄｅｌ ｏｆ ａ ｒｉｇｉｄ ｈｏｍｏｇｅｎｅｏｕｓ ｒｅｃｔａｎｇｕｌａｒ ｐｌａｔｅ ｕｎｄｅｒ ｂａｓｅ—ｔｒａｎｓＪｆｅｒｒｅｄ ｓｈｏｃｋ ｅｘｃｉｔａｔｉｏｎｓ ２．１．３带弹性支撑的刚性均质矩形薄板 对于带弹性支撑的刚性均质薄板的扩展传递矩阵 的推导，与一般刚体不同之处在于式（５）右边还应包含 弹性力，即： ｆｍ（ ｃ＋２０ｌＹｃ＋Ｏ１Ｙｃ＋ １Ｙｃ＋ ）＝ｇ ，，一ｑ ，０一ｋｘｘ， ｔｍ（ ｃ一２０１ｘｃ一 １ ｃ— ｌｘｃ＋ ）＝ｇｙ，一ｇｙ．０一ｋｙｙ， ．这里，ｋ 、ｋ 分别为弹性支撑在 、Ｙ两个正交方向上的 弹性系数， ，、Ｙ，为刚体相对基础的位移量即弹性支撑 的变形量。 另外，由刚性均质矩形薄板的定义，在图２所示相 同坐标系中有６＝０，即： ｃ ｌ， ，Ｃ ，ｌ， ｘｏ ｘｌ， ｏ ，ｌ 其余的推导原理和步骤同上小节中的一般情况， 由于篇幅问题，这里不再赘叙。在这里，引入带弹性支 撑的刚性均质矩形薄板是为了转子（轴系）系统弹性支 撑位置建模的需要，这在本文后面的计算实例中将会 看到。 振动与冲击 ２０１１年第３０卷 ２．２平面运动柔体 ２．２．１一般情况 柔体状态矢量的定义中，除要考虑刚体状态矢量 中各个状态变量外，还必须考虑反映物体变形相关的 物理量。应用模态方法描述物体变形，将使柔体上任 一点变形随时间的变化可用特征矢量对应的广义坐标 描述。因此，只要在刚体状态矢量定义的状态变量中 加入描述柔体变形对应的广义坐标即可作为柔体状态 矢量的定义。定义平面运动柔体元件的状态矢量为： ［ ，Ｙ，０，Ｍ，ｇ ，ｑ ，ｇ１，ｑ２，…，ｇ ，１］ 其中： ，Ｙ，０，Ｍ，ｇ ，ｑ 的含义与刚体状态变量中的相 同，ｑ ，ｑ ，…，ｇ 为用模态方法描述变形的广义坐标。 采用和图１平面运动刚体相同的坐标系，不同的 是，如果将柔体运动分解为第一阶段的类刚体运动和 第二阶段的纯变形运动，则柔体上任一点． 相对基础的 运动即在坐标系Ｘ２０ Ｙ 中的位矢可表示为： ｒ１ ｒＩ＋Ａｒ。 式中：ｒ，为输入端，在坐标系 ：ｏ：Ｙ：中的位矢，ｒ 为在 连体坐标系Ｘ４０ Ｙ４中柔体上 点相对于连体系原点的 位矢，Ａ为Ｘ４ｏ Ｙ 向Ｘ３ｏ Ｙ，的坐标转换矩阵，ｒ。 可分解 为未变形（变形前）位矢２ 与变形矢量 ，的和，即： ０＝１＂ｌ＋Ａｒ。 ＝ｒ，＋Ａ（ｆ。４，＋ ） 写成对应矩阵形式： ＝ ＋Ａ ］＋ Ａ＝ＣＯＳ０ ２…： 因０：为小量，故上式化简为： ［ ］ ＝［ 】，＋［ ：＋ ６ｕ ＋０ ２ （ｂ＋２ ６＋ ｖ＋） ］ｃ 。 式中：｛ｂ ，ｂ ｝ ’和｛“， ｝ 分别为柔体连体坐标系 Ｘ４０ Ｙ 中 点位置未变形时的坐标向量和变形量，两者 相加即为变形后的坐标向量。 在全局绝对坐标系ｘｏｙ中，柔体上任一点 的坐 标为： ＝蛳＋ 『ｃｏｓ０１ ｓｉｎ０１ １ Ｂ＝ｌ【 一ｓｉｎ０１ ｃｏｓ０１ ＪＩ 将式（１０）代入上式得到： ＝ 『 ，＋０１Ｙ，＋Ｘ ＋ｂ１＋ ＋（０１＋０２）（ｂ２＋ ） １ 【一０１ ，＋Ｙ，＋ｙＢ＋ｂ２＋ 一（０１＋０２）（６ｌ＋“）Ｊ 进而得到： ，＋Ｘ８＋ ＋０１ Ｙ，＋２０１Ｙ，＋０１Ｙ，＋ （０１＋０２）（ｂ２＋ ）＋２（０１＋０２） ＋（０１＋０２） Ｙｌ＋ＹＢ十１）一０ Ｌｘｌ一２０ＬｘＩ一０１ｘＩ一 （ ＋ ）（ｂ ＋Ｍ）一２（ ＋ ）厶一（０ ＋０ ） 取连体糸原点Ｏ 为矩心，利用质点糸干目对十沽动 矩心的绝对动量矩定理： ＋ｍＹｏ＋ｍ ４ｃ ０４ ： （１２） 式中：Ｇ仉＝ｆｒｏ￣×ｒ‘ｏ４ｉｐｄ ｄｙ４是平面运动柔体相对于 点０ 的绝对动量矩，ｒ。 ｃ为柔体质心Ｃ相对于０ 的绝 对矢径，ｍ为柔体质量，０。 为点０ 的绝对加速度， 为 输入输出端对矩心的力矩矢量和。在连体系０４Ｘ Ｙ 中， 有： ｄＧｏ４：未 ＋ ） ＋（ ］ ｐ ＋ ｙ ＋ 一（ ＋ｕ）ｖ￣ｐｄ ｄｙ （１３） 整理得： ：２ ｆ Ｅ（ ）厶叱 ＋ Ｊ ｒ『。 ｌ ２ｐｄ ４ｄｙ４＋Ｊ［（Ｙ４＋ ） 一（ ＋ ）ｖ］ｐｄ ｄｙ ２．２．２ Ｅｕｌｅｒ．Ｂｅｒｎｏｕｎｉ梁 考虑基础冲击运动的Ｅｕｌｅｒ—Ｂｅｒｎｏｕｎｉ梁的力学 方程： —警 ：ｄｘ —＝ ｍｌ Ａ．一 一 （．Ｉ ＋十 ：一 ）Ｊ ｌ］ ｄ 一＝ ｌ：ｍ［ Ａ 【 （ ，＋ ．Ｊ十 ．１ ＋ ）＋ ］ （１４）（１斗Ｊ ｄ ＾ ‘ 其中：ｇ表示微元 的 受力，ｍ为梁的线质量 密度，Ｘ 和 的含义同 式（１１）。 用模态方法可以 直接得到构件小变形。 运动时满足边界条件 的模态。用模态方法处 理梁的纵向变形 （横 向变形“等于零），令： （ ４，ｔ）＝ 图３考虑基础冲击运动的 平面运动Ｅｕｌｅｒ—Ｂｅｍｏｕｎｉ梁 ∑ｓｉｎ ） Ｆｉｇ．３ Ａ Ｅｕｌｅｒ—Ｂｅｒｎｏｕｎｉ ｂｅａｍ ｕｎｄｅｒ ＝１ ｂａｓｅ—ｔｒａｎｓｆｅｒｒｅｄ ｓｈｏｃｋ ｅｘｅｉｔａｔｉｏｎｓ （１５） 第２期 贺少华等：一种面向平面多刚柔系统的冲击响应建模和计算方法 如果用有限阶模态来逼近系统的解，例如取／＇ｔ＝ 结合式（１１），得到： ３，则： （ ）＝∑ｓｉｎ ｇ （ ） 所以，由上述变形公式，综合式（１１），沿等截面均 ［ ］＝［Ｘ １，＂一［＂Ｏ ｌ ｙ，Ｉ 一＋２０ １Ｙ，ｌ ＋一Ｏ ｌＹ ｌ，＋Ｘｔ ￣］ 将上式代入式（１９），利用上文中的Ｎｅｗｍａｒｋ１Ｂ线 质梁轴向积分（１４）式得： ｑ４． （１，ｔ）一ｇ４． （０，ｔ）＝ 『 ，＋ ＋０１ ，＋２０ Ｙ‘，＋ Ｙ，一（０ ＋０ ） ，－１ ｍｌ （０ｌ＋０２） ＋２（０１＋０２）１ ，＋ ｌ＋ Ｌ （０１＋０２）ｏｌ ，＋（０ｌ＋０２）（ １＋ ２）ｂ１ Ｊ ｍ（Ｏ ＋ ｚ）（、 １丁 ｊ２ ，＂ＩＴｇ。） ＋ ２ｍ（Ｏ ＋ ｚ）（ ２．，耵ｑ。） （１６） ｇ４， （１，ｔ）一ｑ４， （０，ｔ）＝ 『（０ ＋０２） ，＋（０ ＋０ ） ＋２（０ ＋０２） １夕，＋］ ｌ （０１＋０２）１Ｙ，＋ ，＋ 一０１ｘ，一 ｌ＋ Ｌ ２０１ｘ，一 ｌ ，一（ １＋ ２）ｂ１ ．Ｊ ｍ（Ｏ。＋ ｚ）（ ＋ ：）（＿ｇ２ ＋ ２。）＋ ２ｍ（Ｏ ）（ ＋ ｚ）（ ，２盯ｑ。）＋ ｍ（ ＋ ） （１７） 式中： 、Ｙ 即为 （０，ｔ）、Ｙ（０，ｔ），ｌ为梁的长度，ｍ为梁 的质量。所以根据本文第－－ｄ，节的Ｎｅｗｍａｒｋ 线性化 方法，得到Ｅｕｌｅｒ—Ｂｅｍｏｕｎｉ梁的两端力ｇ（投影到全局绝 对坐标系）ｔ 时刻的传递方程为： ｇ （１，ｔ）＝ ５．１ ，＋ｓｃｓ２Ｙ，＋ ５，，３０２＋ｑ （０，ｔ）＋ ５，７ｑ１＋￣５，９ｑ３＋ ５，１０ ｑ （Ｚ，ｔ）＝ ６，１ ，＋￣６，２Ｙ，＋ ６，３０２＋ｑ （０，ｔ）＋ ６，７ｑ１＋ ６，９ｑ３＋ ６，１０ （１８） 式中：各状态变量的系数 ， 一， ． ， ， 一， ， 和扩 展系数 ． 。、 ． 。为ｔ 时刻系统各参数值以及ｔ 时刻 已知的冲击激励值的函数。同样，相对于传统的平动 冲击激励形式，这里的冲击激励包含了平动激励和旋 转激励两种成分。 由（１２）式得Ｅｕｌｅｒ—Ｂｅｒｎｏｕｎｉ梁的另一动力学方程： ｍ ｚ上２ ｄ ＋ｍ 上（ ＋ ；）ｄ 一 ｍＪｌ０ ＋ｍ 。４ｃ—ｍｘ。 ）， （１９） 式中：｛Ｘ。 ，Ｙｏ４｝ 点０ 的绝对加速度，由： ＝ ＋ ］ 性化方法和变形的模态处理方法，通过积分可以得到 Ｅｕｌｅｒ．Ｂｅｒｎｏｕｎｉ梁ｔ 时刻两端力矩 的传递方程为： （Ｚ，ｔ）＝ ４．１ ，＋ ４．２），，＋ ４．３０２＋Ｍ（０，ｔ）＋ ，５ｇ４， （０，ｔ）＋ ．６ｑ４， （０，ｔ）＋ ．７ｇ１＋ ．８ｑ２＋￣４，９ｑ３＋ ．１ｏ （２０） 另外，容易知道， （ｚ，ｔ）、Ｙ（ｚ，ｔ）与 （０， ）、ｙ（Ｏ，ｔ） 的传递关系为： （Ｚ，ｔ）＝ （０，ｔ）＋／ｃｏｓ０２， Ｙ（Ｚ，ｔ）＝Ｙ（０，ｔ）一／ｓｉｎ０２ （２１） 又：０（Ｚ，ｔ）＝ｏ（ｏ，ｔ）＝０２，ｇ１（Ｚ，ｔ）＝ｇ１（０，ｔ） ９２（ｆ，ｔ）＝ｑ２（０，ｔ），ｇ３（ｆ， ）＝ｇ３（０，ｔ）（２２） 综合式（１８） 式（２０）、●式（２１）、０ ０ 式（２２），得到考虑基础 鼠０ ０ ０ ０ 冲击运动的平面运动Ｅｕｌｅｒ—Ｂｅｒｎｏｕｎｉ梁输入输出端状 态矢量［ ，Ｙ，０，Ｍ， ，０ ｇ ，ｇ。，ｇ：，…，ｑ ，１］０ 缸 之间ｔ０ ０ ０ ０ 时刻 的扩展传递矩阵为： ０一 鼠０ ０ ０ ０ ０ ０ Ｏ ｌ ０ ０ ０ ０ ０ ０ Ｏ ０ ０ ｌ ０ ０ ０ ０ ０ ０ ０ ０ 亭 ０● ０ ０ ０ ０ ｒ＝ ０ ０ ０鼠 缸 ０ Ｏ ０ ０ ０ ０鼠０ ０ ０ ｌ ０ ０ ０ ０ ０ ９ ９ 鼠９ ０ Ｏ●０ ０ ０ 缸０ ０ ０● ３基础冲击平面刚柔体系统冲击响应算法的 一般步骤 采用本文方法在计算机上求解基础冲击平面刚柔 体系统冲击响应算法的一般步骤如下： （１）建立如图所示一类坐标系，包括全局绝对坐 标系、随动基础坐标系（平动和转动坐标系）、系统相对 基础连体坐标系（平动和转动坐标系）； （２）确定描述柔体元件变形的广义坐标个数，确 定系统的初始条件和边界条件，令ｉ＝１； （３）根据初始条件即ｔ 时刻系统各参数值和各 个时刻均已知的冲击激励形式，求得各元件ｔ 时刻的 扩展传递矩阵和系统总传递矩阵； （４）由边界条件和总传递方程，求得 时刻边界 状态矢量中的未知量； （５）由元件的传递方程一次计算得到ｔ 时刻各个 元件的状态矢量； ９８ 振动与冲击 ２０１１年第３０卷 （６）利用ｔ 时刻的状态矢量计算位置坐标、转角 以及描述柔体元件变形的广义坐标在ｔ 时刻对时间ｔ 激励对响应的贡献大小将在下一小节中单独讨论。在 基础平动加速度冲击激励下，采用有限元方法和本文 的一阶导数和２阶导数； （７）用模态方法求解柔体元件所考查位置的 变形； 方法计算得到的各支撑位置处的最大位移如表１所 示，两种计算方法的最大相对误差不超过１０％，误差的 原因主要来自两个方面：一是实际空间轴系模型简化 成平面“刚体一梁”系统带来的误差；二是用本文方法 （８）令ｉ＝ｉ＋１，用上一时刻的计算结果作为初始 条件，从步骤（３）开始，重复上述步骤，直到所要求的时 计算时，用模态方法处理梁的变形，模态截取是有限 间为止。 ４工程应用实例 ４．１ 系统模型和基础冲击激励形式 某型推进轴系系统模型如图４所示，具体的参数 在文献［５］中所列。根据本文研究理论，由于系统严格 轴对称，所以在弘平面内将其简化为图４所示的平面 多刚柔体模型，其中的“刚体”为螺旋桨（位置１），等效 成上文２．１．２节中的刚性均质矩形板，其集中质量和 转动惯量根据其实际值很容易求得。“柔体”为将弹性 轴等效成的Ｅｕｌｅｒ．Ｂｅｒｎｏｕｎｉ梁。另外，对于支撑轴承 （位置２—６）的处理，我们以等效弹簧代替。在本文研 究的平面多刚柔体系统中，将轴承及其支撑处等效成 带弹性支撑的刚性均质矩形薄板（上文２．１．３节所 述）。因为在以往的机械系统和结构冲击响应计算中， 冲击激励的形式多采用加速度的形式，所以本文也采 用加速度形式的冲击激励，如图５所示，其中图５左图 为基础平动加速度时间历程激励，右图为基础转动角 加速度时问历程激励。采用单一的加速度冲击激励形 式，即上文各表达式中的位移激励和速度激励项均为 零，这使得计算变得相对简单。 位置Ｈ立置２位置３位置４ 位置５ 位置６ 位置７ 图４某推进轴系系统及其平面刚柔体模型 Ｆｉｇ．４ Ａ ｐｒｏｐｕｌｓｉｏｎ ｓｈａｆｔ ａｎｄ ｉｔｓ ｅｑｕｉｖａｌｅｎｔ ｐｌａｎａｒ ｍｕｌｔｉ－－ｉｒｇｉｄ－－ｌｆｅｘｉｂｌｅ ｂｏｄｙ ｓｙｓｔｅｍ ｍｏｄｅｌ 一１４０ 图５基础平动冲击激励及旋转冲击激励函数 Ｆｉｇ．５ Ｂａｓｅ—ｔｒａｎｓｆｅｒｒｅｄ ｓｈｏｃｋ ｉｎｃｅｎｔｉｖｅｓ： ｔｒａｎｓｌａｔｉｏｎａｌ ａｎｄ ｒｏｔａｔｉｏｎａｌ ｆｕｎｃｔｉｏｎｓ ４．２本文方法与有限元方法的对比 为计算程序和计算简单性起见，在验证本文计算 理论准确性时，我们先不考虑旋转冲击激励，对于旋转 的，可能漏掉了后面某一阶或几阶对影响贡献大的模 态，使得与有限元时域瞬态响应计算的结果相差偏大。 表１ 本文方法与有限元方法计算结果对比 Ｔａｂ．１ Ａ ｃｏｍｐａｒｉｓｏｎ ｏｆ ｒｅｓｕｌｔｓ ｂｙ ＦＥＭ ａｎｄ ｔｈｅ ｍｅｔｈｏｄ ｉｎ ｔｈｉｓ ｐａｐｅｒ ４．３旋转冲击激励对响应的影响分析 目前，对于机械结构设备的冲击响应计算特别是 舰船设备的非接触水下爆炸冲击响应计算均没有考虑 旋转激励的成分，这是为了计算的简便，并不是旋转激 励成分任何情况下都能忽略，事实上，美国海军舰载设 备动力学分析方法（ＤＤＡＭ）的最新技术 Ｊ一一一旋转 自度谱响应计算就是从考虑冲击旋转激励成分出发 的，另外，文献［７］用分析力学方法证明了转子系统基 础冲击旋转成分对响应的影响是明显的。在这里，运 用本文提出的建模和计算方法对旋转激励对响应的贡 献进行考证，得到如表２所示的结果，结果表明，即使 表２ 考虑旋转激励对响应的影响分析 Ｔａｂ．２ ｒｏｔａｔｉｏｎａｌ ｃｏｍｐｏｎｅｎｔ ｏｆ ｓｈｏｃｋ ｉｎｃｅｎｔｉｖｅｓ’ ｃｏｎｔｒｉｂｕｔｉｏｎ ｔｏ ｔｈｅ ｏｖｅｒａｌｌ ｒｅｓｐｏｎｓｅ 最大响应（位移，单位：ｍ） 位置 不考虑旋转激励 考虑旋转激励 【下转第１０９页） 第２期 周靖等：脉冲地震动缩放水平对结构非线性位移反应影响的分析 １０９ （３）选择合理的ＩＭ参数、选择合适目标的地震记 录、减小缩放地震动记录的水平，能有效减小估算脉冲 地震作用下结构非线性位移反应的偏差。 参考文献 『Ａ］．Ｔｈｅ ８ Ｕ．Ｓ．Ｎａｔｉｏｎａｌ Ｃｏｎｆ．ｏｎ Ｅａｒｔｈｑｕａｋｅ Ｅｎｇｉｎｅｅｒｉｎｇ［ｃ］．Ｓａｎ Ｆｒａｎｃｉｓｃｏ，２００６． ｉｎｇ ｂｉａｓ ｉｎ ｓｔｕｃｔｕｒａｌｒ ｒｅｓｐｏｎｓｅ ｃａｕｓｅｄ ｂｙ ［７］ Ｂａｋｅｒ Ｊ Ｗ．Ｍｅａｓｕｒｇｒｏｕｎｄ ｍｏｔｉｏｎ［Ａ］．Ｔｈｅ ８ Ｐａｃｉｉｆｃ Ｃｏｎｆ．ｏｎ Ｅａｒｔｈｑｕａｋｅ Ｅｎｇｉｎｅｅｉｒｎｇ［Ｃ］．Ｖａｎｃｏｕｖｅｒ，Ｐａｐｅｒ，２００７，０５６． ｏｒｄ Ｐ Ｊ．Ｎｕｍｂｅｒｓ ｏｆ ｓｃａｌｅｄ ａｎｄ ［８］ Ｈａｎｃｏｃｋ Ｊ．Ｂｏｍｍｅｒ Ｊ Ｊ，Ｓｔａｆｆｍａｔｃｈｅｄ ａｃｃｅｌｅｒｏｇｒａｍｓ ｒｅｑｕｉｒｅｄ ｆｏｒ ｉｎｅｌａｓｔｉｃ ｄｙｎａｍｉｃ ａｎａｌｙｓｅｓ ［１］Ｓｈｏｍｅ Ｎ，Ｃｏｒｎｅｌｌ Ｃ Ａ，Ｂａｚｚｕｒｒｏ Ｐ，ｅｔ ａ１．Ｅａｒｔｈｑｕａｋｅｓ， ｒｅｃｏｒｄｓ，ａｎｄ ｎｏｎｌｉｎｅａｒ ｒｅｓｐｏｎｓｅｓ［Ｊ］．Ｅａｒｔｈｑｕａｋｅ Ｓｐｅｃｔｒａ， １９９８，１４（３）：４６９—５００． ［Ｊ］．Ｅａｔｒｈｑｕａｋｅ Ｅｎｇｉｎｅｅｒｉｎｇ ａｎｄ Ｓｔｒｕｃｔｕｒａｌ Ｄｙｎａｍｉｃｓ， ２００８，３７（１４）：１５８５—１６０７． Ｂ．Ｓｅｉｓｍｉｃ ｄｒｉｔｆ ａｎｄ ｄｕｃｔｉｌｉｔｙ ［９］ Ｋｒａｗｉｎｋｌｅｒ Ｈ，Ｍｅｄｉｎａ Ｒ，Ａｌａｖｉ［２］Ｉｅｒｖｏｌｉｎｏ Ｉ，Ｃｏｒｎｅｌｌ Ｃ Ａ．Ｒｅｃｏｒｄ ｓｅｌｅｃｔｉｏｎ ｆｏｒ ｎｏｎｌｉｎｅａｒ ｓｅｉｓｍｉｃ ａｎａｌｙｓｉｓ ｏｆ ｓｔｒｕｃｔｕｒｅｓ［Ｊ］．Ｅａｒｔｈｑｕａｋｅ Ｓｐｅｃｔｒａ， ２００５，２１（３）：６８５—７１３． ｄｅｍａｎｄｓ ａｎｄ ｔｈｅｉｒ ｄｅｐｅｎｄｅｎｃｅ ｏｎ ｇｒｏｕｎｄ ｍｏｔｉｏｎｓ［Ｊ］． ［３］Ｂａｋｅｒ Ｊ Ｗ，Ｃｏｒｕｅｌｌ Ｃ Ａ．Ｖｅｃｔｏｒ—ｖａｌｕｅｄ ｇｒｏｕｎｄ ｍｏｔｉｏｎ Ｅｎｇｉｎｅｅｒｉｎｇ Ｓｔｕｒｃｔｕｒｅｓ，２００３，２５（５）：６３７—６５３． ｉｎｔｅｎｓｉｔｙ ｍｅａｓｕｒｅｓ ｆｏｒ ｐｒｏｂａｂｉｌｉｓｔｉｃ ｓｅｉｓｍｉｃ ｄｅｍａｎｄ ａｎａｌｙｓｉｓ ［１０］ Ａｋｋａｒ Ｓ，Ｙａｚｇａｎ Ｕ，Ｇｕｌｋａｎ Ｐ．Ｄｒｉｔｆ ｅｓｔｉｍａｔｅｓ ｉｎ ｆｒａｍｅ ［Ｒ］．Ｊｏｈｎ Ａ．Ｂｌｕｍｅ Ｅａｒｔｈｑｕａｋｅ Ｅｎｇｉｎｅｅｒｉｎｇ Ｃｅｎｔｅｒ， ｂｕｉｌｄｉｎｇ ｓｕｂｊｅｃｔｅｄ ｔｏ ｎｅａｒ—ｆａｕｌｔ ｇｒｏｕｎｄ ｍｏｔｉｏｎｓ［Ｊ］．Ｊｏｕｒｎａｌ Ｓｔａｎｆｏｒｄ，ＣＡ，Ｒｅｐｏｒｔ，２００５，１５０，３２１． ｏｆ Ｓｔｕｒｃｔｕｒａｌ Ｅｎｇｉｎｅｅｒｉｎｇ，２００５，１２１（７）：１０１４—１０２４． ［４］Ｌｕｃｏ Ｎ，Ｂａｚｚｕｒｒｏ Ｐ．Ｅｆｆｅｃｔｓ ｏｆ ｅａｒｔｈｑｕａｋｅ ｒｅｃｏｒｄ ｓｃａｌｉｎｇ ｏｎ Ｂａｋｅｒ Ｊ Ｗ．Ｑｕａｎｔｉｔａｔｉｖｅ ｃｌａｓｓｉｆｉｃａｔｉｏｎ ｏｆ ｎｅａｒ—ｆａｕｌｔ ｇｒｏｕｎｄ ｎｏｎｌｉｎｅａｒ ｓｔｒｕｃｔｕｒａｌ ｒｅｓｐｏｎｓｅ［Ｒ］．ＰＥＥＲ Ｌｉｆｅｌｉｎｅｓ Ｐｒｏｇｒａｍ， ｍｏｔｉｏｎｓ ｕｓｉｎｇ ｗａｖｅｌｅｔ ａｎａｌｙｓｉｓ Ｊ １． Ｂｕｌｌｅｔｉｎ ｏｆ ｔｈｅ Ｒｅｐｏｒｔ ｏｎ ＰＥＥＲ—ＬＬ Ｐｒｏｇｒａｍ Ｔａｓｋ １ Ｇ００，Ａｄｄｅｎｄｕｍ， Ｓｅｉｓｍｏｌｏｇｉｃａｌ Ｓｏｃｉｅｔｙ ｏｆ Ａｍｅｒｉｃａ， ２００７， ９７（５）： ２００５，６４． １４８６—１５０１． ［５］Ｌｕｃｏ Ｎ，Ｂａｚｚｕｒｒｏ Ｐ．Ｄｏｅｓ ａｍｐｌｉｔｕｄｅ ｓｃａｌｉｎｇ ｏｆ ｇｒｏｕｎｄ ｍｏｔｉｏｎ ［１２］ ｈｔｔｐ：／／ｓｔａｎｆｏｒｄ．ｅｄｕ／一ｂａｋｅｒｊｗ／ｐｕｌｓｅ—ｃｌａｓｓｉｉｆｃａｔｉｏｎ．ｈｔｍｌ， ｒｅｃｏｒｄｓ ｒｅｓｕｌｔ ｉｎ ｂｉａｓｅｄ ｎｏｎｌｉｎｅａｒ ｓｔｕｒｃｔｕｒａｌ ｄｒｉｆｔ ｒｅｓｐｏｎｓｅｓ？ ２００９．０８．２０． ［Ｊ］．Ｅａｒｔｈｑｕａｋｅ Ｅｎｇｉｎｅｅｒｉｎｇ ａｎｄ Ｓｔｕｒｃｔｕｒａｌ Ｄｙｎａｍｉｃｓ， ［１３］ 周靖，陈凯亮，罗高杰．速度脉冲型地震地面运动强度 ２００７，３６（１３）：１８１３一ｌ８３５． 表征参数评估［Ｊ］．振动与冲击，２０１０，２９（７）：１５３—１５８． ［６］Ｂａｚｚｕｒｒｏ Ｐ，Ｌｕｃｏ Ｎ．Ｄｏ ｓｃａｌｅｄ ａｎｄ ｓｐｅｃｔｕｒｍ—ｍａｔｃｈｅｄ ｎｅａｒ— ｓｏｕｒｃｅ ｒｅｃｏｒｄｓ ｐｒｏｄｕｃｅ ｂｉａｓｅｄ ｎｏｎｌｉｎｅａｒ ｓｔｕｒｃｔｕｒａｌ ｒｅｓｐｏｎｓｅｓ？ （上接第９８页） 去除误差因素，旋转激励成分使得总体响应偏大（位置 击环境。经典有限元方法验证了本文理论和算法程序 ６除外），这一结论虽不能断定具有普遍意义，但它至少 的准确性。本文证明了转动冲击激励对响应的贡献不 说明旋转激励成分对响应的影响是不能忽略的，这也 能忽略，这与以往采用不同方法得出的研究结论相一 与文献［７］的结论相一致，虽然两者采用的是两种不同 致，另外，表２的计算显示了本文算法的另一优越性， 的计算方法。 即目前一般的商业有限元软件无法进行旋转冲击激励 ５ 结 论 响应计算，而本文算法却能实现这一点。 参考文献 结合多体系统动力学理论、传递矩阵法理论和数 ［１］Ｒｕｉ Ｘ Ｔ，Ｗａｎｇ Ｇ Ｐ，Ｌｕ Ｙ Ｑ，ｅｔ ａ１．Ｔｒａｎｓｆｅｒ ｍａｔｉｒｘ ｍｅｔｈｏｄ 值算法依据，本文建立了一套完整的复杂（多体）系统 ｏｆｒ ｌｉｎｅａｒ ｍｕｌｔｉｂｏｄｙ ｓｙｓｔｅｍ ｄｙｎａｍｉｃｓ［Ｊ］．Ｍｕｈｉｂｏｄｙ Ｓｙｓｔｅｍ 冲击（瞬态）响应算法理论，该算法主要理论根据是将 Ｄｙｎａｍｉｃｓ，１９（２００８）：１７９—２０７． ［２］Ｌｕ Ｗ Ｊ，Ｒｕｉ Ｘ Ｔ．Ｓｔｕｄｙ ｏｎ Ｄｙｎａｍｉｃ Ｓｉｍｕｌａｔｉｏｎ ｏｆ Ｃｏｎｔａｉｎｅｒ “传递矩阵”概念引入到多体系统动力学理论中，“传递 Ｒｅｐｌｅｎｉｓｈｍｅｎｔ Ｓｙｓｔｅｍｓ．Ｐｒｏｃｅｅｄｉｎｇ ｏｆ ｔｈｅ Ｉｎｔｅｒｎａｔｉｏｎａｌ 矩阵”的引入克服了经典有限元方法、多体系统方法建 Ｃｏｎｆｅｒｅｎｃｅ ｏｎ Ｍｅｃｈａｎｉｃａｌ Ｅｎｇｉｎｅｅｒｉｎｇ ａｎｄ Ｍｅｃｈａｎｉｃｓ［Ａ］． 模求解复杂系统力学问题的“弊端”。我们一旦建立起 Ｎｅｗ Ｙｏｒｋ：Ｓｃｉｅｎｃｅ Ｐｒｅｓｓ ＵＳＡ Ｉｎｃ．，２００５：７４２—７４５． 了多体系统元件响应求解需要的扩展传递矩阵，通过 ［３］Ｒｕｉ Ｘ Ｔ，Ｈｅ Ｂ．Ａｄｖａｎｃｅ ｉｎ Ｄｉｓｃｒｅｔｅ Ｔｉｍｅ Ｔｒａｎｓｆｅｒ Ｍａｔｉｒｘ Ｍｅｔｈｏｄ ｏｆ Ｍｕｈｉｂｏｄｙ Ｓｙｓｔｅｍ Ｄｙｎａｍｉｃｓ［Ｍ］．Ｔｈｅ Ｓｅｃｏｎｄ 组装系统总的传递矩阵，结合边界条件经过数值积分 Ｉｎｔｅｒｎａｔｉｏｎａｌ Ｃｏｎｆｅｒｅｎｃｅ ｏｎ Ｄｙｎａｍｉｃｓ Ｖｉｂｒａｔｉｏｎ ａｎｄ Ｃｏｎｔｒｏｌ， 就能得到响应结果，并且，这些典型元件的冲击扩展传 Ｂｅｉｊｉｎｇ，Ｃｈｉｎａ，２００６． 递矩阵一旦建立好，我们就可以“一劳永逸”的在其他 Ｎｏｉｓｅ ａｎｄ Ｖｉｂｒａｔｉｏｎ Ｃｏｎｔｒｏｌ，４（２００３）：１６—１８． 包含它们的系统的冲击响应计算中直接使用它们，无 ［４］陈志敏．基于刚柔混合建模的柴油机抗冲击特性研究 需再重新计算，所以我们甚至可以建立一个元件冲击 ［Ｄ］．武汉：海军工程大学，２００７． ［５］朱小平，冯奇．船舶主推进系统冲击研究［Ｊ］．船舶力学， 扩展传递矩阵库，直接供计算查找使用，这样大大简便 ２００７，１１（１）：１４３—１５１． 了计算程序，节省了计算资源和时间。此外，本文对于 ［６］汪玉，华宏星．舰船现代冲击理论及应用［Ｍ］．北京：科 冲击激励的考虑，不仅仅局限于传统的平动冲击激励 学出版社，２００５． 形式，还包含了转动冲击激励形式，这使得所提算法理 『７］Ｓｕａｒｅｚ Ｌ Ｅ，Ｓｉｎｇｈ Ｍ Ｐ，Ｒｏｈａｎｉｍａｎｅｓｈ Ｍ Ｓ．Ｓｅｉｓｍｉｃ Ｒｅｓｐｏｎｓｅ Ｏｆ Ｒｏｔａｔｉｎｇ Ｍａｃｈｉｎｅｓ［Ｍ］．Ｅａｒｔｈｑｕａｋｅ Ｅｎｇｉｎｅｅｒｉｎｇ 论更具有广泛的适应性，也可以更真实地仿真一般冲 ａｎｄ Ｓｔｒｕｃｔｕｒａｌ Ｄｙｎａｍｉｃｓ，２１（１９９２）：２１—３６．
