基于弹塑性接触的柔性多体系统碰撞动力学

经 引起工程 界 的重视 ． 例如 ， 受 冲击 载 荷作 用 时 ， 在
的动 力学模 型 ， 假 设模 态 法 进 行 离 散 ， 究 了在 用 研
角速 度 已知 的情 况 下 弹塑 性 应 变 对梁 的振 动 特性
的影 响． 此 基 础 上 ， ｅｓ ｙ 等 进 一 步 研 究作 在 Ｇ ｒｍａｒ ｔ 自由旋 转运 动单摆 的动力学 性态 ， 考虑 弹 塑性 变形
在汽 车碰撞 过程 中 ， 某些 致伤 部位 的损 伤 机理 分析 也是 目前需 要解 决 的问 题． 由此 可见 ， 塑性 多体 弹 系统 动力学 的研究 具有 重要 的工程 应用 前景 ．
性 应 变一 移 关 系 ， 虚 功 原 理建 立 了 单 柔性 梁 的 位 用
动力 学变分 方 程 ， 用假 设 模 态 法 离 散． 此基 础 利 在
性梁系统的刚一 柔耦合动力学方程． 对重 力作用下 的柔性 单摆和 双摆数值仿 真结果表 明， 塑性应 变引起横 向
变形绝对值增大和横 向振动振幅衰减 ， 角加速度 突变 时塑性效 应最 为显 著． 在 关键词 弹塑性梁系统 ， 非线性本构关 系 ， 动力学
引 言 பைடு நூலகம்
目前 ， 机械 系统各 部件 的弹塑性 动力 学 问题 已
１ 梁 变 形 的描 述
平 面柔 性梁 如 图 １所 示 ， 立惯性 基 ｅｅ 建 和柔
性 梁 的浮动基 ｅｅ． ｒ 为 浮动基 原点 关 于惯性 基 ｂｂ设 ０ 。 的位置 向量在 惯性 基 下 的坐标 阵 ， ０为浮 动基 关 于
目前 ， 国内外对 弹塑性 多体 系统 动力 学研 究 才 刚起 步． ｅｓ ａｒ 等建 立 了作 旋转 运动 弹塑 性 梁 Ｇ ｒｍ ｙ ｔ

— 98 —
状态 。 恢复系数方法是一种简化的 处理方法 , 在真 实状况下 影响碰撞行为的参数 很多 , 如相对速度 、物体 的形状 、材料的 特性等 , 恢复系数方法 抛开 这些 因素 , 仅用 接触 刚度和 碰撞 恢复系数两个参 数来 计算碰 撞行 为 。 如 1 图所 示两个 球体 相对运动并发生碰撞 。 假设碰撞过程 在瞬 间完成 , 在 碰撞 瞬间 , 物体 的位 置保 持不变 , 而速 度 发生 跳跃 变 化 。 恢复 系 数的 选择 有 多种 方 法 。 最常用的是 牛顿 恢复 系数 ε N, 假定 碰撞 前后 的接 触点 的法向相对速度可由 恢复系数确定 , 即 :
图 1 两个球体碰撞示意图
() 系数 ; δ 为碰撞前两个物体相对速度 。 由于推导过程中假设 e ≈ 1, 所以得到的表达式只 能适合 ·
如图所示 I 、J 两 个 的 物体 发 生碰 撞 , 定 义 法 向 接触 力 F , 正值表示分离力 ; K 表示刚度系数 ; δ 为碰 撞物体变 形量 , n 或者为侵入深度 ; e 为渗透 深度 的指数 , 简 称刚 性指数 ; D为 阻尼系数 ; δ 为两个物 体的相 对速 度 。 则法 向接 触力的 计算 公式为 :
第 25 卷 第 10期 文章编号 : 1006 9348( 2008) 10 0098 -04
计 算 机 仿 真
2008 年 10 月
多体系统动力学仿真中的接触碰撞模型分析
安雪斌 , 潘尚峰
( 清华大学精密仪器与机械学系 , 北京 100084) 摘要 : 为了明确多体系统动力学仿真时接触碰撞模型的选 择依据 , 对 恢复系数 模型 , I M P A C T 函 数模型 , 迟 滞阻尼模 型三种 模型进行了分析讨论 , 并指出了模型中各参数的物理意义 。 最后以小球与平面 碰撞问题为例 , 通过理论推导与 试验相结合 的方法确定了模型的参数 , 在 A D A M S 软件中进行了仿真计算 , 得到了三种模型的计算结果 。 计算结果表明 , 恢 复系数方法 只能得到碰撞后的速度 ; I M P A C T 函数模型可以得到较好的仿真结果 , 但是参数确定困难 ; 修正 的迟滞阻尼模型 可以解决恢 复系数接近 1的接触问题 。 关键词 : 接触 ; 恢复系数法 ; 迟滞阻尼模型 ; 动力学仿真 中图分类号 :T P 391. 9 文献标识码 : B

1 1 · · K 、 设M 、 Δq 、 Δp 分别是悬臂杆在模态坐标下 —1 —1
4
算例仿真与分析
针对图 1 所示碰撞系统， 设计悬臂杆为铝合金杆 和线 性 强 化 材 料， 长 度 l1 = 0. 3 m， 横 截 面 积 A1 =
引言
杆碰撞现象普遍出现在工程机械和工程结构 中
［1 ］
验研究
［68 ］
。
动态子结构法通过缩聚系统自由度， 大幅度提高 计算效率， 已经被广泛运用于振动分析中。Wu S C 和 Haug E J［9］390-398 采用约束删除—添加技术， 将动态子 结构方法首次应用于弹性单杆—单梁碰撞问题。郭安 萍、 洪嘉振和刘锦阳等 春等
图4 Fig． 4
固定界面子结构
Fixed interface substructure
第 36 卷第 1 期
孔德平等： 双柔性杆弹塑性碰撞动态子结构方法的研究
003
（ s） 力控制方程。子结构交界面上的位移 a J ，实际是子
结构之间保证满足位移协调条件的公共坐标 ， 在拼装 需要消除。 过程中会产生冗余，
＊ ＊ ＊ 孔德平 ＊ 尹晓春 ＊ 骞朋波 沈煜年 ( 南京理工大学 理学院，南京 210094 )
KONG DePing
YIN XiaoChun
QIAN PengBo
SHEN YuNian
( School of Science，Nanjing University of Science and Technology，Nanjing 210094 ，China)
［1315 ］ ［1012 ］ ［2 ］
。例如冲击锤机械系统
中软锤和硬杆之间的
［4 ］
撞击回弹现象、 撞击铆接机

多体系统碰撞动力学中接触力模型的研究进展王庚祥;马道林;刘洋;刘才山【期刊名称】《力学学报》【年(卷),期】2022(54)12【摘要】接触碰撞行为作为大自然与多体系统中的常见现象,其接触力模型对于多体系统的碰撞行为机理研究与性能预测至关重要.静态弹塑性接触模型与考虑能量耗散的连续接触力模型是研究接触碰撞行为的两类不同方法,在多体系统碰撞动力学中存在诸多共性与差异.本文分别从上述两类接触模型的发展历程入手,详细介绍了两类模型的区别与联系.首先,根据阻尼项分母中是否含有初始碰撞速度将连续接触力模型分为黏性接触力模型与迟滞接触力模型,讨论了能量指数与Hertz接触刚度之间的关系,阐述了现有连续接触力模型在计算弹塑性材料接触碰撞行为时存在的问题.其次,着重介绍了分段连续的准静态弹塑性接触力模型(可连续从完全弹性转换到完全塑性接触阶段),分析了利用此类弹塑性接触力模型计算碰撞行为的技术特点.同时,以恢复系数为桥梁和借助线性化的弹塑性接触刚度,避免了Hertz刚度对弹塑性接触刚度的计算误差,根据碰撞前后多体系统的能量与动能守恒推导了弹塑性接触模型等效的迟滞阻尼因子.探索了连续接触力模型与准静态弹塑性接触力模型之间的内在联系,数值计算结果定量说明了人为阻尼项代表的能量耗散与弹塑性接触力模型中加卸载路径代表的能量耗散具有等效性.另外,为了避免阻尼项分母中初始碰撞速度在计算颗粒物质动态性能时导致的数值奇异问题,通过求解等效的线性单自由度欠阻尼非受迫振动方程获得了阻尼项分母中不含初始碰撞速度的连续接触力模型,并以一维球链为例,证明了该模型相比EDEM软件使用的连续接触力模型具有更高的精度.最后,本文分析了当前多体系统碰撞动力学的研究现状,并简要展望了多体系统碰撞动力学中接触力模型的发展趋势与面临的挑战.【总页数】28页(P3239-3266)【作者】王庚祥;马道林;刘洋;刘才山【作者单位】西安建筑科技大学机电工程学院;北京大学工学院;上海交通大学船舶海洋与建筑工程学院;埃克塞特大学工程、数学和物理科学学院 QF【正文语种】中文【中图分类】TH113.1【相关文献】1.汽车-行人碰撞多体动力学仿真与抛距模型研究2.计及侵入体作用的人体系统碰撞动力学模型研究3.可变形多体系统动力学及接触力学在凸轮机构学中的应用述评4.浮筏系统中的碰撞动力学离散模型(2):随机模型5.浮筏系统中的碰撞动力学离散模型(1):确定性模型因版权原因，仅展示原文概要，查看原文内容请购买。

第25卷第6期高压物理学报V o l.25,N o.6 2011年12月C H I N E S EJ O U R N A L O F H I G H P R E S S U R E P H Y S I C S D e c.,2011文章编号:1000-5773(2011)06-0508-06三维多物质弹塑性流体动力学E u l e r方法的并行算法研究及程序测试*马天宝1,费广磊1,张文耀2(1.北京理工大学爆炸科学与技术国家重点实验室,北京100081;2.北京理工大学计算机学院智能信息技术北京市重点实验室,北京100081)摘要:并行计算是解决爆炸与冲击问题大规模数值模拟最有效的手段之一㊂针对E u l e r方法并行程序设计的复杂性,阐述了三维多物质弹塑性流体动力学程序MM I C-3D并行设计的总体策略,基于消息传递接口(M P I)设计出相应的P MM I C-3D并行程序,并提出了一套实用的程序测试方案㊂结合聚能射流形成过程的数值模拟算例,在八节点的集群上测试了加速比㊁并行效率及可扩放性,分析了影响并行性能的因素㊂关键词:爆炸与冲击;E u l e r方法;并行计算;消息传递接口;程序测试中图分类号:O347;O358文献标识码:A1引言近年来,随着大规模科学与工程计算的需求,许多计算问题已经超出单机所能承受的能力范围,并行计算对于大规模科学与工程计算越来越重要[1]㊂目前,并行计算机的基本存储方式主要有共享存储与分布式存储两种㊂M P I(M e s s a g eP a s s i n g I n t e r f a c e) 消息传递接口是消息传递函数库的标准规范,是目前广泛使用的并行编程工具[2],M P I基于分布式存储,但同样适应于共享存储,具有移植性好㊁功能强大㊁效率高等多种优点㊂MM I C-3D(M u l t i-M a t e r i a l i nC e l l f o r3D)[3]是基于E u l e r型有限差分方法的用于爆炸与冲击问题仿真计算的三维多物质弹塑性流体动力学程序,能处理3种及3种以上物质混合格的界面计算问题,解决了三维E u l e r方法中多物质计算的难题,实现了对空中爆炸㊁密实介质中爆炸及聚能射流形成等典型的爆炸与冲击问题的数值模拟㊂在单机32位操作系统下,受限于计算机内存和计算速度,原有MM I C-3D 串行程序最多只能计算200万网格的问题,远远达不到工程和科研的需求,迫切需要将原有的串行程序改造为并行程序,以扩大计算规模,加快求解速度㊂并行程序设计不但包含了串行程序设计,而且还包含了更多富有挑战性的问题[4],由于并行程序需要通信㊁同步等操作,使得并行程序设计远比串行程序设计复杂得多㊂L a g r a n g e方法,无论采用显式或隐式差分格式,最后都会归结为代数或矩阵运算,对于这类运算都有现成的并行算法可用㊂采用E u l e r 方法,由于需要处理物质在网格间的输运问题,隐含在输运算法中的数据相关性及子区域间的关联性不易发现㊂本文针对E u l e r方法并行设计的复杂性,阐述了处理E u l e r输运算法中处理数据相关性的一种方法,以及挖掘子区域间关联性的过程㊂程序测试是程序开发的一个重要环节,考虑到并行程序的测试复杂性,将测试分为两个阶段,提出了一些实用的测试策略,缩短了整个程序的开发周期,并在八节点集群上测试了程序的并行性能㊂*收稿日期:2010-08-20;修回日期:2010-11-26基金项目:国家重点基础研究发展计划(2010C B832706);爆炸科学与技术国家重点实验室自主课题(Z D K T10-03b);国家自然科学基金(10972041)作者简介:马天宝(1981 ),男,博士,副教授,主要从事计算爆炸力学研究.E-m a i l:m a d a b a l@b i t.e d u.c n2 MM I C -3D 数学模型及数值方法介绍MM I C -3D 串行程序采用不考虑外力㊁外源和热传导,非守恒形式的E u l e r 流体弹塑性动力学偏微分方程组[5]㊂数值计算采用算子分裂格式,将上述非守恒形式的方程组分为L a g r a n g e 步和E u l e r 步进行计算,并在计算中按空间x ㊁y ㊁z 3个方向进行分裂㊂MM I C -3D 采用模糊界面方法,该方法是指在一个混合网格中,不区分物质界面;根据模糊方法计算体积比,把体积比作为模糊权重系数;对网格进行分类,不同类网格之间视为物质界面;对介质进行模糊排序,决定输运优先权和模糊输运表;根据模糊权重计算输运量;按模糊输运表进行输运;在建模和计算中应用模糊方法,故称为 模糊界面方法[6]㊂3 P MM I C -3D 并行算法实现的关键性问题P MM I C -3D (P a r a l l e l i z a t i o n f o rM u l t i -M a t e r i a l i nC e l l -3D )采用域分解并行策略,域分解是指把计算域分成若干子区域,一个处理器处理一个或多个子区域㊂采用域分解并行策略,P MM I C -3D 的并行算法设计主要应考虑如下问题:(1)如果原有串行算法有数据相关性,如何解除相关性;(2)采用计算域分区的并行方式,各子区域间关联性如何,即需要多少额外的网格储存临近子区域的信息㊂3.1 数据相关性分析MM I C -3D 程序采用欧拉型有限差分算法,分L a g r a n ge 步和E u l e r 输运步分别进行计算,使得当前计算网格最多受到周围26个网格的影响,这是能采用计算域分区方式对MM I C -3D 并行的有利条件,同时必须分析计算程序所有语句间的依赖关系,称之为相关分析(D e p e n d e n c y A n a l y s i s )㊂对于采用域分解策略的三维多物质流体弹塑性并行程序来说,在1个时间步内,在同1个空间三重循环下若满足:当前网格物理量的更新依赖于周围网格相关变量更新后的值,则该算法会存在数据相关性㊂用表达式表示为a ᶄi ,j ,k =f (a i ,j ,k ,a ᶄi ʃ1,jʃ1,k ʃ1)(1)其中a 的更新直接或间接依赖于周边网格a 更新后的值㊂采用E u l e r 方法的数值模拟,数据相关性一般出现在E u l e r 输运步上㊂模糊界面方法采用方向分裂 输运,3个方向的输运在1个空间三重循环下进行,当前网格的物理量的更新依赖于周边相关物理量更新后的值,更新后网格物理量影响周围网格相关量的更新,因此存在数据依赖性,且这种依赖关系影响到整个计算域㊂如何解除数据相关性,且保持原有的计算精度,要与具体算法结合起来㊂为消除P MM I C -3D 中E u l e r 输运步的数据相关性,当前网格输运采用相邻网格更新前的值,同时为了消除因此而带来过量输运的问题,将原有在1个空间三重循环下完成的3个方向的输运改为在3个空间三重循环下完成,即每个空间三重循环只进行1个方向的输运,3个方向的输运顺序随时间步交替变换㊂为考核因解除数据相关性而带来的精度影响,采用H a r i v e 和F l e t c h e [7-8]设计的测试精度的方法㊂该方法定义了L 1误差,其公式如下E =ði ,j ,k A (i ,j ,k )ΩF (i ,j ,k )t =T-F (i ,j ,k)e (2)式中:A(i ,j ,k )Ω为网格的体积,F (i ,j ,k )t =T为数值计算的网格介质体积分数,F (i ,j ,k )e为准确的网格介质体积分数㊂数值算例计算参数见表1,其结果为计算到400步的结果㊂由表1可以看出,并行算法的精度和原有串行算法的精度基本一致,说明所采用的解除数据相关性的策略是可行的㊂表1 计算参数和结果T a b l e 1 C a l c u l a t i o n p a r a m e t e r s a n d r e s u l t sC e l l n u m b e r S p a t i a l s t e pT i m e s t e p F l u i dv e l o c i t yEP a r a l l e l a l go r i t h m S e r i a l a l go r i t h m 50ˑ50ˑ500.10.1u x =0.5;u y =0.50.1640.162905 第6期 马天宝等:三维多物质弹塑性流体动力学E u l e r 方法的并行算法研究及程序测试3.2 子区域边界网格(层)数量的确定 子区域边界网格是用来储存临近子区域相关变量的额外网格㊂采用E u l e r 数值方法,子区域间的关联性,即子区域边界网格的数量取决于两个因素:(1)因算法本身所固有的因素,当前网格物理量依赖于周边网格物理量的更新;(2)当前网格物理量的更新需要周边网格的物理量㊂考虑到通信所占用的开销,子区域边界网格(层)的数量应该考虑以上两个因素所增加的网格(层)数的最小值㊂以一维E u l e r 输运步为例,说明子区域边界网格数量的确定过程㊂如图1(a )所示,其中u 为当前k 网格的速度值㊂k 网格的更新(输运)影响k +1网格的物理量,也就是说k +1网格的物理量依赖于k 网格和k +1网格物理量的更新(输运),即k 网格为k +1网格的依赖网格,k +1网格为k 网格的影响网格㊂因此若想保证图1(a )左端1网格更新正确,需要在计算域左侧增加1层网格,并参与计算更新(输运),如图1(a )的灰色虚网格㊂由于在实际的更新运算中,往往需要周边网格物理量的信息,因此考虑如下的差分格式u n +1k =(1-2r )u n k +r (u n k +1+u n k -1)(3)则k 网格更新需要k ㊁k +1㊁k -1网格的物理量㊂综合以上两种因素,对于一维E u l e r 输运步的运算来说,需要在计算域增加3层子区域边界网格,左边增加2层,右边增加1层,如图1(b)所示,其中左边灰色虚网格参与运算,两端的虚网格不参与运算更新㊂(a )E x t r a c e l l a d d e d c o n s i d e r i n g t h e f i r s t f a c t o r (b )E x t r a c e l l s a d d e d c o n s i d e r i n g tw o f a c t o r s 图1 一维E u l e r 输运步子区域边界虚网格确定过程F i g.1 T h e p r o c e s s o f e x t r a c e l l s a d d e d i nE u l e rm o d e l 4 P MM I C -3D 并行程序测试G.J .M ye r s 在他的名著‘软件的测试技巧“一书中给出测试的定义: 程序测试是为了发现错误而执行程序的过程 [9]㊂基于M P I 并行程序的测试与调试的主要困难在于:除了串行程序的所有问题之外,并行程序还会有一些其它的问题存在㊂具体体现在:(1)并行算法设计的复杂性,对于一个串行程序很容易实现的问题,并行实现起来可能就困难重重,潜在的算法逻辑错误不易发现;(2)多进程之间协同作业(任务的映射与分发),导致并行程序比串行程序难于驾驭;(3)基于M P I 的并行程序增加了通信㊁进程同步等操作,增加了程序潜在的风险,一些不当的操作可能会导致程序异常中断或死锁㊂考虑到基于M P I 并行程序的测试与调试的困难,把对P MM I C -3D 并行程序的测试分为2个阶段:阶段1,不依赖于实际的物理模型㊁测试分区对计算结果的影响;阶段2,设计合理的物理模型,测试计算结果是否符合物理规律㊂阶段1是阶段2的基础,只有在阶段1测试正确的基础上才能进行阶段2的测试㊂(a )2D g e o m e t r i cm o d e l (b )3D m o d e l 图2 聚能装药二维几何模型及三维模型图F i g .2 2D g e o m e t r i cm o d e l a n d 3D m o d e l o f s h a p e d c h a r ge 聚能射流算例的计算几何模型如图2所示,其中图2(a )为二维尺寸结构图,图2(b )为三维模型图㊂药柱直径为60mm ,高为90mm ㊂依据网格步长设计了4个算例,见表2㊂聚能射流的数值模拟涉及到的材料包括炸药㊁金属药型罩及空气3种介质㊂炸药采用B 炸药,金属罩采用45钢㊂爆轰产物采用J W L 状态方程;空气采用理想气体状态方程;对于金属材料,考虑其在高温㊁高压㊁高应变率下表现的动态行为,采用M i e -G r ün e i s e n 状态方程描述㊂起爆方式采用点起爆,并采用简单的燃烧模型模拟爆轰波在炸药中的传播过程㊂015 高 压 物 理 学 报 第25卷表2 射流算例模型T a b l e 2 M o d e l o f s h a p e d c h a r g e je t M o d e lC h a r g em a s s /(g )S p a t i a l s t e p/(m mˑm mˑm m )C e l ln u m b e r M o d e lC h a r g em a s s /(g )S p a t i a l s t e p/(m mˑm mˑm m )C e l ln u m b e rE x a m p l e 13401.00ˑ1.00ˑ1.0080ˑ80ˑ180E x a m p l e 33420.04ˑ0.04ˑ0.04227ˑ227ˑ480E x a m pl e 23440.06ˑ0.06ˑ0.06151ˑ151ˑ300E x a m pl e 43420.03ˑ0.03ˑ0.03301ˑ301ˑ600对于射流算例,阶段1的测试采用射流头部速度作为监测变量,监测精度为10-11m /s,算例1的分区方式为:3ˑ1ˑ1,1ˑ3ˑ1,1ˑ1ˑ3;算例2的分区方式为:1ˑ4ˑ1,1ˑ1ˑ4;算例3的分区方式为:2ˑ2ˑ2,1ˑ2ˑ4;算例4的分区方式为:2ˑ2ˑ8,2ˑ2ˑ6㊂测试表明,不同分区方式下的计算结果保持一致,表明P MM I C -3D 并行程序在网格数较多的情况下,分区方式不影响计算结果(算例4网格数为5436万)㊂图3展示了4种算例的三维射流图,从图3中可以看出,算例2㊁算例3及算例4的射流形状明显优于算例1,算例3和算例4射流形状差别不大㊂图4为4种算例的射流头部速度随时间的变化曲线,从图4中可以看出,随着网格数的增多,射流最大头部速度略有增加,算例2㊁算例3及算例4比较接近,三者最大头部速度分别为4.4㊁4.5和4.7k m /s ,比实验结果略小,而算例1最大头部速度仅为4.0k m /s㊂图3 32.51μs 时4种算例的三维数值模拟图F i g .3 T h r e e -d i m e n s i o n a l s i m u l a t i o n g r a p h s a t 32.51μs o f f o u r e x a m pl e s 图4 4种算例的射流头部速度随时间变化曲线F i g .4 T h e j e t t i p v e l o c i t y c h a n ge sw i t h t i m e 5 P MM I C -3D 并行性能测试并行程序除了满足分区方式不影响计算结果的最基本要求外,加速比㊁并行效率及可扩放性也是衡量并行程序质量的主要性能指标㊂P MM I C -3D 程序并行性能测试在自主定制的八节点集群上进行,每个节点包含两颗I n t e l 四核E 5620C P U ,主节点24G 内存,其余节点12G 内存㊂5.1 加速比及并行效率测试并行系统的加速比(S p e e d u p)是指对于一个给定的应用,并行程序的执行速度相对于串行程序的执行速度加快了多少倍,也称为 绝对加速 (A b s o l u t eS p e e d u p );对于给定问题,同一程序在单C P U 的运行时间除以在多个C P U 运行的时间,称为 相对加速 (R e l a t i v eS p e e d u p)㊂加速比除以处理机个数,称为并行效率㊂加速比及并行效率测试采用第4节中的射流测试算例中的算例2㊁算例3作为测试算例,二者的网格数分别为684万和2473万㊂射流测试中的算例2㊁算例3的加速比及并行效率如图5所示㊂从图5中可以看出:加速比随着进程数的增加而增加,用64进程时算例2和算例3的加速比可达到16倍和20倍;并行效率随着进程数的增多而降低,用64进程时算例2和算例3的并行效率只有25%和30%,这是由于进程数增多,通信开销所占的比重也会增加;一样的进程数下,加速比及并行效率随着网格数的增多而增大㊂115 第6期 马天宝等:三维多物质弹塑性流体动力学E u l e r 方法的并行算法研究及程序测试图5 加速比及并行效率F i g .5 T h e s p e e d u p a n de f f i c i e n c y ve r s u s t h en u m b e r of p r o c e s s o r s 5.2可扩放性测试图6 可扩放性折线F i g .6 S c a l a b i l i t y可扩放性是指在确定的应用背景下,计算机系统(或算法或编程等)的性能随处理器的增加而按比例提高的能力㊂可扩放性是和并行算法以及并行计算机体系结构放在一起讨论的,某个算法在某个机器上的可扩放性反映了该算法是否能有效利用不断增加的C P U 的能力㊂采用射流算例作为测试算例,其结果如图6所示,计算域为10c mˑ10c mˑ20c m ㊂对于不同的进程数,每个进程分配的网格数固定为100万㊂测试的进程数如图6横轴所示,纵轴为计算1000个时间步的总时间㊂理想的可扩放性折线应该是一条平行于横轴的直线㊂由图6可知:随着进程数的增多,平均每个进程处理固定网格数的时间逐渐增加,64进程下处理时间是单进程下处理时间的2.8倍㊂5.3 影响并行性能的因素图7 通信时间在总时间的百分比F i g .7 T h e p e r c e n t a ge of c o mm u n i c a t i o n t i m e t o t h e t o t a l t i m e影响并行性能的因素主要包括通信所占用的开销以及负载不均衡所造成的同步等待时间等㊂图7为射流测试算例3在不同进程数下通信时间占总时间的百分比曲线㊂从图7中可以看出,通信时间在计算总时间的百分比随着进程数的增多而快速增大,64进程时达到48%㊂这是由于:(1)P MM I C -3D 中每一个子模块中都有诸如密度㊁质量㊁速度㊁能量及动量等变量更新,而这些变量的变化会影响下一个子模块的计算,因此需要数据通信;(2)每个节点采用千兆网络及千兆交换机连接,实测数据通信速度为35~45M ,网络通信延迟及带宽也影响通信㊂通信开销比重过大是影响并行性能的最重要因素㊂6 结 论(1)E u l e r 数值方法由于需要处理物质在网格间的输运问题,隐含在输运算法中的数据相关性及子区域间的关联性不易发现,因此并行算法设计中潜在的逻辑性错误不易发觉㊂215 高 压 物 理 学 报 第25卷(2)由于P MM I C -3D 中通信过多及自主定制的集群通信延迟及带宽的限制,通信开销比重过大,影响了P MM I C -3D 的并行性能㊂优化程序结构㊁减少通信量是改善P MM I C -3D 并行性能的重要手段㊂(3)从射流算例的模拟结果看,随着网格数的增多计算精度越来越高,数值模拟三维图片的分辨率越来越高,表明所采用的并行算法是合理的;P MM I C -3D 并行程序增大了计算规模,加快了计算速度,达到了并行程序设计的目的㊂R e f e r e n c e s:[1] D o n g a r r a J ,F o s t e r I ,F o xG ,e t a l .S o u r c e b o o k o f P a r a l l e l C o m p u t i n g [M ].T r a n s l a t e d b y MoZY ,C h e n J ,C a oXL .B e i j i n g :P u b l i s h i n g H o u s e o fE l e c t r o n i c s I n d u s t r y ,2005.(i nC h i n e s e )D o n ga r r a J ,F o s t e r I ,F o xG ,等.并行计算综论[M ].莫则尧,陈 军,曹小林,译.北京:电子工业出版社,2005.[2] M a l a r d J .M P I :A M e s s a g e -P a s s i n g I n t e r f a c eS t a n d a r d [R ].U K :T h eU n i v e r s i t y o fE d i nb u r g h ,1994.[3] X u eM Y ,N i n g JG.R e s e a rc ho n 3DE u l e r i a nN u m e r i c a l S i m u l a t i o n f o rE x p l o s i o n [J ].A c t aA r m a m e n t a r i i ,2006,27(6A ):129-132.(i nC h i n e s e)薛妙轶,宁建国.三维爆炸问题的E u l e r 数值方法研究[J ].兵工学报,2006,27(6A ):129-132.[4] C h e nGL ,A nH ,C h e nL .T h eA r c h i t e c t u r e o f P a r a l l e l C o m p u t e [M ].B e i j i n g :H i gh e rE d u c a t i o nP r e s s ,2004.(i nC h i n e s e )陈国良,安 虹,陈 崚,等.并行算法实践[M ].北京:高等教育出版社,2004.[5] M aTB ,W a n g C ,N i n g JG.M u l t i -M a t e r i a l E u l e r i a nF o r m u l a t i o n s a n dH y d r o c o d e f o r t h e S i m u l a t i o n o f E x p l o s i o n s [J ].C M E S :C o m p u t e rM o d e l i n g i nE n g i n e e r i n g &Sc i e n c e s ,2008,33(2):155-178.[6] N i n g JG ,C h e nL W.F u z z y I n t e r f a c eT r e a t m e n t i nE u l e r i a n M e t h od [J ].S c iC h i n aSe rE-E n g M a t e rS c i ,2004,47(5):550-568.[7] H a r v i eDJE ,F l e c t c h e rD F .A N e w V o l u m eo fF l u i dA d v e c t i o nA l g o r i t h m :T h eS t r e a m S c h e m e [J ].JC o m p u t P h y s ,2000,162:1-32.[8] H a r v i eDJE ,F l e c t c h e rDF .A N e w V o l u m e o fF l u i dA d v e c t i o nA l g o r i t h m :T h eD if i n e dD o n a t i ng R e g i o nS ch e m e [J ].I n t JN u m e rM e t hF l ,2001,35:151-172.[9] S hi JM.S o f t w a r eE n g i n e e r i n g -P r i n c i p l e ,M e t h o d a n dA p p l y [M ].B e ij i n g :H i gh e rE d u c a t i o nP r e s s ,2002.(i nC h i n e s e )史济民.软件工程 原理㊁方法与应用[M ].北京:高等教育出版社,2002.S t u d y o nP a r a l l e lA l g o r i t h mo fE u l e r i a n M e t h o d f o rT h r e e -D i m e n s i o n a l M u l t i -M a t e r i a l P l a s t i c -E l a s t i cH yd r o k i ne t i c s MA T i a n -B a o 1,F E IG u a n g -L e i 1,Z H A N G W e n -Y a o 2(1.S t a t eK e y L a b o r a t o r y o f E x p l o s i o nS c i e n c e a n dT e c h n o l o g y ,B e i j i n g I n s t i t u t e o f T e c h n o l o g y ,B e i j i n g 100081,C h i n a ;2.B e i j i n g L a b o r a t o r y o f I n t e l l i g e n t I n f o r m a t i o nT e c h n o l o g y ,S c h o o l o f Co m p u t e rS c i e n c e a n dT e c h n o l o g y ,B e i j i n g I n s t i t u t e o f T e c h n o l o g y ,B e i j i n g 100081,C h i n a )A b s t r a c t :P a r a l l e l c o m p u t i n g o f 3De x p l o s i o na n d s h o c k p r o c e s s e s o n t h e p a r a l l e l c o m pu t e r i s e f f e c t i v e m e a n s f o r t h e l a r g e -s c a l en u m e r i c a l s t u d y o f e x p l o s i o na n d s h o c k p r o c e s s .C o n s i d e r i n g t h e c o m p l e x i t y o f t h e p a r a l l e l p r o g r a mm i n g ,t h e o v e r a l l s t r a t e g y f o r p a r a l l e l p r o g r a mm i n g o f 3D m u t i l -m a t e r i a l h y d r o -e l a s t o p l a s t i ch y d r o c o d eMM I C -3Dw a s d i s c u s s e d ,a n d t h e P MM I C -3D p a r a l l e l h y d r o c o d ew a s d e s i gn e d b a s e do nM P I (M e s s a g e P a s s i n g I n t e r f a c e ).I n a d d i t i o n ,a p r a c t i c a l p l a n o f p r o g r a mt e s t i n g w a s p r e s e n -t e d .T h e s p e e d u p ,e f f i c i e n c y a n ds c a l a b i l i t y o f t h eP MM I C -3D p a r a l l e lh yd r o c o d ewe r e t e s t e do nt h e c l u s t e r c o n s i s t i n g of 8n o d e s b a s e d o n t h e n u m e r i c a l e x a m p l e o f s h a p e d c h a rg e j e t ,a n d th e e f f e c t o f t h e b o t t l e n e c k s o f P MM I C -3D p a r a l l e l h yd r o c o d ew a s d i s c u s se d .K e y wo r d s :e x p l o s i o na n d s h o c k ;p a r a l l e l c o m p u t i n g ;M e s s a g eP a s s i n g I n t e r f a c e (M P I );p r o g r a mt e s t 315 第6期 马天宝等:三维多物质弹塑性流体动力学E u l e r 方法的并行算法研究及程序测试。

0 前言
随着机器人技术不断向高速度、高精度、轻型 化和大跨度方向发展，柔性机器人动力学的研究越 来越受到重视。单机器人的操作比较有限，机器人 之间的协调操作可以弥补其不足。协调操作不可避 免地要求人们研究机器人与工作对象之间、机器人 与机器人之间、机器人内部各个杆件之间的碰撞及 其动力学响应问题，从而为机器人轨迹规划与实时 控制提供依据。碰撞是机器人操作过程中不可避免 的动载条件，柔性机器人的碰撞将产生系统速度的 突变， 碰撞过程产生的巨大冲击力将激发高阶模态， 从而影响系统的动态特性。此外，机器人杆件柔性 程度的不同，也对碰撞系统的动力学效应产生不同 的影响。总之，碰撞对机器人轨迹规划和实时控制 等带来严峻的挑战，国内外学者对柔性机器人碰撞 动力学给予了高度重视[1-5]。KHULIEF 等[1]运用广 义冲量－动量定理研究了柔性体的碰撞，求出碰撞 后系统的动力学响应。 CHAPNIK 等[2]对单臂机器人 碰撞动力学进行了研究，采用有限元法建模，建模 时考虑了多种阻尼效应，其理论计算结果与文中自 行设计的试验结果吻合得很好。刘锦阳等[3]采用子 结构法[4]， 由 JOURDAIN 变分原理导出了柔性机械 臂在碰撞前、碰撞阶段以及碰撞后的动力学方程， 它区别于文献 [1] 的显著特色在于能获得碰撞过程 中的动力学性态，并能求出撞击力。刘才山等[5]针 对做大范围回转运动的柔性梁与固定斜面碰撞的情
=− m ∆z
ˆ ∂V Ι ∂z
(9)
式中，m 为柔性机器人系统的广义质量矩阵，可以 为碰撞前后广义速度之 由现有的文献[6]得到， ∆z =z −z 0 。 差，即 ∆z 此处的关键问题是式(9)右端函数的计算。将式 (7)代入式(9)的右端，得 = RΙ m ∆z (10) Ι 式中 R 为广义冲量阵，其形式为

导弹发射动力学基础理论与数值解法一、弹架系统动力学基本任务和发展现状1.发射动力学基本任务研究弹架系统动态优化设计的基本理论和方法，其目的在于寻求合理和实用的分析计算方法，以保证动载作用下结构的安全、经济及使用性能，使火箭导弹的发射精度和可靠性符合要求。

2.发射动力学发展现状近年来，多体系统动力学经过了长足发展，从建模方法的研究进入更深层次问题的分析，比如结构柔性、间隙、接触碰撞和摩擦等影响动态分析精度的方面，为发射动力学问题的解决提供了更为准确有效的理论方法。

同时随着计算机技术的进步和有限元学科的发展，不仅刚体，弹塑性体的动力学分析也变得十分容易，利用计算机模拟进行发射动力学研究将成为必然趋势。

虚拟样机技术的使用将传统的产品设计循环过程以数字化方式进行，可以避免样机的重复构建，不仅能够缩短武器系统研发时间周期和降低开发成本，并且非常利于开展协同工作。

现在，Newton－Euler法、Lagrange法、Roberson－Wittenberg 法、Kane法和变分法是用来建立火箭导弹发射动力学模型的常用方法。

使用计算机对建立的发射动力学方程进行数值求解可以解决许多以前难以解决的问题。

美国等发达国家早已开发成熟的动力学软件，只需先建立动力学系统的三维实体模型，然后施加系统各部件的约束和外部激励，就可以直接计算并且输出系统的动态响应，避免了大量烦琐的计算机编程。

多学科间、多仿真软件间的联合仿真对发射动力学的进一步发展有着非常重要的意义。

目前，国内学者在发射动力学的理论基础和虚拟仿真方法上进行了深入研究，取得了极大的成果。

3.弹架系统动力学模型分类简介弹架系统的动力学模型可以按参数类型或参数的分布规律分类。

其按参数类型有物理模型、模态模型和响应模型三类。

物理模型是描述结构的物理参数（例如惯性、刚度、阻尼）之间的关系，同时用数学模型描述结构动态特性的数学方程，用空间模型，即几何模型，来描述结构动态特性的空间关系；模态模型是描述动态特性的模态参数（包括固有频率、振型和阻尼）；响应模型是描述结构输入和输出的传递关系，一般用频率响应函数矩阵表示。

械
工
程
学
报
第 48 卷第 19 期期
学进行了动力学建模理论、数值计算方法、接触碰 撞试验等方面的大量研究，并取得了一系列成果。 但由于碰撞过程的强非线性、高度耦合、数值计算 困难等复杂特性，目前对这些问题的研究还远未成 熟。对碰撞过程的正确处理是解决柔性多体系统碰 撞动力学问题的关键， 按照对碰撞过程假设的不同， 可以将多体系统碰撞动力学建模方法分为冲量动量 法、连续接触力法、接触约束法等几种方法，上述 几种方法各有优势和局限性[5]，目前国内外学者在 柔性多体系统碰撞动力学建模方面还没有形成统一 的意见。YIGIT 等[6]用动量平衡法结合恢复系数研 究了定轴转动柔性梁的碰撞问题；华卫江等[7]使用 广义冲量—动量方程和碰撞恢复系数方程推导出两 个柔性机器人系统发生碰撞的动力学方程。 KHULIEF 等 [8]提出了用线性弹簧阻尼模型近似模 拟多体系统碰撞过程的方法； 方建士等[9]采用 Hertz 接触理论和非线性阻尼理论建立接触—碰撞模型； 刘才山等[10]根据 Hertz 接触理论和线性切线接触刚 度建立碰撞接触模型和系统含碰撞的一致线性化的 动力学方程。WU 等[11]用子结构离散的方法研究了 两个柔性体的接触碰撞问题；刘锦阳等[12]在子结构 离散的柔性多体系统动力学基础上导出了描述柔性 机械臂三个不同阶段的动力学方程；董富祥等[13]建 立了基于变约束的柔性体碰撞动力学方程，结合连 续介质力学间断面理论给出碰撞瞬时速度跳跃公 式。 KIM 等[14]在处理两个柔性多体系统的接触问题 时采用接触区域局部有限元和其余区域传统建模相 结合的方法。石明全等[15]针对多接触动力学系统计 算效率问题，提出了一种提高仿真效率的方法。 本文基于柔性多体系统刚柔耦合动力学理论， 建立了作大范围回转运动柔性悬臂梁无碰撞时的刚 柔耦合动力学方程。提出了一种碰撞求解方法，即 初始碰撞时刻首先利用冲量动量法进行一次求解以 得到碰撞初始条件，之后使用接触约束法求解，导 出了系统的碰撞动力学方程。给出了接触、分离判 据，实现了系统在碰撞前、碰撞过程、碰撞后各阶 段的动力学转换与求解。编制仿真软件，对实例进 行了动力学仿真，得到了系统的全局动力学响应。

西北工业大学博士学位论文柔性撞击系统的建模、精细算法及控制研究姓名：赵玉立申请学位级别：博士专业：一般力学与力学基础指导教师：邓子辰20031106摘要本文针对柔性撞击系统的建模、数值模拟、控制及精细积分理论在其中的应用进行了研究。

对于柔性撞击系统的建模问题，分别采用有限元思想和模态法建立了作回转运动的柔性梁动力学模型，然后又分别利用约束方程和Ｈｅｒｔｚ接触定律给出了其与固定斜面发生接触的撞击模型，指出它们的动力学方程分别为微分／代数混合方程及复杂的非线性动力学方程。

同时利用Ｈｅｒｔｚ接触定律建立了球一弹性体撞击模型，其方程为含积分项的非线性动力学方程。

运用指数矩阵的２”类算法，构造了微分方程的精细积分算法，并分析了计算精度。

将之应用于刚性方程及非线性方程的求解中，并对非线性方程的求解从显式和隐式两个方面给出了改进措施，得出：使用精细积分法可以给出某些微分方程的精确解，且绝对稳定。

对于数值模拟问题，针对不同的动力学模型，分别利用精细Ｂａｕｍｇａｒｔｅ违约修正法及状态方程的精缅积分法对系统进行了数值模拟，得到了满意的数值解。

对于控制问题，建立了简化柔性梁与固定斜面发生撞击时的动力学控制方程，并利用线性二次最优（ＬＱ）控制策略对之进行了控制研究。

将撞击看作为系统的一种扰动，利用Ｈ。

控制理论框架下的扰动抑制问题理论，建立了简单柔性撞击系统的控制模型。

两种控制方法的要点是求解矩阵黎卡提方程，所以还重点讨论了基于结构力学与最优控制模拟关系的黎卡提方程２”类精细算法。

数值模拟结果表明，文中的控制策略是可行的，而且可通过文中的方法可得到矩阵黎卡提方程的精确解（在计算机有效精度范围内）。

关键词：柔性撞击系统，精细积分法，ＬＱ控制，黎卡提方程，Ⅳ。

控制ＡＢＳＴＲＡＣＴＩｎｔｈｅｐｒｅｓｅｎｔｐａｐｅｒ，ｍｏｄｅｌｉｎｇ，ｎｕｍｅｒｉｃａｌｓｉｍｕｌａｔｉｏｎ，ｃｏｎｔｒｏｌｏｆｆｌｅｘｉｂｌｅｉｍｐａｃｔｓｙｓｔｅｍａｎｄａｐｐｌｉｃａｔｉｏｎｓｏｆｔｈｅｐｒｅｃｉｓｅｉｎｔｅｇｒａｔｉｏｎｍｅｔｈｏｄｗｅｒｅｓｔｕｄｉｅｄＴｏｍｏｄｅｌｔｈｅｆｌｅｘｉｂｌｅｉｍｐａｃｔｓｙｓｔｅｍ，ｆｉｎｉｔｅｅｌｅｍｅｎｔｔｈｅｏｒｙａｎｄｍｏｄｅｍｅｔｈｏｄｗｅｒｅｒｅｓｐｅｃｔｉｖｅｌｙｕｓｅｄｔｏｅｓｔａｂｌｉｓｈｔｈｅｄｙｕａｍｉｃｅｑｕａｔｉｏｎｏｆｔｈｅｒｏｔａｔｉｎｇｆｌｅｘｉｂｌｅｂｅａｍ，ｔｈｅｎｉｍｐａｃｔｍｏｄｅｌｓｗｅ∞ｂｕｉｌｔａｎｄｕｓｉｎｇｃｏｎｓｔｒａｉｎｔｅｑｕａｔｉｏｎａｎｄＨｅｒｔｚｃｏｎｔａｃｔｔｈｅｏｒｙ，ｔｈｅｙｗｅｒｅｄｉｆｆｅｒｅｎｔｉａｌ／ａｌｇｅｂｒａｉｃｅｑｕａｔｉｏｎｃｏｍｐｌｉｃａｔｅｄｎｏｎｌｉｎｅａｒｅｑｕａｔｉｏｎｒｅｓｐｅｃｔｉｖｅｌｙ．Ａｔｔｈｅｓｓｎｌｅｔｉｍｅ，ｔｈｅｍｏｄｅｌｏｆｂａｌｌ—ｆｌｅｘｉｂｌｅｂｏｄｙｉｍｐａｃｔｓｙｓｔｅｍｗａｓｂｕｉｌｔｕｓｉｎｇＨｅｒｔｚｃｏｎｔａｃｔｔｈｅｏｒｙ，ｉｔｓｅｑｕａｔｉｏｎｗａｓｎｏｎｌｉｎｅａｒｉｎｔｅｇｒａｌｅｑｕａｔｉｏｎ．Ｔｈｅｐｒｅｃｉｓｅｉｎｔｅｇｒａｔｉｏｎｍｅｔｈｏｄｏｆｄｉｆｆｅｒｅｎｔｉａｌｅｑｕａｔｉｏｎｗａｓｐｒｅｓｅｎｔｅｄｂｙｕｓｉｎｇ２“ａｌｇｏｒｉｔｈｍ，ｉｔｓｎｕｍｅｒｉｃａｌｐｒｅｃｉｓｉｏｎｗａｓａｌｓｏａｎａｌｙｚｅｄ。

Ｍａ ｒ ．２ ０１ ３
变 拓 扑 柔 性 多体 系统 接 触 碰 撞 动 力 学研 究 木
洪嘉振 刘铸永
２ ０ ０ ２ ４ ０ ） （ 上海交通大学工程力学系 ， 上海
摘要
在 实际工程领域 中存在着大量接触碰撞 等非连续动力学 问题 ， 现有 的解决柔 性多体 系统连续动力 学
过程 的建模理论 与方法 ， 已经无 法解决或无法很好解决这些 问题． 本 文基于变 拓扑思想 ， 提 出了附加接触 约 束的柔性多体 系统 碰撞 动力学建模理论 ； 通过设 计柔 性 圆柱 杆接 触碰撞 实验 ， 验证 了所提 出附加 约束接触 碰撞模型 的有效性 ； 针对柔性 多体 系统全局动力 学仿 真面 临时间 和空间 的多尺度 问题 ， 提 出多变量 的离散 方法 ， 从而提 高了柔性多体系统非连续 动力 学的仿真效率． 关键词 柔性多体 系统 ， 接触碰撞 ， 变拓扑 ， 数值仿真 ， 实验研究
２ ０ １ ２ － ０ ４ — １ １ 收到第 １ 稿， ２ ０ １ ２ － ０ ６ － １ ５收到修改稿． 国家 自然科学基金重点项 目（ １ １ １ ３ ２ ０ ０ ７ ， １ １ ２ ０ ２ １ ２ ６ ） ， 博士点基金 （ ２ ０ ０ ９ ０ ０ ７ ３ １ １ ０ ０ ０ ９ ） 、 上海市 自然科学基金（ １ １ Ｚ Ｒ １ ４ １ ７ ０ ０ ０ ） 资助 十 通讯作者 Ｅ — ｍ ａ ｉ ｌ ： ｚ ｈ ｕ ｙ ｏ ｎ ｇ ｌ ｉ ｕ ＠ｓ ｊ ｔ ｕ ． ｅ ｄ ｕ ． ｃ ｎ
动
力
学
与
控
ห้องสมุดไป่ตู้
制
学
报
２ ０ １ ３年第 １ １ 卷
建模 ， 但柔 性体 的大 范 围空间运 动及 弹性 变形 不 加 区分地采 用惯 性坐 标系下 的节 点坐标 表 示 ， 由此 得 到 的动力 学方 程 的广义 坐标 数 目非 常 庞大 ， 计算 效 率低 ， 难 以应 用 到 复杂 的柔 性 系统 动力 学 分 析 中． 同时 由于该 方法无 法 区分 刚体运 动和 弹性 变 形 ， 给 控制 系统 的设 计带 来 难 以逾 越 的 困难 ． 因此 ， 迫 切 需要 提 出一种 柔 性 多体 系统 碰撞 动 力 学 过 程 全 局 仿 真方法 ， 以满 足柔 性多体 系 统碰撞 动力 学 对通 用

总之，考虑刚柔耦合效应的柔性多体系统碰撞动力学研究具有重要的理论和 实践意义，将在机器人技术、航天技术等多个领域发挥重要作用。
感谢观看
然而，本次演示的研究仍存在一些不足之处。首先，实验对象仅为特定类型 的柔性机器人，其结论是否适用于其他类型的柔性多体系统仍需进一步验证。其 次，本次演示的理论模型简化了一些复杂因素，如摩擦力、气体阻力等，这些因 素在某些应用场景下可能对碰撞过程产生重要影响。
展望未来，我们提出以下研究方向：
1、研究不同类型和结构的柔性多体系统撞动力学的影响；
刚柔耦合效应主要通过引入弹性变形、阻尼等因素，考虑物体的变形和能量 耗散对碰撞过程的影响。在刚柔耦合模型中，碰撞过程中物体的速度、位移和加 速度不仅与物体的质量、惯性和碰撞速度有关，还与物体的变形、阻尼等因素有 关。这些因素的存在增加了碰撞过程的复杂性和不确定性，使得对柔性多体系统 碰撞动力的研究更具挑战性。
柔性多体系统碰撞动力的实验研 究
为了深入了解刚柔耦合效应对柔性多体系统碰撞动力学的影响，我们设计了 一系列实验进行研究。实验中使用了多个柔性关节机器人作为研究对象，通过高 速摄像机和力传感器等设备采集了碰撞过程中物体的速度、位移和作用力等数据。
实验结果表明，在碰撞过程中，柔性关节机器人的碰撞速度、碰撞时间和作 用力均小于刚性机器人。这是由于柔性机器人在碰撞过程中发生了明显的变形和 阻尼耗散，使得其碰撞过程更加缓慢和柔和。此外，实验结果还显示，随着碰撞 速度的增加，柔性机器人的变形量和阻尼耗散也相应增加，从而减小了碰撞冲击 力。
2、考虑更复杂的外部条件和因素，如不同温度、湿度和气压等环境条件对 柔性多体系统碰撞过程的影响；
3、研究具有更复杂刚柔耦合特性的柔性多体系统碰撞动力学模型和方法， 以提高模型的准确性和适用性；

力学响应的主流方法 [２]．
１９９７ 年,美国 多 体 动 力 学 专 家 Shabana 教 授
为了实现柔性体大范围运动和大变形的统一描述,
提出了绝 对 节 点 坐 标 方 法 (
Abs
o
l
u
t
eNoda
lCoo
r
Ｇ
[]
该方法基于连续介
d
i
na
t
eFo
rmu
l
a
t
i
on,ANCF)３ ．
质力学的构型描述,在惯性坐标系中定义单元节点
随着航天技术的发展,轻质、高精度、柔性零部
件被广泛应用于众多工程领域,如空间站柔性操作
机械臂、卫星环形桁架天线、航天器太阳帆板、空间
格生成和细化的过程,提高数值仿真计算效率和精
度,实现模型设计和分析一体化 ．
机器人等 ．
这些系统中存在 刚 性 部 件 与 柔 性 部 件 的
大范围运动相互耦合,而且柔性部件的大范围运动
王悦等 [１７]提出了一种基于 T
技术解决了短 粗 梁 结 构 分 析 中 常 见 的 剪 切 和 薄 膜
对受冲击柔性薄板的 动 力 学 分 析 表 明 所 提 出 T 样
集成到现有的 多 体 系 统 动 力 学 分 析 软 件 时 遇 到 不
少困难 ．
网格模型采用同一表达方式,避免了网格细化过程
限元中网格 划 分 的 复 杂 过 程 ．
此 外,和 绝 对 节 点 坐
标方法相比,等 几 何 分 析 的 结 果 精 度 更 高,其 分 析
在柔性多体系统动力学研究中,采用绝对节点
Hughe
s 等 ５ 为 了 解 决 这 一 问 题,提
点转角坐标,得到的柔性多体系统微分 － 代数方程

首先，本文利用第一类 Lagrange 方程建立了柔性多体系统的动力学方 程，并建立了速度意义下的约束方程，方程最终表现形式为微分/代数混合 方程组。在推导的过程中未对方程作任何简化以体现系统得耦合非线性特 征。本文以后的工作都是在此基础上展开的。
其次，本文基于 ADAMS 详细的建立了包含同步展开机构和展开驱动 机构的 DFH-4 卫星系统的动力学模型，并基于 ADAMS 利用建立的卫星系 统动力学模型对卫星太阳翼展开过程进行了详细的动力学仿真研究。
Keywords spacecraft; flexible multi-body dynamics; flexible appendages; deployment process; attitude motion; fault mode
- II -
哈尔滨工业大学工学硕士学位论文
目录
摘要 ............................................................................................................................... I Abstract ....................................................................................................................... II
答 辩 日 期： 2006 年 6 月
授予学位单位： 哈尔滨工业大学
Classified Index：V411.8 U.D.C.: 629.7
Dissertation for the Master Degree in Engineering

参 数（ 变量改为斜体 ）
节距 Ｐ
依据实际需要 ， 合理选用 的链节型号为 ３ Ｂ 整 ２，
尺寸／ ｍ ｍ
分 度 圆 ：＝ ／ｎ ＇Ｚ）３６９ ； ｄｐｓ （ｒ ＝ ５．６ ｉ ｒ ／ 齿 顶 圆 ：ａ ｄ １ Ｐ ｄ， ｄ ：＋ ． － １ ２
ｄ ｒ＝ ＋ １ １ ） ｄ． ａ ｄ （－ ． １ ￣ ６
模型． Ａ Ａ Ｓ环境下 ， 在 ＤＭ 使用宏命 令功能完成虚拟样机 中链 节 、 柔性 铰和接触碰撞的 自动化建 模 ； 合理设置 弹链 系统 的仿 真参数 ， 完成 弹链 的动力学仿 真． 合实 际对 得到 的相关 曲线做简要 分析 ， 出相应 的仿真结论 ， 结 得 验证 了基于柔性 铰弹链 模 型的合理性与正确性 ， 以对 弹链 在运动过程 中的振动 和冲击问题进行很好 的预测 和分析． 可 关键词 ： 弹链 ； 柔性铰 ； 宏命令 ； 动力学仿真
Ａ ｓ ａ ｔＣ ｍｐｅ ｅ ｓ ｅ ｏ ｓ ｅｅ ｉ ｅｌ ｅｒ ｉ ｌｃｍｅｔ ｎ ｎ ｕａ ｉ ｌｃｍｅ ｔ ｆｈ ｉｔ ｂｔ ｅｎ ｂ ｔ ｃ ： ｏ ｒｈ ｎ ｉ ｌ ｃｎ ｉ ｒｄｗｔ ｔ ｎ ａ ｓ ａ ｅ ｎ ｄａ ｇ ｌｒ ｓ ａ ｅ ｎ ｏ ｔｅ ｏｎｓ ｅ ｅ ｒ ｖｙ ｄ ｈｈ ｉ ｄｐ ａ ｄｐ ｊ ｗ
ＨＵＳ ｅｇ ａ， Ｈ ＮＧＸ ａｄ ｎ ， Ｉ ａ ｈ ｎ Ｌ ｉ， Ｉ ｉ ｅ ｈｎ ｈｉ Ｚ Ａ ｉｏｏ ｇ ＬＵＺ ｉ ｕ ， Ｉ ｎ Ｌ ｙ ｃ Ｌ Ｑ ｕ
（ ｏｌｇ ｆＭｅ ｈ ｎｃ ｌａ ｄ Ｅ ｅｔ ｃ ｌＥｎ ｉｅ ｒ ｇ Ｈａｂｎ Ｅ ｇｎ ｅｉｇ Ｕｎｖ ｒｉ ，Ｈａｂｎ ５ ０ ， ｈｎ ） Ｃ ｌ ｅ ｏ ｃ ａ ｉａ ｎ ｌｃｒ ａ ｇｎ ｅ ｎ ， ｒ ｉ ｎ ｉｅ ｒ ｉｅｓ ｙ ｅ ｉ ｉ ｎ ｔ ｒｉ １０ ０１ Ｃ ｉａ

关键词：火箭弹，定向管，接触碰撞，动力学仿真
中图分类号：ＴＰ３９１．９
文献标识码：Ａ
Ｒｏｃｋｅｔ－ ｔｕｂｅ Ｃｏｕｐｌｉｎｇ Ｐｒｏｂｌｅｍ Ｂａｓｅｄ ｏｎ Ｃｏｎｔａｃｔ－ ｉｍｐａｃｔ Ｍｏｄｅｌ
ＴＡＮＧ Ｈｏｎｇ，ＬＥ Ｇｕｉ－ｇａｏ，ＪＩＡＮＧ Ｙｏｎｇ，ＣＨＥＮ Ｘｕｅ－ｙｉ （School of Mechanical Engineering，Nanjing University of Science and Technology，Nanjing ２１００９４，China）
Ａｂｓｔｒ ａｃｔ：Based on elastic plastic mechanics, multi body dynamics and finite element method of contact impact, two contact impact finite element models of rocket directional tube coupling are established. One model is rigid constraint, other is elastic plastic constraint. The dynamic response process of the two contact impact models are simulated by nonlinear dynamics. The Influence of rigid constraint on the initial disturbance is greater by analyzing the simulation results of two models. The contact impact process between rocket and directional tube is described with the elastic plastic constraint rocket directional tube system. The results provide valuable preferences for design and analysis of rocket weapon.

多柔体系统碰撞动力学研究综述
刘才山;陈滨
【期刊名称】《力学进展》
【年(卷),期】2000(030)001
【摘要】多柔体系统碰撞动力学研究具有重要的研究价值和工程实际意义,本文针对多柔体系统碰撞动力学研究中的几个基本问题进行了全面的分析和评述,其中包括多柔体系统动力学方程的描述、碰撞模型的建立、铰接间隙引起的碰撞问题、数值算法、实验研究、控制等几个方面,并根据目前的发展现状和研究中存在的问题,指出了今后多柔体系统动力学碰撞研究中的发展方向.
【总页数】8页(P7-14)
【作者】刘才山;陈滨
【作者单位】北京大学力学与工程科学系,北京,100871;北京大学力学与工程科学系,北京,100871
【正文语种】中文
【中图分类】O3
【相关文献】
1.多柔体系统动力学碰撞研究中的若干基本问题 [J], 刘才山;陈滨
2.多柔体系统动力学中的间隙接触内碰撞 [J], 陈鹿民;阎绍泽;郭峰;金德闻
3.多柔体系统动力学中的间隙接触内碰撞 [J], 陈鹿民;阎绍泽;郭峰;金德闻
4.多体系统动力学碰撞问题研究综述 [J], 董富祥;洪嘉振
5.基于刚-柔耦合多柔体动力学的悬架系统分析及优化 [J], 何通俊;纪玉国
因版权原因，仅展示原文概要，查看原文内容请购买。

变拓扑柔性多体系统接触碰撞动力学研究洪嘉振;刘铸永【摘要】There is an increasing number of discontinuous dynamics problems,such as contact and impact,in the engineering field.The known modelling theories and methods,which are used to solve the continuous dynamics of flexible multibody systems,are not able to solve those problems or to solve those problems perfectly.In this paper,an appending constraint method was proposed to solve the contact/impact dynamics of flexible multibody system based on the idea of topology variable.The impact experiment of flexible bars was designed to verify the validity of the proposed method.At last,a multi-variable method was proposed to solve the multi-scale problems of time and space when the global simulation of the discontinuous dynamics of flexible multibody systems was carried out,and the computational efficiency of the global simulation of the discontinuous dynamics of flexible multibody systems was improved by using this method.%在实际工程领域中存在着大量接触碰撞等非连续动力学问题,现有的解决柔性多体系统连续动力学过程的建模理论与方法,已经无法解决或无法很好解决这些问题.本文基于变拓扑思想,提出了附加接触约束的柔性多体系统碰撞动力学建模理论；通过设计柔性圆柱杆接触碰撞实验,验证了所提出附加约束接触碰撞模型的有效性；针对柔性多体系统全局动力学仿真面临时间和空间的多尺度问题,提出多变量的离散方法,从而提高了柔性多体系统非连续动力学的仿真效率.【期刊名称】《动力学与控制学报》【年(卷),期】2013(011)001【总页数】7页(P5-11)【关键词】柔性多体系统;接触碰撞;变拓扑;数值仿真;实验研究【作者】洪嘉振;刘铸永【作者单位】上海交通大学工程力学系,上海200240;上海交通大学工程力学系,上海200240【正文语种】中文近三十年来国内外学者在多体系统连续过程建模理论、计算方法和实验研究等方面做了大量的工作，取得了一系列重要科研成果，相继解决了动力刚化、微分代数方程数值求解、最优控制等重要问题，建立了程式化的多体系统刚柔耦合动力学建模理论和有效的数值求解算法，并形成了诸如ADAMS、DADS和RECDYN这样的大型工程软件，这些软件已经在众多工程问题的实践中经受了考验，大大提高了新产品的设计开发的效率．然而在现代航空航天、机械制造等领域，高速、轻质的系统面临的动力学环境与过程越来越复杂．柔性多体系统动力学中，现有的解决连续动力学过程的建模理论与方法，已经无法解决或无法很好解决这些问题［1－6］．例如，目前我国载人航天与未来空间站的航天器交会对接;现代航天器已经成为大型空间的机械系统，各类卫星太阳电池阵的展开，到位锁定时的接触碰撞;空间大型柔性部件(机械臂、太阳电池阵、大型天线阵等)的非连续驱动;铰链间隙引起的撞击等等．在这些非连续动力学过程中，系统的拓扑会发生突变，以其作用时间短、撞击力大给系统造成剧烈扰动，引起系统动力学性态的突变，对航天器机构动力学性态的灾难性影响已经引起了关注．柔性多体系统非连续动力学过程的实质是若干不同拓扑系统连续过程的不断切换．也就是具有某个自由度的某种连接关系的一种机构的运动，突然改变为具有不同自由度的另一种机构的运动．当前最典型的非连续问题是接触碰撞．多体系统接触碰撞问题在力学本质上是一种非光滑、变边界的高瞬态、高非线性动力学过程．对这类问题的研究既要考虑柔性体变形运动与大范围运动之间的耦合，又要考虑碰撞动力学过程的合理建模，其中对碰撞过程的正确处理是解决多体系统碰撞动力学问题的关键．现有的接触碰撞模型都还有不足之处．基于刚体假设的冲量动量方法虽然具有计算效率高、直观等优点，且对碰撞过程中出现的能量不协调等问题进行了很好的解决，仍然难以适应以上工程领域对碰撞力及系统碰撞动力学响应求解提出的需求．基于弹簧阻尼的连续碰撞力公式的出现在一定程度上克服了碰撞力难以求解的困难，并在工程领域中得到了相当广泛的应用，长期以来也为两个问题所困扰:碰撞力模型合理形式的选择和碰撞力方程中参数的合理辨识．这两个问题直接决定了仿真结果的可信性．目前存在的实验研究成果和方案往往只针对具体的个案，缺乏普适性．非线性有限元法、无网格等方法的优点是可以对碰撞动力学过程进行建模，但柔性体的大范围空间运动及弹性变形不加区分地采用惯性坐标系下的节点坐标表示，由此得到的动力学方程的广义坐标数目非常庞大，计算效率低，难以应用到复杂的柔性系统动力学分析中．同时由于该方法无法区分刚体运动和弹性变形，给控制系统的设计带来难以逾越的困难．因此，迫切需要提出一种柔性多体系统碰撞动力学过程全局仿真方法，以满足柔性多体系统碰撞动力学对通用性、高效率、准确性和可靠性的要求．本文基于非连续动力学过程通过合理切换可分解为若干连续过程的思路，将复杂机械系统变拓扑过程统一在通用的柔性多体系统模型下，研究通用与高效的柔性多体系统刚柔耦合动力学模型与接触碰撞动力学模型，解决变拓扑柔性动力学全局仿真系统切换的关键问题．采用附加约束的方法能够真正从物理上对接触碰撞问题进行准确描述，从而实现柔性多体动力学全局仿真．通过数值仿真和实验研究验证本文所提出的附加约束方法的有效性．针对变拓扑柔性多体系统动力学的非连续、非光滑特性和存在时间及空间的多尺度问题，提出多变量方法，研究了柔性多体系统非连续动力学全局仿真计算方法．基于弹簧阻尼的多体系统连续碰撞力模型目前在工程领域得到广泛应用．然而该方法将碰撞力作为物体间嵌入量和嵌入速度的函数来计算，且需要人工设置函数中的参数．仿真结果过分依赖于人工参数的选择，况且嵌入量本身也不符合物理实际［1－2］．如图1所示，本文将柔性多体碰撞分为碰撞前、碰撞阶段、碰撞后三个过程，采用附加接触约束描述碰撞阶段的动力学过程，将上述三个过程的动力学方程形式统一．通过引入碰撞检测机制，实现拓扑模型的自动切换．将此方法称为附加接触约束方法．柔性多体系统碰撞前和碰撞后动力学方程为［3］:这里M是柔性多体系统质量矩阵，Φq是多体系统铰约束的雅克比矩阵，λ是与铰约束对应的拉格朗日乘子矢量，w，fo和fu分别是柔性体系统受到的惯性力矢量，广义外力矢量和广义变形力矢量，γ是加速度约束方程的右端项．碰撞阶段，接触位置位移相同，故附加接触位移约束方程为:将其写成加速度约束方程的形式，即结合式(1)，柔性多体系统碰撞阶段动力学方程为:和γC是碰撞过程中附加接触约束方程的左端项和右端项矩阵．如图2所示，水平面上柔性杆B1与柔性杆B2的中心处发生碰撞，B1和B2的轴线相互垂直．B1为自由杆，B2两端受固定端约束．碰撞发生时B1的速度为0．2m/s，B2处于静止状态．杆的几何和材料参数见表1，B1和B2扁平等厚，故该问题可简化为平面问题处理．对此算例分别用ACCM法和LS－DYNA罚函数法进行仿真．使用计算机的配置参数为 AMD Athlon(tm)64 processor 3200+，2．01GHz，1GB 的内存．在利用LS－DYNA罚函数方法仿真时，需要人工设定罚因子．对同一个罚因子，计算步长不同，计算结果也不同．随着计算步长的减小，计算结果将趋于一个稳定解．对于本算例，表2列出了利用LS－DYNA罚函数方法仿真时，对于不同的罚因子获得稳定解所需要的最大时间步长，以及相应的CPU仿真时间．利用本文提出的方法，克服了这种对罚因子的依赖．此外，从图3－6的罚函数方法计算结果(虚线)发现，利用罚函数方法，对于不同的罚因子(依次分别对应罚因子为 0．01、0．1、1．0 和 5．0 时的情况)，所得到的稳定解是不同的．过小的罚因子(图3)会造成仿真结果严重失真．当罚因子增大到一定程度，仿真结果将变化不大(图5－图6)，可以认为仿真结果收敛于“精确解”．比较计算结果，本算例能得到“精确解”的罚因子应大于1．0．理论上罚因子，得到的解为精确解．但实际计算中，罚因子取值越大，为了得到稳定结果，则计算步长要取得非常小，严重影响了计算效率．可见，为了兼顾计算精度和计算效率，罚因子取值不能无限制大，必须要经过多次试算来选择合适的罚因子．从图3－图6的ACM方法计算结果(实线)发现，这种计算“精确解”依赖于罚因子的不足，对于本文提出的方法不存在．为了验证本文提出的碰撞动力学过程建模方法的正确性，设计了两圆柱铝杆碰撞实验．如图7所示为实验装置示意图．该实验装置包含以下部件:两圆柱截面铝杆，Polytech Gmbh公司制造的OFV－300F型激光测振仪(LDV)、动态应变仪、1．5V直流电源和金属实验台架．碰撞过程中杆2的速度由激光测振仪进行测量，它由激光扫描头、控制器、连接箱和工作站组成．如图7所示，一开始将杆B1从限制板上某一确定位置处释放．然后杆B1自由下摆，当下摆到水平位置处时，B1杆达到最大速度，即碰撞前初始速度，碰撞前B1的速度通过激光测振仪进行测量．经过25次测量，得到碰撞时刻杆 B1的速度为0．21m/s，其中速度平均误差小于0．5%．在碰撞实验中碰撞力的时间历程是一个关注点．在一些实验中，直接通过压电传感器测量碰撞力变化历程然而在实验过程中压电传感器必须固定在物体碰撞表面上，这样原有的碰撞界面就会受到损坏，导致压电传感器测量得到的碰撞力并非原来的碰撞力．因此，只能通过速度或者应变间接地验证出仿真结果中碰撞力的正确性．在本部分，将采用不用的建模方法对两杆碰撞动力学过程进行数值仿真，其中碰撞过程中杆的几何和运动学参数如表3所示．如图8－11所示，可以看出，本文提出的附加约束方法仿真结果与实验结果较为吻合．而弹簧阻尼模型仿真结果与实验结果之间存在着很大的误差．如图10所示，为ACM方法仿真结果与实验结果在杆B2上不同点处速度曲线比较，其中杆B2上11个测量点沿杆轴向方向均匀分布．由于图中存在22条曲线，限于空间，没有给出每条曲线的图标．从曲线的比较可以看出，几乎所有仿真结果与实验结果都吻合地很好，这意味着由ACM方法仿真得到的碰撞过程与真实的碰撞过程相吻合．从图11的仿真结果上看，可以看出在碰撞过程中不同位置处的速度并不同时改变，这说明刚体碰撞假定中各点速度同时突然变化的假定与实际情况之间存在着一定得差异．根据如图7所示实验结果，同一速度条件下100cm处和0cm处的时间差为1．9561e－004s，通过实验数据可以得到弹性波在铝杆表面的传播速度为5．1122e+003m/s，而一维情况下弹性波的理论传播速度是5．0918e+003，相对误差是0．4006%．这说明碰撞引起的速度扰动在弹性介质中的确以弹性波的形式进行传播，同时验证了碰撞动力学方程的速度跳跃条件的解释符合真实情况．如图8所示为弹簧阻尼模型和ACM方法仿真得到碰撞力时间历程比较．由于碰撞力难以通过实验的方法直接精确测量，因此只能通过仿真计算得到碰撞力的时间历程．从图8－11中ACM仿真结果和实验结果的对比上看，可以得出由ACM 方法仿真得到碰撞过程与实际的物理碰撞过程一致，这样就可以得出由ACM方法仿真得到的碰撞力与真实碰撞力一致．此外，从图中可以看出这样一个复杂的碰撞力变化历程很难通过简单的弹簧阻尼模型进行建模．为提高全局仿真的计算效率，提出了多变量离散方法．如图12所示，两多体系统1和2发生接触/碰撞．其中，系统1包含m个物体，系统2包含n个物体，发生碰撞的物体为系统1中的第i个物体Bi与系统2中的第j个物体Bj．对发生碰撞的物体Bi、Bj进行分区，接触局部区域采用非线性有限元方法，其他区域用浮动参考系方法．对于接触区域，以系统1上的接触域为从物体，系统2上的接触域为主物体．图13所示铅垂面内双摆由B1、B2两个物体组成．B1为刚体，长为0．5m，质量为5kg，对质心的转动惯量为0．15kgm2．B2为平面弹性体，表4给出了B2的材料参数，同时给出了双摆的初始位置和速度．双摆自水平位置已一定初速度下落，撞击刚性地面后又反弹．其中，表面与地面静摩擦系数为0．7，动摩擦系数为0．6．B2物体可能发生碰撞局部的区域采用有限元节点坐标，其它区域采用浮动基位形坐标和模态坐标，两部分交界如图13所示，浮动基原点取在O2点，横轴附着于B2两节点上．本算例中，对于碰撞过程采用显式法求解，对于未接触和接触分离后的阶段采用隐式法求解．隐式法时间步长取h=5E－3s，未接触阶段和接触分离后阶段分别只用到了19和30个时间步．显式法时间步长取0．2E－6s，碰撞阶段用到了3000个时间步．本例约有120个自由度(B1有3个自由度，B2浮动基3个自由度+20阶模态+47有限元节点)，三个阶段的仿真计算均在半分钟至一分钟内完成(CPU:AMD Athlon(tm)64 Processor 3200+，2．01GHz)．说明多变量方法具有较高的计算效率．本文基于非连续动力学过程通过合理切换可分解为若干连续过程的思路，将复杂机械系统变拓扑过程统一在通用的柔性多体系统模型下，研究通用与高效的柔性多体系统刚柔耦合动力学模型与接触碰撞动力学模型，解决变拓扑柔性动力学全局仿真系统切换的关键问题．基于变拓扑思想，本文提出了附加接触约束的柔性多体系统碰撞动力学建模理论;通过设计柔性圆柱杆接触碰撞实验，验证了所提出附加约束接触碰撞模型的有效性;针对柔性多体系统全局动力学仿真面临时间和空间的多尺度问题，提出多变量的离散方法，从而提高了柔性多体系统非连续动力学的仿真效率．2012－04－11 收到第 1 稿，2012－06－15 收到修改稿．【相关文献】1 Diolaiti N，Melchiorri C，Stramigioli S．Contact impedance estimation for robotic systems．IEEE Transactions on Robotics，2005，21(5):925～9352 Zhang Y，Sharf I．Experimental validation of nonlinear compliant contact force models．Proceedings of the ASME 2007 International Design Engineering Technical Conferences ＆ Computers and Information in Engineering Conference，September 4–7，2007，Las Vegas，Nevada，USA3 洪嘉振．计算多体系统动力学．北京:高等教育出版社，1999(Hong Jiazhen．Computational dynamics of multibody systems．Beijing:Higher Education Press，1999(in Chinese))4 王礼立．应力波基础．北京:国防工业出版社，2005(Wang L L．Foundation of stress waves．Beijing:National Defence Industry Press，2005(in Chinese))5 Thomas Klisch．Contact mechanics in multibody system．Multibody System Dynamics，1998，2(4):335～3546 黄筑平．连续介质力学基础．北京:高等教育出版社，2003(Huang Z P．Fundamentals of continuum mechanics．Beijing:Higher Education Press，2003(in Chinese))7 董富祥．刚柔耦合多体系统碰撞动力学建模理论与实验研究［博士学位论文］．上海:上海交通大学，2010(Dong Fuxiang．Study of modeling theory and experiment for impact dynamicsof rigid－flexible coupling multibody systems．［PhD Thesis］．Shanghai:Shanghai Jiao Tong U－niversity，2010(in Chinese))*The project supported by the National Natural Science Foundation of China(11132007，11202126)，Foundation of Doctoral Program of Ministry of Education ofChina(20090073110009)and Natural Science Foundation of Shanghai(11ZR1417000) † Corresponding author E－mail:zhuyongliu@sjtu．edu．cn。