（机械制造及其自动化专业论文）基于adams平台的柔性体仿真理论的若干研究.pdf

基于a d a m s 平台的柔性体仿真理论的若干研究 摘要 本文以多体系统动力学理论为基础，以a d a m s 2 0 0 5 为平台，研究了柔性 体仿真的相关问题。本文的研究重点主要集中于柔性体的建模、求解等相关理 论。最后，以简单系统为例，在a d a m s 平台上进行了仿真。通过和刚性体进 行对比分析，说明了柔性体对系统性能的影响。 关键词：动力学，多体系统，柔性体，刚性体，建模，仿真 r e s e a r c ho nt h et h e o r yo ff l e x i b l eb o d ys i m u l a t i o nb a s e do i l t h ep l a t f o r mo fa d a m s a b s t r a c t b a s e do nt h ed y n a m i c so fm u l t i b o d ys y s t e m ，s o m ei s s u e sa b o u tf l e x i b l e b o d yk i n e m a t i c s & d y n a m i c ss i m u l a t i o na r er e s e a r c h e d o nt h ep l a t f o r mo f a d a m s2 0 0 5 t h i sp a p e ri sm a i n l ya b o u tt h et h e o r e t i e a ls t u d yo nm o d e l i n ga n d s o l v i n gf o rf l e x i b l eb o d y b yt a k i n ga ne x a m p l eo fas i m p l es y s t e m ，t h es i m u l a t i o n p r o c e s si si l l u s t r a t e do nt h ep l a t f o r mo fa d a m s b y c o n t r a s tw i t hr i g i db o d y ，t h e i n f l u e n c eo nt h ep e r f o r m a n c eo ft h es y s t e mc a u s e db yf l e x i b l eb o d yi sr e v e a l e d k e y w o r d s ：d y n a m i c s ，m u l t i b o d ys y s t e m ，f l e x i b l eb o d y ，r i g i db o d y ，m o d e l i n g ， s i m u l a t i o n 插图清单 图1 1 建模与求解一般过程3 图2 1 柔性体上节点的位置5 图2 2 柔性体变形模型7 图2 3 迭代流程图1 2 图2 4 使用a d a m s f l e x 的步骤1 9 图3 1 威尔逊9 法2 4 图4 1 平面曲柄滑块系统3 5 图4 2 未变形时的平面曲柄滑块系统3 6 图4 3 简支梁的振动3 6 图4 4 刚性连杆o a 4 6 图4 5 状态一，g s t i f f & i n d e x 3 4 9 图4 6 状态二，g s t i f f & i n d e x 3 5 0 图4 7 状态三，g s t i f f & i n d e x 3 5 0 图4 8 状态四，g s t i f f & i n d e x 3 5 1 图4 9 状态二，w s t i f f & i n d e x 3 5 1 图4 1 0 状态三，w s t i f f & i n d e x 3 5 1 图4 1 l 状态四，w s t i f f & i n d e x 3 5 2 表格清单 表2 1 惯性时不变矩阵列表9 表3 1 三种积分方式的比较3 3 表4 1 柔性连杆的自由态正则模态4 9 独创性声明 本人声明所呈交的学位论文是本人在导师指导下进行的研究工作及取得的研究成 果。据我所知，除了文中特别加以标注和致谢的地方外，论文中不包含其他人已经发表 或撰写过的研究成果，也不包含为获得 金魍王些太堂 或其他教育机构的学位或 证书而使用过的材料。与我一同工作的同志对本研究所做的任何贡献均已在论文中作了 明确的说明并表示谢意。 学位论文作者签名： 贸杉签字日期：知汐年g 月g 日 学位论文版权使用授权书 本学位论文作者完全了解金目巴王些太堂有关保留、使用学位论文的规定，有权保 留并向国家有关部门或机构送交论文的复印件和磁盘，允许论文被查阅和借阅。本人授 权金壁王些太堂可以将学位论文的全部或部分内容编入有关数据库进行检索，可以采 用影印、缩印或扫描等复制手段保存、汇编学位论文。 ( 保密的学位论文在解密后适用本授权书) 学位论文作者签名： 泼乳 签字日期：乙名年6 月g 日 学位论文作者毕业后去向： 工作单位： 通讯地址： 导师签名： 签字日期：汐孑年易月fd 日 电话： 邮编： 致谢 本文是在尊敬的导师王健强教授的精心指导下完成的。近三年来，王老师 在学习上给予我很多的关怀和帮助。王老师为学生提供了很好的学习机会，并 且精心指导，使学生充分了解和学习了工业机器人技术及其自动化解决方案。 王老师渊博的学识，开阔的思维，严谨的治学态度，对我产生了深刻的影响! 这三年的研究生生活，更是学生一生中宝贵的财富。学生所取得的所有成绩和 进步，都包含着恩师的心血。王老师渊博的学识、开阔的视野、严谨的治学态 度、求实创新的工作作风，永远引导和激励着学生在科学技术的殿堂里探索前 进。在此论文完成之际，谨向王老师致以最崇高的敬意和衷心的感谢。 同时，感谢实验室的全体同学在生活和学业上给我提供的支持和帮助。 感谢我的家人。他们的支持和理解是我完成学业的前提和动力。 最后向在百忙之中评审本文的各位专家老师表示衷心的感谢。 作者：戴洪光 2 0 0 8 年5 月 第一章绪论 1 1 多体系统动力学简介 机械系统动力学分析与仿真是随着计算机技术的发展而不断成熟的，多体 系统动力学是其理论基础。多体系统是指由多个物体通过运动副连接的复杂机 械系统。多体系统动力学的根本目的是应用计算机技术进行复杂机械系统的动 力学分析与仿真，它是在经典力学基础上产生的新学科分支。在经典刚体系统 动力学的基础上，它经历了多刚体系统动力学和计算多体系统动力学两个发展 阶段，目前已趋于成熟。 计算多体系统动力学中所研究的多体系统，根据系统中物体的力学特性可 分为多刚体系统、多柔体系统和刚柔混合多柔体系统。多刚体系统是指可以忽 略系统中物体的弹性变形而将其当作刚体来处理的系统。该类系统常处于低速 运动状态。多柔体系统是指系统在运动过程中会出现物体的大范围运动与物体 的弹性变形的耦合，从而必须把物体当作柔性体处理的系统。大型、轻质而高 速的机械系统常属于此类。如果多柔体系统中有部分物体可以当作刚体来处理， 那么该系统就是刚柔混合多体系统，这是多体系统中最一般的模型。 1 2 多体系统动力学研究状况 多体系统动力学的核心问题是建模和求解问题，其系统研究开始于2 0 世纪 6 0 年代。从2 0 世纪6 0 年代到8 0 年代，侧重于多刚体系统的研究，主要是研 究多刚体系统的自动建模和数值求解；多体系统动力学在2 0 世纪7 0 年代逐渐 引起人们的注意，在一些系统( 如高速车辆、机器人、航天器、高速机构、精 密机械等) 中，柔性体的变形对系统的动力学行为产生很大影响；到了2 0 世纪 8 0 年代中期，多刚体系统动力学的研究已经取得一系列成果，尤其是建模理论 趋于成熟，但更稳定、更有效的数值求解方法仍然是研究的热点；h a u g 等人在 当时确立了“计算机多柔体系统动力学 这门新的学科；2 0 世纪8 0 年代之后， 多体系统动力学的研究更偏重于多柔体系统动力学，这个领域也正式被称为计 算多体系统动力学。它至今仍然是力学研究中最有活力的分支之一，但已经远 远的超过一般力学的涵义。 对于多刚体系统的建模方法，从2 0 世纪6 0 年代到8 0 年代，在航天和机械 两个领域形成了两类不同的数学建模方法，分别称为拉格朗日方法和笛卡尔方 法。2 0 世纪9 0 年代，在笛卡尔方法的基础上，又形成了完全笛卡尔方法。这 几种建模方法的主要区别在于对刚体位形描述的不同。 航天领域形成的拉格朗日方法，是一种相对坐标方法，是以系统每个铰的 一对邻接刚体为单元，以一个刚体为参考物，另一个刚体相对该刚体的位置由 铰的广义坐标( 又称拉格朗日坐标) 来描述。广义坐标通常为邻接刚体之间的 相对转角或位移。这样开环系统的位置完全可以由所有铰的拉格朗日坐标阵目 所确定。其动力学方程的形式为拉格朗日坐标阵的二阶微分方程组，即 a ( q ，f ) 口= b ( q ，香，f ) ( 卜1 ) 机械领域形成的笛卡尔方法是一种绝对坐标方法，以系统中每一个物体为 单元，建立固结在刚体上的坐标系，刚体的位置相对于一个公共参考基进行定 义，其位置坐标( 也可称为广义坐标) 统一为刚体坐标系基点的笛卡尔坐标与 坐标系的方位( 也称姿态) 坐标，方位坐标可以选用欧拉角或欧拉参数。单个 物体位置坐标在二维系统中为3 个，三维系统中为6 个( 如果采用欧拉参数为 7 个) 。对于由个刚体组成的系统，位置坐标阵q 中的坐标个数为3n ( 二维) 或6n ( 或7n ) ( 三维) ，由于铰约束的存在，这些位置坐标不独立。系统动 力学模型的一般形式可表示为 么口+ ；允= b ( 1 2 ) 【o ( q ，f ) = 0 ( 卜3 ) 式中：为位置坐标阵g 的约束方程，。为约束方程的雅可比矩阵，兄为 拉格朗日乘子。 对于多刚体系统的求解问题，通常采用符号一数值相结合的方法或者全数 值的方法。符号一数值方法是先采用基于计算代数的符号计算方法，进行符号 推导，得到多刚体系统拉格朗日模型系数矩阵简化的数学模型，再用数值方法 求解二阶常微分方程( o d e ) 问题。由于计算机技术的发展，目前全数值方法也 较为流行，就是将多刚体系统拉格朗日数学模型当作一般o d e 问题进行求解， 这方面的技术已经较为成熟。 多柔体系统的动力学数学模型，其形式与多刚体系统相同，可以借鉴多刚 体系统数学模型的求解方法。只是混合坐标中描述浮动坐标系的刚体坐标通常 是慢变大幅值的变量，而描述相对于浮动坐标系弹性变形的坐标却为快变微幅 的变量，两类变量出现在严重非线性与时变的耦合动力学方程中，其数值计算 呈病态，将出现多刚体系统中见不到的数值计算困难。 1 3 计算多体系统动力学建模与求解的一般过程 一个机械系统，从初始的几何模型，到动力学模型的建立，经过对模型的 数值求解，最后得到分析结果，其流程如图1 - 1 所示【1 1 。计算多体系统动力学 分析的整个流程，主要包括建模和求解两个阶段。建模分为物理建模和数学建 模。物理建模是指由几何模型建立物理模型。数学建模是指从物理模型生成数 学模型。几何模型可以由动力学分析系统几何造型模块所构造，或者从通用几 何造型软件导入。对几何模型施加运动学约束、驱动约束、力元和外力或外力 矩等物理模型要素，形成表达系统力学特性的物理模型。物理建模过程中，有 时候需要根据运动学约束和初试位置条件对几何模型进行装配。由物理模型， 采用笛卡尔坐标或拉格朗日坐标建模方法，应用自动建模技术，组装系统运动 2 方程中的各系数矩阵，得到系统数学模型。对系统数学模型，根据情况应用求 解器中的运动学、动力学、静平衡或逆向动力学分析算法，迭代求解，得到所 需的分析结果。联系设计目标，对求解结果再进行分析，从而反馈到物理建模 过程，或者几何模型的选择。如此反复，直到得到最优的设计结果。 i 。 建模 ；i 求解 ；j ： i 一一- i l 一一一。 i ；! 图卜1 建模与求解一般过程 在建模和求解过程中，涉及到几种类型的运算和求解。首先是物理建模过 程中的几何模型装配，图1 - 1 中称为“初始条件计算 ，这是根据运动学约束和 初始位置条件进行的，是非线性方程的求解问题；再就是数学建模，是系统运 动方程中的各系数矩阵自动组装过程，涉及大型矩阵的填充和组装问题；最后 是数值求解，包括多种类型的分析计算，如运动学分析、动力学分析、静平衡 分析、逆向动力学分析等。运动学分析是非线性的位置方程和线性的速度、加 速度方程的求解。动力学分析是二阶微分方程或二阶微分方程和代数方程混合 问题的求解。静平衡分析从理论上讲是一个线性方程组的求解问题，但实际上 往往采用能量的方法。逆向动力学分析是一个线性代数方程组求解问题。这里 面，最复杂的是动力学微分代数方程的求解问题，它是多体系统动力学的核心 问题。 1 4 本课题的研究意义和内容 本课题研究的主要内容如下： 多柔体系统运动学和动力学的基本理论； a d a m s 平台下，柔体系统建模和求解理论； a d a m s 平台下对于曲柄滑块机构的仿真； 通过对比刚性机构和柔性机构的运动学和动力学的不同表现，研究机构 柔性化给系统带来的影响； 本文通过将刚性系统柔性化，增加了仿真的真实性，使得系统的运行过程 更接近于真实的表现。此外，通过对比不同情况下，系统的运动学和动力学方 面的表现，加深了对柔性系统本身的特殊性的认识，为含有柔性机构的机械系 统的设计起到了一定的启示作用。 4 第二章多柔体系统的建模 2 1 多柔体系统建模基本理论 2 1 1 多柔体系统运动学方程的建立 1 柔性体系统中的坐标系 柔性体系统中的坐标系如图2 - 1 所示，包括惯性坐标系) 和动坐标系( e 6 ) 。 前者不随时间而变化，后者是建立在柔性体上，用于描述柔性体的运动。动坐 标系可以相对惯性坐标系进行有限的移动和转动。动坐标系在惯性坐标系中的 坐标( 移动、转动) 称为参考坐标。 e 7 图2 1 柔性体上节点的位置 与刚体不同，柔性体是变形体，体内各点的相对位置时时刻刻都在变化， 只靠动坐标系不能准确描述该柔性体在惯性坐标系中的位置。因此，引入弹性 坐标来描述柔性体上各点相对坐标系统的变形。这样柔性体上任一点的运动就 是动坐标系的“刚性”运动与弹性变形的合成运动。由于柔性体上各点之间有 相对运动，所以动坐标系的选择不是采用连体坐标系，而需要采用随着柔性体 形变而变化的坐标系，即“浮动坐标系 。 在研究多柔体系统时，合适的坐标系是非常重要的。在确定浮动坐标系时， 有两个准则：便于方程建立求解；柔性体的刚体运动与变形运动的耦合尽 量小。目前常见的浮动坐标系大致有如下5 种：局部附着框架、中心惯性主轴 框架、蒂斯拉德框架、巴克凯恩斯框架以及刚体模态框架。采用何种需要据实 际情况而定。 2 柔性体上任意点的位置、速度和加速度矢量【2 】【1 2 】【1 5 1 在分析刚体平面运动的时候，把复杂的刚体平面运动分解为几种简单的运 动。n 3 在对柔性体的运动，尤其是在小变形的情况下，也可以采用类似的方法。 如某柔性体从位置厶运动到位置厶，其间运动可以分解为：刚性体移动一刚性 转动一变形运动。对于柔性体上任意一点p ，其位置向量为： ，= r o + a ( s p + 甜p ) ( 2 1 ) ，为点在惯性坐标系中的向量；r o 为浮动坐标系原点在惯性坐标系中的向 量；a 为方向余弦矩阵；s 。为柔性体为变形时尸点在浮动坐标系中的向量；u 。为 相对变形向量。对于尸点，该单元的变形采用模态坐标来描述，有： u p = p q ( 2 - 2 ) 式( 2 - 2 ) 中，。为点尸满足里兹基向量要求的假设变形模态矩阵，g ，为 变形的广义坐标。 柔性体上任一点的速度向量及加速度向量可以通过式( 2 1 ) 求时间一阶导 数和二阶导数得到： 户，= 矗+ a ( s p + “p ) + a r b p 口r ( 2 3 ) 芦，= 芘+ a ( s p + “p ) + 2 彳p 以+ 彳p 弼 ( 9 - 4 ) 2 1 2 多柔体系统动力学方程的建立 1 外加载荷4 1 【6 1 【1 3 l 【1 4 l 在a d a m s 软件中，外加载荷包括单点力与扭矩、分布式载荷以及残余载 荷三部分。 ( 1 ) 单点力载荷 施加于柔性体上某一标记点的单点力和扭矩必须投影到系统的广义坐标上 才能起作用。力和扭矩以矩阵形式写出，在标记点k 的局部坐标系下表示为： 最= 眈，工，f z 2 ，t k = i t , ，f ，t z 7 ( 2 5 ) 广义力q 由广义平动力、广义扭矩( 以欧拉角表示的广义力) 和广义模态 力组成，可表示为： q = q ，繇，a m r ( 2 - 6 ) 平动坐标下的广义力可以通过转换单点力最到全局坐标基e 7 下来获得， 即： q ，= 彳最 ( 2 7 ) 其中，彳为标记点k 上的坐标系相对于全局坐标的欧拉角变换矩阵，可以 表示为如= 冬锄厶如图2 2 所示，如为局部坐标系b 点相对于全局坐标系 的转换矩阵，即方向余弦阵：4 。为因节点尸的小变形引起的标记坐标的方位变 化而引入的转换矩阵；4 破为定义在柔性体上标记点处的坐标系相对于尸点坐标 系的常值变换矩阵。 6 g 图2 2 柔性体变形模型 作用在柔性体标记点处的合力矩，可用相对于全局坐标的矢量矩阵表达为： t r o t = 如砭+ p x 如砟 ( 2 8 ) p 为从局部坐标基e 6 到受力点处的位置矢量，用斜方阵表示，上式可改写 为： t r o t = 如& + 弛最 ( 2 9 ) 将物理坐标系下的力矩向基于欧拉角的坐标系转换，并利用式= 印，可 得广义合力矩为： g = 么b 】r 【如砭+ 五缸瓦】 ( 2 1 0 ) 通过投影单点力合单点力矩到模态坐标上，可得到p 点处的广义模态力。 如在标记点k 处施加力最和力矩砭，通过转换到全局坐标上，可得广义模 态力为： q f = e + o t r ( 2 1 1 ) 其中，和；对应于节点p 处的平动和转动自由度的模态斜方阵。由于 模态矩阵只定义在节点，故单点力和单点力矩只能施加于节点处。 ( 2 ) 分布式载荷 在a d a m s 软件上，分布式载荷可以通过m f o r c e 的方式来创造。通常在 f e m ( f i n i te l e m e n tm e t h o d s ) 软件中有运动方程： 磊缉+ k x = f ( 2 - 1 2 ) 其中，m 和足为柔性体上有限单元的质量和刚度矩阵，而x 和f 为物理节 点的自由度矢量和载荷矢量。 利用模态矩阵将式( 2 - 1 2 ) 转换到模态坐标g 下，有： f m o q + ，k o q = 7 f ( 2 - 1 3 ) 上式可简化为： 7 聊+ 幻= f ( 2 1 4 ) 其中，m 和k 分别为广义质量矩阵和广义刚度矩阵，厂为模态载荷矢量。 节点力矢量在模态坐标上的投影为： 厂= o r f ( 2 - 1 5 ) 上式中，如果f 是时间的函数，则求解时计算费用太大。一种替代的方法 是假设空间依赖性和时间依赖性可以分开，把载荷f 看成是一系列依赖于时间 的静态载荷的线性组合，即： f ( t ) = s , ( t ) f 1 + + o ) e ( 2 1 6 ) 载荷在模态坐标上的投影计算可在有限元文件m n f 过程中完成，而不必在 a d a m s 仿真中重复进行。如果定义一系列静载荷为载荷矢量，并使其与系统 响应显性相关，即表示成f ( q ，t ) 的形式，则模态力又可表示为： f ( q ，r ) = s l ( q ，t ) f + ( g ，f ) 五 ( 2 1 7 ) ( 3 ) 残余载荷 先前的假设是施加的载荷向模态坐标的投影已实现，即式( 2 一1 5 ) 成立。 然而由于模态截断，很多情况下施加的力并不能进行投影。这些力被称为残余 载荷。可以将其看为已投影到了被截断的高阶模态坐标上。残余载荷可表示为： a f = f p r 】。1 f ( 2 1 8 ) 与残余载荷相关的是残余矢量，可看成是把残余载荷施加于柔性体时产生 的变形。残余矢量可被当作模态振型加入到c r a i g b a m p t o n 模态基中，增强的 模态基完全能够捕捉残余载荷。否则，残余力的丢失不可避免。 2 多柔体系统的能量 ( 1 ) 动能和质量矩阵心町 考虑节点p 变形前后的位置、方向和模态，柔性体的广义坐标可以表示为： 善= 【x ，y ，z ，y ，0 ，矽，g ，( f = 1 ，m ) 】。= ，j c ，g 】1 ( 2 1 9 ) 速度表达式( 2 - 3 ) 在系统广义坐标式( 2 - 1 9 ) 的时间导数手中表示为： 咋= i , - 彳( 跏干“尸) b ，a o 尸】f ( 2 2 0 ) 柔性体的动能为： t = 强1 v d y l em e v r v e + 妒l p 带 其中，m 尸和分别在节点p 的节点质量和节点惯性张量； 点b 相对于全局坐标基的角速度在局部坐标基中的斜方阵表示。 和关系式绵= 昂矽代入式( 2 2 1 ) ，得到动能的广义表达式： ( 2 - 2 1 ) g p b = b 婶，为 将式( 2 - 2 1 ) 丁= m ( 善) 善 ( 2 啦) 上式中的质量矩阵m ( f ) 为3 x 3 维的方阵，表示为： 8 蜒，= 鞋 ( 2 - 2 3 ) 其中下标r ，m 分别表示平动、旋转和模态自由度。质量矩阵的6 个独立 分量分别表示为： m t t = v l m 睁= 一a ( t 2 - t - 3 ) b 峨：rp，allt(m巧=+a才z#3丁)一弓g，b(2-24)b q峨= 7 p 7 一( 巧+ 才丁) 一弓g ， b 帆= 矿( 丁4 + 巧劬) m m m = f 6 其中9 个惯性时不变矩阵见表2 - 1 。 表2 - 1 惯性时不变矩阵列表 惯性时不变矩阵 模态坐标数维数 n f 1 = 脚 标量 p = l r 2 = m 尸郎 ( 3x1 ) p = i ( 3 m ) 才= 坼p j = 1 ，m 尸= l j v ( 3 m ) f 4 = 朋p 昂| p + o p p i l | v ( 3xm ) 巧= 聊p 面彤p j = i ，m p = i ， ( m m ) f 6 = 朋尸p t p + ；。p p z l j v ( 3 3 ) = m p ? i p + i p j 口= l j = 1 ，m ( 3 x3 ) 巧8 = 历p 昂面坷 | d = l j ，k = i ，m ( 3x3 ) = 面肚 p 皇l 从式( 2 2 4 ) 中可以明显看出质量矩阵与模态坐标显性相关，而且由于引 入转换矩阵a 和b ，质量矩阵与系统的方向坐标也显性相关。质量矩阵中的9 9 个惯性时不变矩阵f 1 f 9 可通过计算有限元模型的个节点信息在预处理过程 中一次性得到，从而简化运动微分方程的求解。节点信息包括：每个节点的质 量m ，未变形时的位置矢量s p 以及模态矩阵p 。 ( 2 ) 势能和刚度矩阵 势能一般分为重力势能和弹性势能两部分，可用下列二次项表示： 1 w = 睨( 孝) + 孝2 膨 z ( 2 2 5 ) 在弹性势能中，k 是对应于模态坐标q 的结构部件的广义刚度矩阵，通常 为常量。重力势能睨表示为： 嵫2l 所g a z e2 l p 阮+ a ( s e + 啡g ) 】7 鲥矽 ( 2 - 2 6 ) 出： 其中，g 表示重力加速度， 厂：堡： 。带。 重力对求导得： ( l p d w ) g 嚣【l p ( s e + o e q 煳 彳( f p ；d 形) g ( 2 - 2 7 ) ( 3 ) 能量损失和阻尼矩阵 阻尼力依赖于广义模态速度并可以从能量损耗函数i 的二次项中推导得 r ：扫0 0 二 ( 2 - 2 8 ) 式( 2 - 2 8 ) 称为r a y l e i g h 能量损耗函数。矩阵d 包含阻尼系数以，它是常 值对称阵。当引入正交模态振型时，阻尼矩阵可用对角线为模态阻尼率g 的对 角阵来表示。对于每一个正交模态的阻尼率都可以去不同值，而且还能以该模 态的临界阻尼，的比值形式给出。 3 多柔体动力学方程 柔性体的运动方程从下列拉格朗日方程导出： 瓦d ( 嘴o l 一琵o l + 瓦o f + ( 参见一q = o ( 2 _ 2 9 ) j 瓦嘴- 琵+ 瓦+ 毒) _ 一_ u ( 2 一 【 莎：o 。 其中，杪为约束方程；允为对应于约束方程的拉氏乘子；f 为如式( 2 1 9 ) 定义的广义坐标；q 为投影到孝上的广义力；l 为拉格朗日项，定义为三= t - w ， 丁和w 分别表示动能和势能；f 表示能量损耗函数。 将求得的t ，w ，r 代入式( 2 - 2 9 ) ，得到最终的运动微分方程为： m 芋+ 竹罅i 瓦o m 毋7 善+ 聪+ 疋+ 彤+ 够五= q 2 - 3 0 1 0 其中，善，善为柔性体的广义坐标及其时间导数；m ，肪为柔性体的质量矩阵 及其对时间的导数；娑为质量矩阵对柔性体广义坐标的偏导数，它是一个 d c ( m + 6 ) ( m + 6 ) ( m + 6 ) 维张量，m 为模态数。 2 1 3 多柔体系统的逆动力学 多柔体系统逆动力学问题与多刚体系统逆动力学的相同之处在于：给定末 端的情况下，求关节变量和驱动力( 矩) 。 v l 但是也有不同之处，即多柔体系统 无法通过纯运动几何学关系建立末端物体的位置、姿态变量与各关节变量问的 一一对应关系。它们之间的关系必须通过动力学控制方程的求解才能建立起来， 即需要通过求解方程： q 名一级= a , + q k 旮一 缉(2一)1 3 1 其中，口一约束雅可比矩阵； 五一拉格朗日乘子： q ，一广义力向量： q 一速度二次项广义力向量； k 一刚度矩阵： m 一质量矩阵； 口一广义坐标向量。 但是由于方程左端是待求的，而在方程右端，除末端物体的位置、姿态以 及速度和角速度给定外，其余的广义坐标都是未知的。又m 、k 及q 也都是广 义坐标的时变函数，因此很难求解这组方程。本节主要介绍一种以多刚体系统 逆动力学方法为基础的近似迭代算法。 这种迭代算法是充分利用柔性体的变形在通常情况下是线弹性小变形的特 点，把柔性体的变形运动看作在未变形状态附近的微振动。于是，如下的逆运 动学分析的迭代过程，可以认为能够得到足够精确的解。 i ) 在把柔性体看作刚体的情况下，由给定的末端物体的位置、姿态以及速 度和角速度，根据多刚体逆运动学分析求得各关节变量的值，再由多刚体动力 学方程求得所需的关节驱动力( 矩) ； i i ) 用多刚体动力学分析所得的关节驱动力( 矩) 作为第一次近似，在给 定系统各物体的刚性运动的条件下，求得其变形运动； i i i ) 用已知的变形运动和关节变量的值，求得二次近似的关节驱动力( 矩) 。 反复迭代直至达到所要求的驱动力( 矩) 的精度为止。 迭代流程框图如图2 。3 所示。 图2 3 迭代流程图 1 2 按这种方法进行分析的时候，首先要忽略弹性变形对整体刚性运动的影响， 然后将这些量代入弹性坐标的动力学方程中，求得弹性坐标值。最后应将已知 的刚性及弹性坐标代入柔性模型的动力学方程中的那些相应于刚性坐标的方程 中，从而求出驱动力矩的修正值。注意在随后的迭代分析中，刚性坐标总是被 取为最初从刚性模型逆运动学所求得的值，因此这种方法在很大程度上忽略了 弹性变形对整体刚性运动的影响。在大多数机械系统中，这种影响的确是很小 的，应用非常简化的动力学模型，就能得到令人满意的结果。 2 1 4 多体系统动力学中的刚性( s t i f f ) 问题 刚性问题存在于多刚体系统动力学的某些情形，更普遍的存在于多柔体系 统动力学中，是多体系统动力学的一个重要问题。【2 4 】在多刚体系统运动过程中， 可能会由于系统构件之间过大的差异，如不同物体特性参数的差异、力元( 如 弹簧、阻尼器) 参数的差异等，致使系统中构件运动速度差别很大，从而使描 述系统运动的微分方程呈现刚性特性。在多柔体( 或刚一柔混合) 系统动力学 中，由于柔性体空间大范围运动和其本身小副值弹性变形的耦合，更容易出现 刚性问题。 为描述刚性方程性质，先考虑线性系统： 譬：砂i j f ) + 舯训 ( 2 。3 2 其中，y o ) = y lo ) ，y o ) 】r r 珊为解向量函数，0 ) - - 匦o ) ，九o 妒r m 为已 知向量函数，么灭”为常系数矩阵。设其特征值为2 ，= 口，+ 帮，j = l ，掇， 相应的特征向量为孝，。 定义1 ：如果在线性系统( 2 3 2 ) 中，a 的特征值z ，( ，= l ，o l * 9 朋) 满足： r e 协，j 1 则称式( 2 - 3 2 ) 为刚性方程式，比值s 称为刚性比。通常刚性比s 达到 0 0 0 尸p 1 ) 就认为是刚性的。 对于非线性系统： 鲁= 令歹o ) 为式( 2 3 3 ) 满足初始条件y ( 口) = 儿的精确解， 方程式( 2 - 3 3 ) 作线性逼近： 譬= ，o 炒，歹例i i + o ，夕 t 、，l ，、，v、7 ，、， “ 或 詈= ，协+ f o ) ( 2 - 3 3 ) 在解歹o ) 的邻域内对 其中，p ) 是在点o ，歹o ) ) 处f ( t ，) ，) 的雅可比矩阵8 f ( t ，y ) i ) y 的值。如果矩阵 似 治 卜 j ( t ) 的特征值a i = 彳形) ，= 1 ，m 满足： r e ( a jj 1 j l ，矩阵么或，0 ) 是病态的，故刚性方程也称 为病态方程或坏条件方程。刚性方程数值积分过程中存在快变分量m a x l r e ( a ，) i 和慢变分量m i nr e ( a ) i 的差别。快变分量要求积分步长很小，而慢变分量则使 得在该步长条件下计算步数很多，舍入误差大。这就是刚性方程数值解的困难 所在。考虑到实际问题中可能会出现单个方程情形，或者矩阵彳的特征值有实 部为零或实部为很小正数的情形，可以给出与实际问题更为接近的刚性方程的 定义。 定义2 ：若线性系统( 2 3 2 ) 满足条件： 矩阵彳的所有特征值实部小于不大的正数； 彳至少有一个特征值的实部是很大的负数： 对应于具有最大负实部的特征值的解分量变化是缓慢的。 则称( 2 3 2 ) 是刚性方程。 对于刚性方程数值稳定性的讨论，一般是针对试验方程： 进行的。 鲁吻，y ( o ) = 儿 ( 2 _ 3 6 ) 定义3 ：一个数值方法以定步长h 解试验方程式( 2 3 6 ) ，得到线性差分方 程的解。当n o d 时，若专0 ，则称该方法对步长h 是绝对稳定的。 定义4 ：一个数值方法称为么稳定的。如果它的绝对稳定域包含整个左半 平面r e ( h a ) 0 。 2 2a d a m s f l e x 柔性体建模 a d a m s f l e x 是a d a m s 软件包中的一个集成可选模块，它提供a d a m s 与有限元分析软件a n s y s 、n a s t r a n 、a b a q u s 、i - d e a s 之间的双向数据 交换接口。利用此模块可以考虑物体的弹性，在模型中引入柔性体，从而提高 系统仿真的精度。 2 2 1 柔性体的表示 a d a m s f l e x 是用模态柔性来描述物体弹性的【5 】【引。它基于物体的弹性变 形是相对于连体动坐标系( 或称物体坐标系) 的弹性小变形。同时，物体坐标 系又经历大的非线性整体移动和转动这个假设建立的。 a d a m s f l e x 中的柔性体是用离散化的若干个单元的有限个结点自由度来 表示物体的无限多个自由度的。这些单元结点的弹性变形可近似的用少量模态 1 4 的线性组合来表示。如果物体坐标系的位置用它在惯性参考系中的笛卡尔坐标 x = ( x ，y ，z ) 和反映刚体方位角的欧拉角y = 似，0 ，矽) 来表示，模态坐标用 g = 9 1 ，g ：，鳓) 7 ( m 为模态坐标数) 来表示，则柔性体的广义坐标可选为： h 乒 x y z 缈 目 矽 乃，j = 1 ，m ( 2 - 3 7 ) 那么，柔性体上任一结点( 如第f 点) 的位置向量可表示为： i = 石+ 4 ( 量+ 办g ) ( 2 3 8 ) 其中，4 一物体坐标系到惯性参考系的转换矩阵： 最一结点f 在物体坐标系中未变形时的位置； 谚一对应于结点f 的移动自由度的模态矩阵子块。 将式( 2 3 8 ) 对时间求导，得到该结点的移动速度为： m = 去= 去+ 警( 薯+ 细) + 彳亟盟 = 戈一a ( 霹+ f 6 t q ) c o + 彳谚口 = i e a 墨+ 旃g ) b + 么谚l 孝 ( 2 3 9 ) 其中，c o 一物体坐标系的角速度向量； b e u l e r 角的时间导数与角速度向量之间的转换矩阵； “”一向量对应的对称矩阵。 结点f 的角速度也可以用物体的刚体角速度与变形角速度之和来表示： q = 国+ 矽口 ( 2 4 0 ) 其中，一对应于结点珀勺转动自由度的模态矩阵子块。 2 2 2a d a m s 运动学方程 利用a d a m s 建立机械系统仿真模型时，系统中构件与地面或构件之间存 在运动副的联接。这些运动副可以用系统广义坐标表示为代数方程。【9 】【1 7 】这里 仅考虑完整约束。设表示运动副的约束方程数为砌，则用系统广义坐标矢量表 示的运动学约束方程组为： k ( g ) = f ( g ) ，参( g ) ，轰( g ) r = 0 ( 2 4 1 ) 考虑运动学分析，为使系统具有确定运动，要使系统实际自由度为零，为 系统施加等于自由度( 加一n h ) 的驱动约束： d ( g ，f ) = 0 ( 2 - 4 2 ) 在一般情况下，驱动约束是系统广义坐标和时间的函数。驱动约束在其集 合内部及其与运动学约束合集中必须是独立和相容的。在这种条件下，驱动系 统运动学是确定的，将作确定运动。 由式( 2 4 1 ) 表示的系统运动学约束和式( 2 4 2 ) 表示的驱动约束可以统 一表示为： ( g ，) 2 【- o 。r ( ( g q ，, f t ) ) j 1 = 。( 2 - 4 3 ) 式( 2 4 3 ) 为，2 c 个广义坐标的n c 个非线性方程组，其构成了系统位置方程。 对式( 2 4 3 ) 求导，得到速度约束方程： o ( q ，圣，f ) = 中。( 窖，f ) 圣+ ，( g ，f ) = 0 ( 2 - 4 4 ) 若令v = 一，( q ，r ) ，则速度方程为： o ( q ，口，f ) = 。( g ，f ) 口一v = 0 ( 2 4 5 ) 对式( 2 4 4 ) 求导，可得加速度方程： o ( q ，口，口，f ) = g ( g ，f ) 牙+ ( g ( g ，f ) 口) g 口+ 2 0 讲( g ，f ) 口+ 盯( g ，t ) = 0 ( 2 - 4 6 ) 若令r = - ( o 口口) g q - 2 a , 班口一中盯，则加速度方程为： ( g ，口，牙，r ) = 。( g ，f ) 百- r ( q ，口，f ) = 0 ( 2 4 7 ) 矩阵。为雅可比矩阵。如果的维数为m ，q 的维数为刀，那么。维数为 m x n 矩阵，其定于为( q ) ( “) = 锄f a q j 。在这里西g 为r l c x n c ( 珂厅个运动学约束， ( n c n h ) 个驱动约束，n c 个广义坐标) 的方阵。 坐标【1 0 1 ”1 ，即q = i x ，y ，z ，0 ，伊r ，令r = i x ，少，z 】r ，7 = 眇，0 ，矽r ，q = r r ，7 r r 。构 is i n o s i n 4 0 c o s ? f b2 l s m c 口o s c 。0 s 矽。1 一s o m 口j 2 _ 4 9 li 1 6 吃= 尹 ( 2 5 1 ) 考虑约束方程，a d a m s 利用带拉格朗日乘子的拉格朗日第一类方程的能 量形式得到如下方程： 昙c 若卿酗n 羞 5 2 ， t 为系统广义坐标表达的动能，g ，为广义坐标，q ，为广义坐标g ，方向的广 义力，最后一项涉及约束方程和拉格朗日乘子表达了在广义坐标q j 方向的约束 反力。 a d a m s 中进一步引入广义动量： 弓2 吾 简化表达约束反力为： q = 喜五景 这样方程( 2 4 8 ) 可以简化为： 动能可以进一步表达为： e 一o r 2 g q 卜1 2k r g k + 妒召7 历户 其中m 为构件的质量阵，为构件在质心坐标系下的惯量阵。 将式( 2 5 5 ) 分别表达为移动方向和转动方向有： 幺一罢= 绋一g ( 2 5 7 ) e 一嚣2 q ，一c ， ( 2 - 5 8 ) 其中，幺= 丢( 薏) - 裂d 脚o ) = 肺，罢= 。= e 薏) - 占丁聊。由于b 中 包含欧拉角，为了简化推导，a d a m s 中并没有迸一步推导立，而是将其作一 个变量求解。 式( 2 5 7 ) 可以简化为： m v = 级一c r ( 2 5 9 ) 这样，a d a m s 中每个构件具有如下1 5 个变量( 而非1 2 个) 和1 5 个方程 ( 而非1 2 个) 。 捌 一 一 一 ( ( ( ( 变量： 方程： v = 【圪，巧，圪r r = k 夕，z r 0 = 【0 ，尼，乞r 魄= q ，r 7 = 【，p ，矽】r ( 2 6 0 ) 椰= 级一g v ：厌 寥一薏2q ，一q(2-61) 0 = b 7 j b c o , 眈= 户 集成约束方程a d a m s 可自动建立系统的动力学方程，即微分一代数方程： 户一娶+ 叱日r f ：o 。 p ：望 ( 2 - 6 2 ) u = 口 c b ( q ，f ) = 0 f = f ( u ，q ，f ) 其中，尸为系统的广义动量；日为外力的坐标转换矩阵。 2 2 4 在a d a m s v i e w 中使用柔性体 在a d a m s v i e w 中使用a d a m s f l e x 建立柔性体并进行仿真的步骤如图 2 - 4 所示【1 1 1 ： 第一步，引入检查中性文件( m n f ) ：就是运用a d a m s f l e x 将有限元分 析生成的柔性体描述文件一模态中性文件，输入a d a m s v i e w 中，然后运用 a d a m s v i e w 中提供的约束将柔性体与模型中其它零件连接起来，或者将模型 中已有的该零件的刚体模型替换为柔性体模型。将柔性体引入a d a m s v i e w 后，可以通过查看其信息或使用a d a m s l i n e a r 对其进行检查。 第二步，对柔性体进行设置：a d a m s f l e x 提供了指定柔性体属性的选项， 如：修改柔性体的阻尼率，选择用于仿真计算的模态以减少仿真时间等。 第三步，对模型进行仿真：对包括柔性体的模型进行仿真，与对刚体模型 进行仿真的操作过程类似。 第四步，观察仿真结果：完成仿真之后，可以通过动画显示或绘图窗口来 观察仿真的结果，a d a m s v i e w 会自动记录下零件的位置、速度、加速度等结 1 8 果。还可以观察每一阶模态变形的动画显示。 图2 4 使用a d a m s f l e x 的步骤 2 3 本章小结 本章主要介绍了多柔体系统运动学和动力学的建模理论，以及在 a d a m s f l e x 中建模的主要步骤。此外，还简要介绍了多柔体系统逆运动学的 相关理论。经由对研究对象的数学模型的建立，为接下来的分析求解做好了必 要的准备。 1 9 第三章多体系统的求解 3 1 多体系统求解的数值方法 在建模和求解的过程中，涉及到几种类型的运算和求解：初始条件计算、 运动学分析、动力学分析、逆动力学分析和静平衡分析。初始条件计算是非线 性位置方程的求解；数学建模是系数矩阵操作；运动学分析是非线性的位置方 程和线性的速度、加速度方程的求解；动力学分析是二阶微分方程或二阶微分 方程和代数方程混合的问题；逆动力学分析是线性代数方程组求解；静平衡分 析从理论上是线性方程组的求解。总的来说，计算多体系统动力学涉及的基本 运算包括线性方程组的求解、非线性方程组的求解、常微分方程组( o d e ) 求 解和微分代数方程组( d a e ) 求解。由于对应的解法众多，本节将就其中的几 种常见类型的问题，选择一种方法进行阐述。 3 。1 1 求解非线性方程组的牛顿一拉斐尔( n r ) 法 1 单变量情况