（固体力学专业论文）大挠度空间梁的静、动力学建模、分析与计算.pdf

南京航空航天大学硕士学位论文 摘要 本文从三维连续体的应变能出发，推导了空间大挠度粱的静、动力学有限元方程。 在静力学方程中既考虑了几何非线性项，也考虑了弯扭耦合项，可以描述空f 目粱发生 的大挠度变形。文中建立了描述大挠度梁变形的t l 法和u l 法有限元方程，以悬臂 梁为研究对象t 将其计算结果与理论解进行了比较，验证了本文所建立方程的有效性。 文中推导的动力学有限元方程中的一致质量矩阵吖和弯扭耦合矩阵0 都考虑了几何 非线性项和非线性惯性项，因此是与系统变形有关的时变矩阵。计算表明，在非线性 情况非常弱的情况下本文的方法与线性解几乎吻合；在非线性因素较强的情况下， 与其他相类似的文献进行了比较，得出了一致的结论。 关键词：大挠度；空间梁；几何非线性；t l 法；u l 法；非线性动力学 大挠瘦空间融的静、动力学建模、分析l j 计算 a b s t r a c t b a s e do nt h es t r a i n e n e r g ye q u a t i o nf o rt h e t h r e e d i m e n s i o n a lc o n t i n u u m ，t h ef i n i t e e l e m e n tf o r m u l a t i o nf o rl a r g ed e f l e c t i o n o ft h r e e - d i m e n s i o n a lb e a mi s d e v e l o p e d t h e g e o m e t r i cn o n l i n e a r i t ya n dt h ec o u p l i n ga c t i o no fb e n d i n ga n dt o r s i o na r ec o m p r i s e di n t h es t a t i cf o r m u l a t i o n ，w h i c hc a nd e s c r i b et h el a r g ed e f l e c t i o no ft h et h r e e d i m e n s i o n a l b e a m s f o c u s e do nac a n t i l e v e rb e a m ，t h et o t a ll a g r a n g ef o r m u l a t i o n sa n dt h eu p d a t e d l a g r a n g ef o r m u l a t i o n sa r ef o u n d e d ；t h ev a l i d i t yo ft h ep r e s e u tm e t h o dk sc o n f i r m e db y c o m p a r i n gt h e c a l c u l a t i o nr e s u l t sw i t ht h et h e o r e t i cr e s u l t s t h em a s sm a t r i xa n dt h e c o u p l i n gm a t r i xo fb e n d i n g a n dt o r s i o nv a r yw i 氆t h et i m ec o t l c e r nt ot h ed e f l e c t i o n so f t h e s y s t e m w h e nt h en o n l i n e a ri n e r t i a lt e r ma n dt h e g e o m e t r i cn o n l i n e a r i t y t e r ma r e c o m p r i s e di nt h ed y n a m i cf o r m u l a t i o n t 沁c a l c u l a t i o nr e s u l t ss h o wt h a tt h es o l u t i o ni s v e r yc l o s et o t h es o l u t i o no fl i n e a rs y s t e mw h e nt h en o n l i n e ；a rt e r mi sw e a k ，a st h e i n f l u e n c eo ft h en o n l i n e a r i t yi so b v i o u s ，t h ec o n c l u s i o nd r a w nf r o mt h ep r e s e n tm e t h o di s t h es a m et ot h ec o n c l u s i o no f t h eo t h e rc o r r e l a t i v el i t e r a t u r e s 。 k e y w o r d s ：l a r g ed e f l e c t i o n ；t h r e e - d i m e n s i o n a lb e a m ；g e o m e t r i cn o n l i n e a r i t y ；t o t a l l a g r a n g ef o r m u l a t i o n ；u p d a t e dl a g r a n g ef o r m u l a t i o n ；n o n l i n e a rd y n a m i c s 2 南京兢空靛多乏太学硕 学稳论交 第一章绪论 疆蓑琨找科学技术瓣飞速发鼹， # 线蠖力学润题越寒邈被人翻萎襁。茏荬建霞 航空、躺天、船舶以及举攀等领域，许多结构在大载精、高速、蠢溺_ 摹疆离压豹条臀 下工露，冀受力元薛在以上这些祭舞下的运囊形式袭埂窭一系瓢复杂戆嚣线性税 浆，这时采用线性系统模烈来分析、模拟结构的力学行为就会带来不可忽略的瀑麓。 本文以靛空舷天班及工犍蓊太变影涔瞧结搀尧磷究鹜豢。霹为粱彗搦在簸空魅 天，车辆，潞艟，秘器入以及桥粱等结购中被广泛庭髑，所以在舆俸处璎上造撵在 工程中很多绒构郡采用的紧结构作为分祈对象，胬在通过对粱结构的动、静力学研 究搽采箕内在麴复杂毒 线馁力学行为帮娥簿。 l ，l 建模方法 旱蘩蠢美梁豹缝穆葫力学建骥蹩基予经典豹线嚣援理论，主蘩照巢麓了毒交形 假设。认为粱的横截面在变形以精仍然保持与中线垂点，而且假蹴粱的几何特性及 专葶料性裁缝合麓鼹的线谯芙系。巍于这嶷线牲系统建模方法其蠢方程式麓擎且基于 瓷羧元涟鞠篡髂实疆等魏患t 线瞧系统建模方法蕊王携孛蔽广泛靛应瑟。渡然，受 冀实结籀瓣交澎劳不是镄大瓣翦撵下，这耱线毪系统建按方法憨够缀好豹臻述囊寅 累统的变形情况，满足工程需要。但是对飞行器以及船舶这类巨型结构丽言，奠构 梅中嚣粱、巍结掇鬻露受剽很夫载蕊趣撵慝，筵德鬻鬻念发生较大瓣变形，繇德楚 谯莳屈曲阶段其几何非线燃现象融经十分严重。送时，根据线性祭统建模方法得如 的计算绪巢与实际绪栗桐院具有明显的误差。这魁由于程经典的线烛模型中采嗣f 簟 小变黟鬏浚忽貉了结稿在变澎孛纛移与应变懿 # 线强关系，遮翠筑漾麓了线链纯簸 设。 多数臻究粱黪势力掌a 嚣婆线性翘熬跫建立在张簸k a r m a n 夫挠凄方疆麴蒸艘 。蠹f u , 谈方程中豹簇应交绦帮了蒹些露缝秘变形彰蛹黧蓉戆菲笺餐磺，麓弯蓝应蹙 倪采蘑致避线性疆论中豹簿强庭变方程。分羁舔tl 法( 全拉格鞠 = 法) 和u ，l 法 ( 熙新的掇格朗同法) 求解时，随辫载荷的增加两者的缩祟渐渐出现技大的差别2 l ， 查婆霪皇塑篓堕鳌：塑查兰壅蒸：坌篓羔茎簦一 只有当；掣j 1 的时候，两种方法才能得到同一结果，当挠度较大的时候t l 滋是 l 舔 没有意义的。c o o k ，r d 冽在1 9 7 4 年藏论证了在求解y o n k a r m a n 意义下的大挠度 滴邃啐 ，u 毛法浅子鳖拣黉强黧爨薪套潮毽含了膜癌变葶【l 弯盏蘧攀懿菲线经项e 疆薅关予大挠菠粱的动力学建摸大雾是在爨黎弯趣噩| l 率方程墨稳辈线性项携固 时考虑戮非线性惯性顼。然话根籀h a m i l t o n 原理娥者牛顿第二运动定律褥出的运动 微分方程。 e v e n s e n n 零in 3 y f e 魏秘餐在建立粱戆大挠废运动微分方程嚣孪考惫了梁涎霉镤馁 碾，瞧在搬分方瑷中对粱熊弩兹魏事疆秘然终了线性化歉邂。c r e s p o 妇s i i v a _ ( 6 f 7 与g l y n n 8 m 】，f i o 1 基于h a m i l t o n 原理，程保留弯硅扭曲率的非线性的情况下导出了保 留三爨雾线敬磺静犬撬发浆在任鬈邃赛条释下露夸簸一毒蘧一亵转摄凄翁 线程编 微分方程，褥出了非线性惯性项和非线性曲率项对两端距离不可疑的大挠度梁的动 力学响应有魄较明显的辩晌。d 壬h o d g e s 羊嚣e 嚣。d o w e u 在磷究宣秀辊旋冀静 气凌凌震特梭l l 亨势鬟裁麓h a m i l t o n 嚣潦鞠牟羲繁二邋动定律箍簿露了囊翼翦羧韵 拙伸一弯曲一弯曲一扭转耦台的旋冀的动力学微分方程。在以前工作晌基础上，c r e s p o d as i l v a 对粱豹鸳基一弯靛一掇转一级向投传款动力学觏题又数了进步的秘究。能把 凌孛瑟变彭，鸾熬夔率器骥谴壤熬葵线魏嚣素懿彰穗餐雩l 天r 动力学方纛，势辍鬻 了对菲线犍影响堤为显蓑豹蘸三狳菲线性颈【1 3 j t i l 4 。滋些磺究残暴鞭鹱，在线牲簇 统中被忽略了的由非线饿弯曲曲率所引起的非线性项与非线性惯憔项是同阶次的， 嚣显会霹系绫鹣豌痉产蔓鹋患鹃影蕊，嚣豆系笺抟稳凑璃应蔹藏予裙始条 孛。 综合以上文献的研究埘以看出，建模是基于h a m i l t o n 原理或者牛顿第二运动定 律，在考虑到电子大变形掰雩 超稿位移一应变韵鞯线镶菸盛保留撵瞧粱在弯曲变形 黠熬夸熬辫率瓣嚣绞篷磺翡肇疆上接导爨寒豹粱懿丈挠度运动凌力学徽分方程。 系列实验和实践表明，保留到二阶或者量阶非线性项的影响所建立的梁的大挠度送 嬲微分方程是瀵怒工程或麓装求豹。 1 2 难点以及求解方法 凡翅非线性阙题与线性系统鞠题豹缀岑嚣裂戡夜予线蛀系绞忽瞧了弹性毒枣受 南京航空航天大学硕士学位论文 力后在空间位置的改变。几何非线性问题由于其结构的刚度与变形之间的耦合对问 题的求解带来困难。 目前求解大挠度梁的动、静力学的方法大致来说有三类。一类是解析方法，另 一类是数值方法，最后一类是半解析方法。解析法有摄动法、平均法、k b m 渐进法 以及目前使用广泛的多尺度法等。这些解析方法的优点在于能够定性的给出一些局 部动力学行为，但主要只适用于求解低维、弱非线性动力系统的近似解，对比较复 杂的结构、强非线性等情况没有办法给出定量的分析。工程结构在实际中会遇到各 种各样的问题，然而，在许多情况下，解析解只能求解简单结构的问题，精确的理 论解不容易得到，这时就必须依赖数值方法来描述。数值方法有很多种，有非线性 迦辽金法、配点法、非线性有限元法等。现代计算机技术的迅猛发展为数值解法提 供了广阔的发展空间。 本文研究的对象是一空间大挠度梁结构，其非线性程度取决于其发生的非线性 位移的大小。很难具体的将其划分为强非线性系统或是弱非线性系统。而且由于梁 结构的连续性，其实际系统的惯性、弹性和阻尼等都是连续分布的，因此属于连续 系统或分布参数系统。要确定连续系统中无数个质点的运动形态需要无限多个广义 坐标，显然是不能将其作为简单的低维系统来处理的。在这种情况下，要获得梁的 动、静力学行为特性是很困难的，到目前为止还没有见到有关文献求解几何非线性 梁结构的精确理论解的方法。 由于人们在工程实践中遇到了大量的非线性问题，按照传统的线性有限元方法 计算出来的结果已经不满足生活生产的需要，非线性有限元就是在这种情况下，从 线性有限元的基础上作为计算非线性结构问题的一种数值方法提出来的。非线性有 限元的优点在于：原理简单、易于计算机实现以及计算效率高等。随着现代计算机 技术的飞速发展，非线性有限元方法也得到了迅速的应用和推广。 1 2 1 静力学方程的求解方法概述 对于大挠度梁结构来说，位移和应变之间的已经不再表示为线性关系。而且由 于梁结构的大挠度变形，位移和应变之间呈非线性关系，刚度矩阵也不再是常数矩 阵，而是单元位移的矩阵函数。几何非线性有限元是在传统方法的基础上，引入了 大挠度空间梁的静、动力学建模、分析与计算 应变与位移之间的几何非线性关系，然后再将系统动力学方程中的非线性项作近 似变换，非线性项就可以表示为与节点位移有关的几何刚度矩阵，即 k = k ( d 1 ，d 表示单元的节点位移。平衡方程可写为： k ( d ) d = f ( d )( 11 ) 只有当位移无限小、线弹性材料以及加载过程中边界条件性质不变时，上式才 可以被看成是线性的，一般情况下均是非线性情况。求解上面的非线性方程组的常 用方法是迭代求解。主要有全量法和增量法两大类。全量法又可以分为 n e w t o n r a p h s o n 法，拟n e w t o n r a p h s o n 法和修正的n e w t o n r a p h s o n 法等。增量法 可以分为增量加载法，线性加载法等。 全量法就是在总的载荷作用下进行线性迭代求解计算的方法。这种方法是从求 解非线性方程的根提出来的。在数学上，对于没有直接方法求根的非线性方程，往 往采用线性化序列来逼近。为了表示方便，令 庐( d ) = 庐( d ) 一p 式中 ( 12 ) 庐( d ) = k ( d ) d( 1 3 ) 如果d 是方程( 1 1 ) 的根，则有 ( d ) = 芗( d ) 一p = 0( 1 4 ) 上面是从数学角度来阐述的，从力学角度来看，芗( d ) 表示结构的内力，而( d ) 则是表示内力与外力不平衡的力，也就是内力与外力作用的残差力。只有在找到在 外力p 作用下的精确的位移d ( 方程的根) 的时候，不平衡力( 残差) 为零，即 庐( d ) = 0 。 全量法就是要构造一个线性逼近序列 一，d 町 ( 1 5 ) 使得当n 斗时，有： d ”斗d ，( d ”) j0 ( 1 6 ) 南京航空航天人学硕十学位论文 d 表示系统的精确解( 方程的根) 。 主要有三中构造序列( i 5 ) 的方法，分别是前面介绍的n e w t o n - r a p h s o n 法，修 f 的n e w t o n - r a p h s o n 法以及拟n e w t o n r a p h s o n 法。 n e w t o n - r a p h s o n 法是对具有一阶导数的连续函数痧在d ”处作一阶t a y l o r 展 开，如果d ”是已知的，那么就可以得到线性逼近的序列公式： k ：m h d 扣k 一妒11 dc * l 】= d + d f ( 1 7 ) ，世：叫分别为d n 处的不平衡力和切线刚度矩阵，切线刚度矩阵为： k ! 一，：丝i 捌j d 。d ( n 显然，我们可以事先给定一个位移d 旧，通过上式的迭代求解获得满足式( 1 6 ) 的位移序列。由式( 1 7 ) 构成的每次迭代计算都要重新形成切线刚度矩阵，因此对于 多自由度结构，其计算量是比较大的，但是由于其始终是沿着切线方向逼近精确解， 故收敛较快。 为了避免重新形成切线刚度矩阵，而提出了修正的n e w t o n ，r a p h s o n 法，即在 每次迭代中，将足j ”改用初始切线刚度矩阵，并且为了提高收敛速度，采用每两次 迭代之后对后一次结果进行修正。 拟n e w t o n r a p h s o n 法提出于上世纪6 0 年代初期，但在非线性有限元方面的应 用却始于7 0 年代末。特别是8 0 年代中期，该方法得到了许多研究者的注意。这种 方法实际上也是一种修正的n e w t o n r a p h s o n 法，但它吸收了n e w t o n - r a p h s o n 法和 修正的n e w t o n r a p h s o n 法的长处。 1 , 2 2 动力学方程的求解方法概述 解析方法的发展主要采用小参数摄动法，即将微分方程的解展开成幂级数，这 种方法源于p o i n c a r d l 9 世纪末研究天体运动所提出的直接摄动法。其主要适用对象 是下面的这类弱非线性自治系统的初值问题： 大挠艘空间梁的静、动力学建校、分析与计算 j 靠( f ) + 嗣“( ，) 2 印( “( f ) ，n ( 嘞f 1 8 ) l “( o ) = a o ，女( o ) 2 0 这哥中方法认为方程的解可以写为“= u ( 1 ，) ，将其在线性孵嫩近展开成5 的幕级数： o ，g ) = z o ( r ) + 捌l o ) + 6 2 2 0 ) 十 ( 1 9 ) 将式( 1 9 ) 代入蓟原微分方程，根耩8 的丽次幂，可以得到一缀可戳依次求解的常徽 分方程初缎问题。 但这种方法随者时间的增加会使系统融现永年项。l i n d s t c d t 对p o i n c a r 4 的直接 摄动法 乍了重要的液避，提出了p * l 法，搜其能获德自治系统属期强动的致有效 展开，使得解一致收敛，p l 法又称为奇舜摄动法。2 0 世纪2 0 年代起，乌克兰基 辘举派熬数学家霹菲线往紫激分方程靛濠避瓣法遴行了深入豹繇褒，缝寒l 激有代表 性的成果就是以学者k r y l o v ，b o g o l i u b o v ，m i t r o p o l s k y 的名字命名的k b m 平均法 1 7 1 ，k b m 平均法建基于微分方程瓣常数交易法，继而与小参数法结合而成的。长 期以来苏联学者在发展单囱由度系统非线性理论方面做出了极大的贡献。上世纪5 0 年代，美围学者s t u r r o c k 擒出来了寻求不附阶次的近似解采用不同时间尺度的思想 的多尺度法 1 s 。 但是，由于非线性系统的复杂性，解析法的分析研究仅仪限于1 - 2 自由度系统， 对予 线瞧攘动方程懿末缨骜蘸还没骞一耱逶爱予不嚣类鳖夔方糕懿逶鼹浆簿叛 法，仅仅只有极少数的方程可以求得其精确解1 9 1 ，而且难以对高维系统做如分析。 嚣魏，人们追切需簧寻求缎震其寄离精度、高效率麓数毽计算方法，适用予高维菲 线性系统长期发展预测及念局非线性特性分析。 长期以来，人们为了援商非线性方程的计算精度，研究了许许多多的数值算法 以及各种修芷格式，这些方法实隧上都是代数方程豹逐步遴挂差分格式。鼹藏常用 的数值方法有e u l e r 法、r u n g e k u t t a 法、中心差分法、h o u b o l t 法、w i l s o n 。0 法、 n e w m a r k - f l 法等。 e u l e r 法属于最为简单的显式藏分格式低精度单步法，时间步长o ( h ) 蓑分法造 成的离散谖差为o ( h 2 ) 。在e u l e r 法的基础上，人们为了简化计算减少计算辫，将其 修燕为隐式积分格式的梯形法，继褥发展了中心差分法，拱彤法积中心差分法熬裹 散误差都疑o ( h 3 ) a 预测校正滋始建于a d a m s 多步法，发展为高阶预测校 南京航空航天人学硕士学位论文 正的单步r u n g e k u t t a 法，离散误差为o ( h 5 ) ，一直沿用至今。 以上这些方法从本质上来晓都属于e u l e r 型有限差分法的范畴。所谓e u l e r 型有 限差分法就是指微分方程在时间步长的端点满足，并且根据实际的需要，假设合适 的试探函数在时间域内拟台真实的曲线，以此来近似的逼近真实系统在时间域内的 系统方程。在实际应用上，人们常常采用取平均、线性、等加速度等各种假设，将 非线性微分方程化为每一时间步长出内的线性方程，然后按照中心差分法等递推算 法以及各种修正形式来计算非线性方程的数值解【2 。 在非线性常微分方程的计算上，从方法的不断改进到误差分析，发表的论文浩 如烟云。近来郑兆昌等提出了考虑高阶余项的非线性系统连续线性化模型及其数值 计算方法f 2 1 1 ：随后，他们又对该方法进行了改进，在状态空间罩通过t a y l o r 级数展 开推导出非线性振动系统的瞬态线性化方程，并将t a y l o r 展开式的高阶余量表述为 系统状态的表达式而引入到线性化的方程中，从而进一步提高了计算效率【2 2 1 。 1 3 本文的主要工作以及各章的内容简介 1 3 1 本文主要工作 本文利用d h h o d g e s 和e h d o w e l l 研究直升机旋翼的气动响应特性时所推 导的旋翼的纵向拉伸一弯曲一弯曲一扭转耦合的旋翼的动力学微分方程中所推导的位 移应变的非线性关系，导出了空间大挠度梁的静、动力学有限元方程组，给出 了一般梁类结构在受外载作用下的静、动力学行为的有限元描述方法。根据本文给 出的有限元方程，计算细长悬臂梁结构在外载作用下的静力学变形：比较了在静平 衡位置处梁的固有频率与线性系统中计算出来的梁的固有频率的差别：计算了计及 几何非线性因素的梁的动力学响应问题。 1 3 2 各章内容安排 第一章介绍了几何非线性问题在近年来的研究发展情况；讨论了空间大挠度 梁的建模思路和方法；简要的介绍了几何非线性方程的分析方法。 第二章系统的介绍了几何非线性问题的建模方法以及建模所需要的应变 大挠度空间粱的静、动力学建模、分析与计算 位移的几何非线性关系。根据h h o d g e s 和eh d o w e l l 的研究结果，给出了一般 空间粱结构在变形过程中的局部坐标系与总体坐标系之间的转换矩阵关系。 第三章采用l a g r a n g e 描述，推导了大挠度梁结构的有限元方程。采用不同的 坐标系作为参考系分别推导了两种不同的有限元描述方法。考虑了全拉格朗日法中 忽略掉的高阶余量，推导出了一种可以适用更大加载步长的全拉格朗同法。 第四章介绍了非线性动力学的分析研究方法，讨论了前人的研究成果。在考 虑几何非线性项以及非线性惯性项的基础上导出了空间大挠度梁的动力学有限元 控制方程组，在此基础上，将得到的动力学有限元方程变换到状态空间里，采用 r u n g e k u t t a 法计算了方程的数值解。 第五章本文总结以及前景展望 南京航空航天大学硕十学位论文 第二章大挠度梁的应变位移基本方程 在固体力学中，往往将物体视为无限多个质点的组合。这些物质的质点在空间 中所占有的位置构成了相应的物体形状。 物体形状的变化是通过发生变化前后的两种物体位置的变化对比来确定的 2 3 1 ， 2 ”。一般认为物体发生了变形就一定会存在位移，而物体发生了位移却不一定 会产生变形。根据h e l m h o l t z 分解原理可知，有限变形的产生可能是因为大位移、 大转动或大应变，或者是兼而有之。应用大家熟悉的无限小位移的c a u c h y 公式所 定义的应变张量到有限变形的情况时可能会导致荒谬的结果。 为了得到在普遍情况下物体变形的一般描述，必须引入相应的参考系来考察物 体变形的情况。通常，人们采用物质坐标系和空间坐标系来描述物体的变形，与上 面两种坐标系项对应的是l a g r a n g e 描述和e u l e r 描述两种确定物体变形的描述方 法。不管是l a g r a n g e 描述还是e u l e r 描述，两种方法都是建立在下面这个假设的基 础上的，即：物体由质点构成，质点在空间上是连续分布的，并且物质中的各个质 点在物质的变形过程中是无裂缝不重叠的，代表该过程的数学表达式是：l d j 0 在l a g r a n g e 描述中，将质点在每一瞬时的位置x ，看成是时间t , n 物质坐标，的 函数。在e u l e r 描述中质点在每一瞬时的位置z ，看成是时间t , u 物质坐标x ，的函 数。两种描述方法上的不同使得相应给出的应变位移关系也不同。 由于本文所涉及的研究仅仅局限于固体力学范畴，因此，在文中以后的描述中 都采用以物质坐标为变量的l a g r a n g e 描述方法。在物体变形过程中，只要参考构型 q 。是已知的，就可以用参考构型来确定当前构型。从而建立起几何非线性的应变 位移关系。 梁的变形问题从本质上来说是非线性的，线性假设只是实际问题在一定条件下 的一种简化。在小变形情况下，应变和位移之间近似表现为简单的线性关系，然而， 弹性梁在受到比较大的载荷作用或者发生较大的转动的时候，按照线性假设所建立 的应变位移关系不在成立。必须按照有限变形理论来建立大挠度粱的动、静力方程。 人挠度空间粱的静、动力学建模、分析与计算 在这一章节，将先从建立有限变形情况的坐标转换关系开始，推导出大挠度梁 的应变一位移关系。 2 1 结构模型与坐标系 这里采用的结构模型为一矩形机翼( 如图2 1 所示) ，机翼沿展向没有预扭 转；沿展向作大挠度变形( 如图2 2 、图2 3 所示) ，忽略机翼横截面翘曲的影自： 翼剖面上的质量中心、弹性中心和气动力中心在弦向不重合；机翼根部没有预安装 角。当该机翼的展弦比大于等于1 0 的时候，可以将其简单的看成为发生大挠度和 扭转相结合的粱结构。 图2 1 矩形机翼结构模型 梁的弯曲变形如图2 2 所示，此时，梁不产生沿轴向的扭转。x ，y ，z 为梁的初 始坐标系。梁在横向( 即z 方向) 的位移为w ，弦向( 即y 方向) 的位移为v ，轴 向( 即x 方向) 的位移为“。“，v ，w 分别平行于j ，j ，尼。变形前中性轴上的任意一点 的位置为( x , o ，0 ) ，在变形后位置为 + “，v ，w ) ，一 图2 2 机翼变形及坐标系 南京航空航天大学硕1 ：学位论文 如图2 2 所示，当梁发生沿着z 方向和沿着y 方向的横向弯曲变形后，固连在 梁截面上的坐标系的坐标方向矢量由原来的i ，_ 7 ，f 变成了i ，_ 7 ，f ，i ，_ 7 ，f 分别平 行于x ，y ，z 。坐标系统。 梁截面如图2 3 所示，叩、亭是建立在垂直弹性轴的截面上的截面坐标系，其 交点就是该横截面的剪切中心。图示为空间梁在发生绕x 轴方向的转动角以后所 处的位置。 如图所示，x ，y ，z 坐标系初始坐标系，为分析大挠度梁变形所采用的初始 构型，x ，y ，z 坐标系截面坐标系，为当前构型，现在所需要做的就是找到在两 种构型下的物体位萱的关系。 ， 图2 3 机翼截面及坐标系 2 2 两种坐标系之间的转换关系矧 坐标变换的目的就是希望通过坐标变换找到初始坐标系同变形以后的坐标系 之间的某种函数关系。 在建立坐标转换关系之前，为了以后叙述方便，并且符合本文讨论的目的，先 引入以下几个假设： 1 梁在横向( 即沿y 轴和z 轴方向) 发生大挠度变形在z 轴方向上的变 形是小变形。 大挠度空间梁的静、动力学建模、分析与计算 2梁的端面不发生翘曲。 3线弹性假设( 即梁在变形过程中构成梁的物质始终保持应力应变 的线性关系，不发生屈服) 。 初始坐标系与截面坐标系之间的关系可以表示如下： 即豳 ( 2 1 ) 式中，t 为初始坐标系与截面坐标系之间的转换矩阵。i ，_ 7 ，云和i ，乃岳分别如上所 述。 r 矩阵对于建立大挠度梁的应变位移关系是至关重要的一个量，其矩阵中的 各项可以由“，v ，w ，来确定。 先将转换矩阵用欧拉角表示如下 c o s ，c o s y s i n o s i n f l c o s f c o s 8 s i n - 一c o s o s i n f l c o se - + s i n 8 s i n - c o s f l s i n - c o s o c o s - 一s i n q s i n f l s i n o s i n o c o s g - 一s i n 5 s i n f l c o s 0 s i n8 c o s 口s i n 臼 c o s p c o s 0 ( 2 2 ) 上面矩阵的各个项次可以由简单的矢量叉乘方法得到。式中的各个欧拉角如下 图所示： 图2 4 微元体变形以及欧拉角 南京航空航天大学硕十学位论文 上面的e u l e r 角定义了梁在变形前和变形后的关系。由图可以很容易的得出v ，w 和p ，f 之间的关系如下： c o s c = ( 23 ) 将式( 2 3 ) 代入式( 2 2 ) q a ，可以得到用梁上任意一点的位移所表示的坐标转换矩阵如 下： t = 隔 竺耋匾竺耋 陌 ：! ! 璺& v 匾竺塞o r 踊 却 西 s i n 护肛歹 c o s 口、；一c 争2 却 o r c o s p 雁f 万一s i n 嗉警 陌 一型! 翌：翌：堡量 ( 2 4 ) 现在需要确定的就是欧拉角o , n 变形之间的关系了。考虑在梁上任意一点的位置矢 哆商警 i l = i i 疗 f 矗 硌 n n c s s 大挠皮空间粱的静、动力学建模、分析与计算 量i 的长度，发生了d r 的变形，那么截面坐标系绕垂直通过截面的轴x 也相应的发 q ：娑( 2 5 ) 詈2 毒并 眨s ， r 。= r 7 将式( 2 4 ) ，( 2 7 ) 代入式( 2 6 ) ，从中消去c o ，( - 0 k ，得 0 = 函 a2 vo w 券o 争2 一c 詈，2 6 ( 27 ) ( 2 8 ) 当梁发生的变形满足大转动，小应变的情况下，这时候有：d rad x ，将这一关系代 入式( 28 ) 并且略去式( 2 8 ) 中二阶以上的高阶小量，得： 州一r 窘芸出 ( 2 9 ) 将上式代入式( 2 4 ) ，同样的保持到二阶小量。得到了用梁在几个方向上的变形表达 的转换矩阵r 如下： 丁= 1 一生一丝。， l 一 v 22 c 。s - w s i n 乒c 。s ( 方“w ) ( 1 一孚) v s i n - - w r c 0 s 方_ s i n ( 乒“w ，) ( 1 一争 上式中：阳一f 窘豢出 。i 。痧( 1 一娑) 。方( 1 一嬖) ( 2 1 0 ) 至此，已得到了用梁在几个方向上的变形表达的转换矩阵r ，这对后面推导出大挠 度空间粱的应变位移关系是很重要的。 南京航空航犬大学硕十学能论文 2 3 大挠度空闻粱的凡何菲线性应交位移关系澄j 厩彤一幌识：2 防却西】k 。 i d r ，j ( 2 _ 1 1 ) l d g j 五= r 歹云】 茗三“ + ， ； ( 2 、t ：， 无= 乱吲：。( 2 1 3 ) 令w = v = w = = 0 ，并将该条件代入式( 2 1 0 ) ，褥 神州斟 强砷 上式中：r 。= v l 舢，氏= f l 。 对式( 2 ，1 2 ) ，( 2 + 1 4 ) 衷皴分，痤爱式( 2 。7 ) ，莠缓竣吒：z 盯：，= 气，= 。，可敬缮爨k j 豹 表达式，该液达式非常复杂，在这里就不将它的各个项次一一表达，只将与文牵有 关的l ，岛2 ，衷达如下： 夫撬凄空翔粱趣棼、勘力学建模、分辨等诗嚣 气；塞+ ；) 2 ；尝2 扣2 萋2 嘞一霸窘嘞多年。窘 ：一善娑 ( 2 15 ) 钆= 墨罢 上式晤蠢，二阶以上的两阶小量。 上嚣瑟瘦燹方程中包含了出于辍游搜棼零致静中辙线变形熬貘盛变： ，* 罢+ 吉( 争十i 1c 争& v 2 由摇转联弓 起靛藏瑟上发生的膜废变： 钿殍= i 1 国2 噤) 2 幽纯弯曲所引起的弯曲应碴= ： a 2 vp8 2 w 5 一印可 可 斑弩睦、攮转联合作用所学 起的鸯盛拯转藕合盼应交： ( 2 1 6 ) 0 1 7 ) ( 2 1 8 ) = 彩窘o - v 彩窘 鳓窜) 将上嚣褥裂熬盘炎张量表示巍受艘靛工程瘟变彭式热下； f 占“端蜀， 2 奄：器+ 2 0 ) l s ，f 2 e a3 至髭，褥弱了大挠爱粱戆建交馥秽关系。肖鬻蓉裂，上露黪鏖变霞移 荚系表达式是一缎非线性方程。将式( 2 1 5 ) 退化为平面方程，即：w ；秘：0 ， “喾0 ，v 喾0 ，得： = 罢+ 吾c 一玎窘 z 一， 式8 1 7 ) 就是大家所熟籍的平磷形式静g r e e n 陂变关系式。日暴粱发生的变形裰 小，也就楚普通妁小挠度情况，这对侯，斌q ，1 1 ) 可以进一步运化为枣变形情况下懿 线性应交位移关系。 南京航空航天大学硕十学位论文 a2 v 占，- 一叩萨 ( 2 2 2 1 2 4 物理方程 方程( 2 1 5 ) 是定义在初始构型上的g r e e n 应变，所对应的应力关系也应当是定义 在初始构型上的k i r c k h o f f 应力。由于在这里讨论的问题仅仅局限于大转动、小应 变情况，不必区分c a u c h y 应力和k i r c h h o f f 应力，而且在整个变形的过程中h o o k e 定律均适用口1 。 因此，应力应变关系可表述如下： 盯捌= 上占n t 口= g 。 ( 2 2 1 ) d 吖= g 。f e 为弹性模量，g 为剪切模量。 2 5 小结 l 本章在l a g r a n g e 描述下建立了在大转动、小应变的情况下空间梁的应变位 移方程。 2 从上面得到的应变位移方程来看，即使是小应变情况下，最简单的梁的应 变方程也十分复杂。方程中既包括简单的纯弯曲所引起的弯曲应变项，又包括 轴力引起的膜应变项，还包括复杂的弯曲扭转耦合所引起的弯扭耦合应变 项。 3 讨论了在l a g r a n g e 描述下所适用的物理方程。 4 由于应变方程非常复杂，这给后续的研究带来了很大的困难。因此，如何去寻 找合理、有效的研究方法以及计算手段是开展后续研究所要考虑的首要问题。 大挠度空间梁的静、动力学建模、分析与计算 第三章t l 法和u l 法 前一章基于文 1 2 】的建模方法，即在l a g r a n g e 描述下系统的建立了大挠度空洲 梁的应变位移，并得到了适用于该方程组的物理方程。然而要求出其在外力作 用下的变形还需要其他的条件，如能量原理、虚功原理等。 本文研究的对象为大挠度空间梁，在具体问题上，将选择悬臂梁作为研究的例 子。悬臂梁结构在工程中随处可见，在第二章中所采用的机翼模型就可以看成是一 端固连于机身一端自由的悬臂梁。 几何非线性问题向来都是国内外研究的一个热点，对于梁、杆等结构的几何非 线性分析人们已经做了大量的研究。 为了研究梁柱的弹性稳定性问题，人们早在6 0 年代初就推导了几何非线性粱 单元刚度矩阵。z a r g h a m e e 和s h a h 2 6 1 以及c a i l n o r 【2 7 】等人提出了t ，l 方法，但假定结 构只有小的转动，且只考虑轴向和弯曲的效应，但他们在方程中忽略了剪切和扭转 的影响。t e z c a n 和m a h a p a t r a t 2 剐c h u 以及r a n p a t s t r e i t e r 2 9 】等人从梁柱理论中推导出 了考虑大位移，小转动的u l 法的几何刚度矩阵，他们仍然没有考虑扭转的影响。 他们得出的u l 列式的转换矩阵是根据计算后的位型不断更新的，这也是u ，l 法与 t l 法在计算上的最大的区别。但是他们并没有给出有效的大转动的转换规律o0 1 。 一般来说，由于方程的复杂性，要求解大挠度梁的静力学问题是非常困难的。 首先，在几何非线性情况下，线性系统中常用的叠加原理不再适用；其次，几何非 线性问题的平衡方程是建立在变形后的位置上的，然而变形后的位置在求解前是未 知的，这就是几何非线性问题的复杂性所在：再次，在求解几何非线性问题的时候 往往要进行多步迭代求解，这将带来巨大的工作量，这种工作量远远大于求解线性 问题时的工作量。当然，梁的几何非线性问题并非没有精确解，2 0 0 2 年，为了解决 宁波甬江口的甬江铁路大桥的对接问题，钱伟长导出了非线性大挠度梁的基本微分 方程，并给出了解析解d “。然而，除了平面梁结构能够得到理论解以外，还没有见 到过有关文献给出空间梁的理论解，因此人们求解梁的大挠度变形问题通常是采用 数值解法。 人们求解几何非线性问题的数值方法也可以大致的分为两大类。第一类是通过 南京航空航天大学硕+ 学位论文 一一。 微分方程，寻求满足微分方程的数值解；另一类则是通过建立系统的非线性有限元 方程来求解。 几何非线性问题的微分方程数值解归根到底就是求解满足一定边界条件的非 线性偏微分方程。加权残数法是求解这类非线性微分方程的一种常用方法。具体来 说又可以分为非线性迦辽金法、最小二乘配点法等。其基本思路是构造一组满足全 部或部分边界条件的函数作为试探函数，将其代入到微分方程，方程的残差就反映 了不平衡力的大小。为了尽量减少残差，非线性迦辽金法的方法是使不平衡力关于 各个试探函数的平均做功为零，从而得到一组关于待定系数的非线性方程组，解该 方程组即可得到各个试探函数的系数，各个试探函数的和就是该非线性微分方程的 解；最b - 乘配点法减少残差的办法是记残差的各个项次的平方和为，由最小二 乘法准则得出各个试探函数的待定系数，从而得到非线性微分方程的解”2 1 。 非线性有限元解是在线性有限元基础上发展出来作为求解非线性问题的一种 数值解法。对于大型结构来说求解量大，并且要求求解效率高，而有限元法的最大 的特点在于表达形式整齐统一、组装规则简单方便，故有限元法的可编程性好，便 于在计算机上实现。在现代计算机飞速发展的今天，有限元法得到了广泛的应用。 要得到大挠度梁的非线性有限元方程必须从能量原理出发，通过能量原理得到 结点位移与刚度矩阵的关系来建立非线性有限元方程。就目前的非线性分析方法来 说，连续体的能量增量方程能够很好的描述物体的连续变形过程，并且精度高，理 论成熟，应用广泛。几乎所有的非线性有限元方程都是从能量增量方程中推导出来 的。其基本思路从已知构型出发，给定一个加载步长，根据已知构型的状态参数： 应力、应变、位移，得出下个构型的状态参数，通过连续不断的加载过程最终达 到求解的目标。增量非线性有限元方程中最常用的列式就是t l 和u l 列式。 t l 法又称全拉格朗日法( t o t a ll a g r a n g i a nf o r m u l a t i o n ) ，t l 法在物体的整个变 形过程中始终以初始构型作为参考构型，所以在变形中应当采用k i r c h h o f f 应力和 g r e e n 应变。 u l 法又称更新的拉格朗日法( u p d a t e d l a g r a n g i a nf o r m u l a t i o n ) ，u l 法在变形 过程中始终以相邻的构型作为参考构型，因此在变形过程中应当采用c a u c h y 应力 和a l m a n s i 应变。 大挠度空间粱的静、动力学建模、分析与计算 采用什么样的构型作为参考构型是tl 法和ul 法的重要差别，而且由于两种 列式所用到的应力和应变都不同，这两种有限元列式方法所采用的本构方程也是不 4 样的。 但两者所使用的参考构型不管是初始构型还是当前构型，其构型始终是以已知的 构型作为参考构型，所不同的是t l 法中参考构型始终是不变初始构型的，而ul 法所采用的参考构型随着求解不断变化，正是由于上面的差别，t l 法更便于理解， 而u l 法的优点在于可以使公式、数学推导更简化。 b a t h e 3 引、m e e k 3 4 1 、n a r a y a n a n 3 ”、吕和祥2 1 等通过对几何非线性梁单元的研究 表明u l 法比t l 法在计算上更有效。 本章先从三维连续体变形的应变能出发，推导了空间梁结构的受轴向拉压，弯 曲，剪切以及扭转等普遍情况下的几何非线性有限元列式的两种列式。并且在t l 法有限元列式过程中保留了传统的t l 法在推导过程中略去的高阶余量，从而获得 了在较大加载步长的情况下的解。 3 1 梁单元的应变能 参考前一章式( 2 15 ) 以及式( 2 ，2 0 ) 可得 = 妻+ ；毫) 2 + ；毫) 2 + ；c 7 2 ，营2 一一劬窘嘞声均万c 3 2 w ：娑( 3 i ) a 击 2 节蔷 式( 3 1 ) 就是空间梁结构的受轴向拉压，弯曲，剪切以及扭转等普遍情况下的 几何非线性应变方程。在前一章里我们已经讨论过该方程中各个项次的物理意义， 在此不再做解释。 假设梁单元的位移函数为： ”= u a d 。 v = v a d 。 w = w a d 。 咖= a d ! ( 3 2 ) 南京航空航天大学硕十学位论文 上式中的删、v a 、w a 、删分别为单元节点位移的位移函数“，v ，w ，的型函 数：d 。为单元的节点位移；矩阵a 为一与型函数系数有关的常数矩阵； 上式中的u 、y 、o 可以根据计算所需要达到的精度以及计算机运算机时 来具体选择。 分别简记。，、占、s ；f 为s l 、s 2 、屯。将( 3 2 ) 代入( 31 ) 得： 铲b l 、d 。+ 专即。枷c ，d e 5 2 = b l 2 d 。( 3 3 ) q = b l 3 d 。 其中 l b l = u a 一7 1 v a q w a b = w a d 。矿_ + w a d ，陟么+ ( 7 72 + f2 ) m _ 以中_( 3 4 ) i 反= 9 v a d 。中一r l w a d 。o a 导：母_(35)b 【3 = r o a “7 上式中各项的意义说明如下： 吼，吼：，眈，