（工程力学专业论文）广义精细积分法及其应用.pdf

中山大学硕士学位论文 中山大学硕士学位论文：广义精细积分法及其应用 专业：工程力学 硕士生：林敬华 指导教师：富明慧教授 摘要 本文提出了一种求解非齐次常微分方程组特解的广义精细积分法，并将其应 用于非线性动力问题、薄板弯曲问题和随机振动问题的分析中。首先，针对非齐 次项为多项式、指数函数、三角函数以及它们乘积的情况，选取一个d u h a m e l 积分形式的特解，并将积分区段f 划分为2 份。然后给出积分在一个微小区段 r 2 上的近似值，并建立一种递推关系，求出积分在区段f 2 。1 上的值，如此 类推，只需迭代次就能精确地计算出积分和，从而得到非齐次方程的一个特 解。该方法将通解的精细积分和特解的精细积分两个过程有机地结合起来，形成 一种广义精细积分法，因为后者能充分利用前者的中间结果，所以这种结合具有 极高的效率。与其他特解精细积分法相比，广义精细积分法保持高精度的同时， 既避免了矩阵求逆又给出了统一的计算公式。由于任意形式的非齐次项通常可拟 合为样条函数或者f o u r i e r 级数，因而广义精细积分法是普遍适用的。算例结果 证明广义精细积分法的有效性。 对于非线性问题，本文的做法是将非线性项纳入到非齐次项，然后对待求量 v 在时间段【f l ，f + ，】进行预估并用三次多项式逼近。这样就能够进一步将非线性项 也表达为三次多项式，再利用广义精细积分法求解后便得到 ，的修正值。重复上 述过程便得到一种迭代格式。算例结果表明，上述迭代方法计算精度高、收敛速 度快。 本文还将广义精细积分法应用于薄板弯曲问题和随机振动问题。薄板弯曲问 题半解析化以后就转化成了一阶非齐次常微分方程组的两点边值问题。原有的两 点边值问题精细积分法结合广义精细积分法，可以给出一种求解非齐次方程组两 点边值问题的稳定、有效的算法。虚拟激励法将随机振动问题转化为求解结构在 确定性荷载作用下的响应问题，结合广义精细积分法，可以使虚拟激励法的计算 效率得到很大的提高。 关键词：精细积分，指数矩阵，非线性，结构动力学 中山大学硕士学位论文 t i t l e ：g e n e r a l i z e dp r e c i s ei n t e g r a t i o nm e t h o da n di t sa p p l i c a t i o n s m a jo r ：e n g i n e e r i n go fm e c h a n i c s n a m e ：l i nj i n g h u a s u p e r v i s o r ：p r o f e s s o rf um i n g h u i a b s t r a c t t 1 1 i sd i s s e r t a t i o np r e s e n t st h eg e n e r a l i z e dp r e c i s ei n t e g r a t i o nm e t h o d ( g p i m ) t o o b t a i nt h ep a r t i c u l a rs o l u t i o nt of i r s t - o r d e rn o n - h o m o g e n e o u so r d i n a r yd i f f e r e n t i a l e q u a t i o n sw i t hc o n s t a n tt o e m c i e n t sw h e nt h en o n - h o m o g e n e o u st e r m sa r eo ft h e f o l l o w i n gf o r m s ：p o l y n o m i a l s ，t r i g o n o m e t r i cf u n c t i o n s ，e x p o n e n t i a lf u n c t i o n sa n d t h e kp r o d u c t so fe a c ho t h e r f i r s t l y , g p i mc h o o s e st h ep a r t i c u l a rs o h t i o nf o r mo f d u h a m e l si n t e g r a la n dp a r t i t i o n st h ei n t e r v a lo fi n t e g r a lri n t o2 “s u b - i n t e r v a l s o v e r t h es h o r t e s ti n t e r v a l 0 , r 2 】，t h ed u h a m e l si n t e g r a lc a nb ec a l c u l a t e da p p r o x i m a t e l y b yt a y l o re x p a n s i o n u s i n gt h ea p p r o x i m a t i o n , t h er e s u l to fd u h a m e l si n t e g r a lo v e r t h ei n t e r v a l 0 ，f 2 。1 c a nb eo b t a i n e d a n dt h ed u h a m e l si n t e g r a lo v e r o ，f 2 n - i 】 c a nb eo b t a i n e d b yu s i n gt h ei n t e g r a lo v e r 【o ，r 2 n q 】s ot h ec o m p u t a t i o no f d u h a m e l si n t e g r a lo v e rt h ei n t e r v a l 0 ，f 】c a nb ec o n v e r t e di n t oar e c u r s i o n i nt h i s d i s s e r t a t i o n , t h ep r o c e s so fp r e c i s ei n t e g r a t i o nm e t h o d ( p i m ) f o rg e n e r a ls o l u t i o na n d t h eo n eo fg p i m 南rp a r t i c u l a rs o l u t i o na r ec l o s e l yi n t e g r m e d f o rt h el a t t e rc a nm a k e 如l lu s eo ft h em i d d l er e s u l to ft h ef o r m e r ，t h i sc o m b i n a t i o ni so fah i g hc o m p u t a t i o n a l e f f i c i e n c y c o m p a r e dw i t ho t h e rm e t h o d sc o m p u t i n gd u h a m e l si n t e g r a l ，g p i mi sn o t o n l ya sp r e c i s ea sp i mb u ta l s oa v o i do fm a t r i xi n v e r s i o n a n dag e n e r a lf o r m u l ao f g p i mi s g i v e n t o g r e a t l y s a v e c o m p u t a t i o n t i m e f o ro t h e rf o r m so f n o n - h o m o g e n e o u st e r mc a l lb ec o n v e r t e di n t oap o l y n o m i a l ( s p l i n ef u n c t i o ni s s u g g e s t e d ) a n dt h et r i a n g u l a rs e r i e s g p i mi ss t i l la v a i l a b l e n u m e r i c a le x a m p l e sa r e g i v e nt od e m o n s t r a t et h ev a l i d i t ya n de 伍c i e n c yo fg p i m a n di nt h ed i s s e r t a t i o n , g p i mh a sa l s ob e e nu s e dt os o l v es o m ep r o b l e m s ，s u c h a sn o n l i n e a ro r d i n a r yd i f f e r e n t i a le q u a t i o n s ，t h i np l a t eb e n d i n ga n dr a n d o mv i b r a t i o n c o m b i n i n gg p i ma n di t e r a t i v em e t h o d ah i g hp r e c i s ea l g o r i t h mi sg i v e nt os o l v e n o n l i n e a re q u a t i o n s f o rn o n l i n e a rp r o b l e m , t h ee q u a t i o n sc a nb e e q u i v a l e n t l y s e p a r a t e di n t ol i n e a ra n dn o n l i n e a rp a r t sa n dt h en o n l i n e a rp a r t sc a nb es e e na s n o n h o m o g e n e o u st e r m s t oa p p r o x i m a t et h eu n k n o w nf u n c t i o ni nt h en o n l i n e a rt e r m s ， ac u b i cp o l y n o m i a li n t e r p ol a t i o nc a nb ec o n s t r u c t e db yu s e do ft h ef u n c t i o nv a l u e s a n dt h ed e r i v a t i v ev a l u e sa tt h eb o t he n d so ft h ei n t e r v a l ( v a l u e sa tt h ef o r m e re n da r e k n o w n , b u tv a l u e sa tt h eo t h e re n da r en e e d e dt ob ep r e d i c t e d ) u s i n gt h e a p p r o x i m a t i n gc u b i cp o l y n o m i a l s ，w ec a nc o n s t r u c ta n o t h e rc u b i cp o l y n o m i a l i n t e r p o l a t i o nt oa p p r o x i m a t et h en o n l i n e a rp a r t s ，t h e n o n h o m o g e n e o u st e r m s s ot h e l 广义精细积分法及其应用 n o n l i n e a rp r o b l e mi sc o n v e r t e dt ol i n e a rp r o b l e mw i t hn o n - h o m o g e n e o u si t e m sw h i c h c a nb es o l v e db yg p i m t i l i si sa ni t e r a t i o np r o c e s s n u m e r i c a le x a m p l e sa r eg i v e n a n ds h o wt h a tf o r t h en o n l i n e a rp r o b l e m , t h em e t h o di sf a s tc o n v e r g e n ta n dh i g h p r e c i s i o n i na d d i t i o nt on o n l i n e a rp r o b l e m , t h i sd i s s e r t a t i o na l s oa p p l i e sg p i mt ot h i np l a t e b e n d i n gp r o b l e ma n dr a n d o mv i b r a t i o np r o b l e m t h i np l a t eb e n d i n gp r o b l e mc a nb e c o n v e r t e di n t oas e m i a n a l y t i c a lf o r mw h i c hi sp r e s e n t e db yt w o p o i i l tb o u n d a r yv a l u e p r o b l e m s ( t p b v p s ) 。c o m b i n i n gg p i ma n dp r e s e n t e dp i mf o rt p b v p s ，t h i s d i s s e r t a t i o n g i v e s a ne f f i c i e n ta n ds t a t i c m e t h o dt os o l v et p b v p so f n o n h o m o g e n e o u se q u a t i o n s p s e u d o - e x c i t a t i o nm e t h o d ( p e m ) c a nt r a n s f o r n lr a n d o m v i b r a t i o np r o b l e mi n t os t r u c t u r a ld y n a m i ca n a l y s i sw i t hi n i t i a lv a l u e s c o m b i n e d g p i ma n dp e m ，a ne f f i c i e n c ys c h e m ei sp r e s e n t e dt og r e a t l yi n c r e a s et h ee f f i c i e n c y k e yw o r d s ：p r e c i s ei n t e g r a t i o nm e t h o d ；e x p o n e n t i a lm a t r i x ；n o n f i n e a r ； s t r u c t u r a ld y n a m i c s ； i v 论文原创性声明内容 本人郑重声明：所呈交的学位论文，是本人在导师的指 导下，独立进行研究工作所取得的成果。除文中已经注明引 用的内容外，本论文不包含任何其他个人或集体已经发表或 撰写过的作品成果。对本文的研究作出重要贡献的个人和集 体，均已在文中以明确方式标明。本人完全意识到本声明的 法律结果由本人承担。 学位论文作者签名： 日期：础年乡月多 曰 b 学位论文使用授权声明 本人完全了解中山大学有关保留、使用学位论文的规 定，即：学校有权保留学位论文并向国家主管部门或其指定 机构送交论文的电子版和纸质版，有权将学位论文用于非赢 利目的的少量复制并允许论文进入学校图书馆、院系资料室 被查阅，有权将学位论文的内容编入有关数据库进行检索， 可以采用复印、缩 学位论文作者签名 日期：矽p 年乡+ 月 法保存 师签名 期：砷年多月蚕如 中山大学硕士学位论文 第一章绪论 1 1 结构动力方程直接积分法简介 任何实际的结构都是连续的，基于连续模型建立的运动方程是偏微分方程。 但在实际计算中，由于模型和边界条件的复杂性，往往需要采用数值的方法，对 实际结构进行离散近似，例如质点弹簧系统；或者采用有限元法对偏微分方程进 行离散。离散的结果是，得到如下形式的方程组： 膨+ 仪+ k x = 厂( 1 - 1 ) 称为结构动力学方程。这是一个二阶常微分方程组，常规可使用r u n g e k u t t a 法求解。但考虑到结构动力学方程自身的特点，人们提出了另外一些具有针对性 的方法。这些方法分为两类：振型叠加法( m o d es u p e r p o s i t i o nm e t h o d ) 和直接 积分法( d i r e c ti n t e g r a t i o nm e t h o d ) 。振型叠加法是先求解原方程的广义特征值问 题，得到频率和振型，通过模态振型对原方程进行解耦，得到n 个独立的方程再 进行求解。振型叠加法求得的解关于时间是解析的，在时间域上均满足方程( 1 1 ) 。 但实际问题中，由于矩阵阶数n 一般比较大，求解全部的频率和模态并不实际， 而高阶模态对结构响应的影响有限，所以往往只需要近似求解前几阶模态；直接 积分法是通过在时间域上的差分离散，得到一系列只在离散的时间点上满足方程 ( 1 1 ) 的解答。直接积分法的格式比较多，有中心差分法、h o u b o l t 法、线性加速 度法、w i l s o n - 9 法、n e w m a r k 法等1 。在采用各种格式时，需注意算法的精度、 稳定性和数值阻尼等问题。 1 2 精细积分法简介及研究现状 为求解一阶常系数线性常微分方程组 多= h v + ，( f ) ( 1 2 ) 必须计算指数矩阵e 村。一般来说，求解该矩阵需要对系数矩阵做相似变换，这 样做效率并不高。而对于现代最优控制，需要快速高效又满足精度要求的算法。 于是v a nl o a n 2 1 、m o l e r 和v a nl o a l l 【3 11 9 7 8 研究了计算指数矩阵的1 9 种算法。 1 9 9 4 年钟万勰提出了精细积分算法【4 】，这是一种2 算法，精度非常高，因而圆 满解决了方程( 1 2 ) 对应齐次方程初值问题的通解问题，并为求解方程( 1 2 ) 提供了 途径。 能够高精度求解方程( 1 2 ) 有很大的意义。热传导问题和现代控制理论本身就 具有方程( 1 2 ) 的形式。对于结构动力问题，实际工程结构都是连续的、复杂的， 必须离散，离散后的方程一般是二阶常系数常微分方程，由勒让德变换则可化成 方程( 1 2 ) 的形式。除此之外，弹性力学的很多问题，例如规则区域的平面问题( 如 矩形域、圆形域和环形域等) 、多层层合板圣维南问题等，在l a g r a n g e 体系下， 广义精细积分法及其应用 基本变量是位移，应变则是位移的偏导数，控制方程是偏微分方程；而在h a m i k o n 体系下，位移和应力是基本变量，若选定一个方向为解析方向，在其他方向上进 行离散，用h a m i l t o n 正则方程表示的控制方程则是一阶常微分方程。 精细积分法自出现以来，在结构动力分析中获得了日益广泛的应用。精细积 分法在求解动力方程的齐次方程时具有突出的优势，一方面它既有极高的精度， 可认为是计算机上的精确解，另一方面该求解过程稳定性好，因此具有较高的效 率。张洪武h 1 研究了迭代次数和初始阶段的截断误差对精细积分法结果精度的影 响，并给出了给定误差的条件下，确定t a y l o r 展开项数和迭代次数的算法。谭述 君等一1 对使用p a d e 展开的精细积分法进行了误差估计，并给出了给定精度时p a d e 展开项数和迭代次数的自适应算法。曾进等1 将精细积分法引入到辛算法中。徐 明毅等一1 给出了精细辛算法的误差估计。 对于非齐次结构动力方程，钟万勰【l l 】将激励在每个时间步中假设为线性变 化，并得到了求解系统瞬态响应的精细积分公式。林家浩，钟万勰等【1 2 】进一步推 广了精细时程积分方法，考虑了非齐次项是三角函数时的情况，并得到了h p d s ， h p d f 计算格式，其中h p d f 计算格式是基于f o u r i e r 级数展开和常微分方程的理 论解建立的精细算法。林家浩又将精细积分方法的思想应用到分层岩层中的平稳 月乍平稳地震波的传播等随机振动领域中“。 上述这些方法求解特解，简单易行，便于实现，但会涉及矩阵求逆的问题 1 4 - 1 5 1 。矩阵求逆不仅计算量大，有时甚至是不稳定或者不可能的，这也在一定程 度上限制了上述精细积分法的应用。为了避免矩阵求逆，任传波1 2 l j 和汪梦甫【2 2 】 提出了利用数值积分求解d u h a m e l 形式特解的方法，这种方法具有简便易行的特 点，其精度会随选取积分点个数的增加而逐渐提高，但与此同时计算量也随之急 剧增加，实际上特解要想达到或者接近通解的精度几乎是不可能的。 针对求逆和精度之间的矛盾，顾元宪 1 6 , 1 7 提出了增维精细积分的方法，该法 的特点是先将非齐次方程转化为齐次方程，然后再进行精细积分，从而避免矩阵 求逆，且具有通解精细积分的精度。但这种转化是有条件的，非齐次项必须是多 项式、三角函数、指数函数或者三者乘积的形式才可以增维，且当非齐次项形式 较为复杂时，增维的技巧性大大提高，缺乏通用的增维方法，而且不同的增维方 案计算量相差很大。谭述君等叫和富明慧等1 通过将通解精细积分的思想引入 至u d u h a m e l 形式特解的求解中，得到了一种同样适用于非齐次项为多项式、三角 函数、指数函数或者三者乘积的高精度求解方法。富明慧等卜w 通过对一类指数矩 阵函数进行研究，对于上述非齐次项，给出了通用的求解公式。 对于非线性动力学问题的研究，有解析法卜和数值法。近年来精细积分法 的发展，为运用高精度数值方法研究非线性动力学问题提供了一个新的计算平 台。求解该问题的思路时是：把非线性动力学方程写成痧：h v + f ( v , f ) 的形式。其 2 中山大学硕士学位论文 中h 为常数矩阵，f ( v ，f ) 为方程的非线性部分。该方程的解可以写成通解加上 d u h a m e l 积分的形式，再运用逐步积分法，在各个时间段内将该方程化为非齐次 线性常微分方程组求解。具体有一下几种做法：赵秋玲卜u 运用分段直接积分法， 将非线性项，( p ，f ) 在k 展开成一次多项式，然后做分部积分。这样在该时刻方程 就化为了一个一阶非齐次线性常微分方程组。因为( ，f ) 包含未知向量p ，所以 用d u h a m e l 积分求特解时必须对步长内的v 进行预估，如何预估成为关键。裘春 航p 峰”等对该法有所发展，可以将( ，f ) 在缸展开成高次多项式，并提出避免求 ( ，r ) 的各阶导数的方法，但该方法仍然需要对矩阵日求逆。李伟东等p w 建议使 用修正刚度的方法来避免h 阵无逆。张素英p ”给出了利用增维方法和在时间步 长内假设y 为常数的求解方法，但该方法在每个步长内都需要重新计算一个传递 矩阵。张洵安畔1 给出了用迭代法求解非线性动力学方程的方法。梅树立p 叫结合 使用同伦摄动法和精细积分法，将非线性动力学方程做同伦摄动展开得到一系列 线性微分方程组，用精细积分法求解这些方程组即可得到各阶近似解。李金华【3 4 l 等通过在时间步长内将指数矩阵用t a y l o r 展开式近似，给出了一种避免矩阵求逆 的方法。向宇等p 叫将非线性项做二阶t a y l o r 展开，并结合增维精细积分法给出 一种避免矩阵求逆的算法。富明慧等w 在文献【1 9 ，2 0 】的基础上给出一种将非线性 项用三次样条近似的预估校正格式。张素英等卜叫提出了一种精细r u n g e k u t t a 法，它结合了精细积分法和r u n g e k u t t a 法的优点，因此是一种好的构思。但由 于该文中仅是将r u n g e k u t t a 法的有关公式直接移植到d u h a m e l 积分中，忽略了 动力学方程的特点，因此精度和效率均不十分理想。富明慧等卜叫根据动力学方程 的特点，重新建立了精细r u n g e k u t t a 法。它不是将r u n g e k u t t a 法的公式简单地 移植到特解的数值积分中，而是利用r u n g e k u t t a 法的几何含义，对数值积分点 处未知的状态参量进行预估，从而得到了精度更高的精细r u n g e k u t t a 法。王海 波等p q 建议，在处理非齐次项时做进一步细分，把线性部分和非线性部分分开处 理。 上面介绍的都是初值问题，对于边值问题，钟万勰p 1 w 在每一个小步长内将 边值问题转化成初值问题的形式，使用合并消元法，得到了与初值精细积分法具 有相同精度的边值问题精细积分法。但应该指出，这种方法只能处理某些特定的 边界条件，并且在求时间域内部的值时常常会出现严重的数值病态问题。 除了上面介绍的常微分问题外，精细积分法在求解偏微分方程中也得到非常 多的重视。近年来，学者们相继将精细积分法应用于求解热传导问题p 6 州1 、水动 力学问题p 1 4 如、随机振动问题h 3 1 、椭圆积分问题1 4 4 乃至求解薛定谔方程h 5 1 。使 用精细积分法求解偏微分方程的个共同特点是使用状态空间法，将原问题半解 析化，变成常微分方程来处理。 广义精细积分法及其应用 1 3 板壳有限元简介 板壳的弯曲问题是一个古老的且到现在还非常活跃的课题。早在十九世纪初 期，力学家们已经开始对板壳结构进行了研究。板和壳可以看成是三维实体的一 种特殊形式，这种结构一个方向的尺寸比其他两个方向要小很多，利用这一特点 可以实现对原问题的消元和降维。薄板理论就是这样一种方法。该理论是以1 8 5 0 年k i r c h h o f f 提出的假设为基础发展起来的m 一。基于该假设，控制方程最后可 约化到一个只关于横向位移( 即挠度) 的非齐次双调和方程。 到了2 0 世纪6 0 年代，有限元发展之初，板壳弯曲问题便成为有限元法最先 应用的领域之一。如果按照完全三维问题来处理，不仅计算量大，而且当网格一 个方向的尺寸比另外两个方向要小很多时，会引起严重的数值病态。因此，基于 薄板理论构造板单元便成为一个自然而然的选择。 薄板理论用挠度w 的一阶导数来近似直法面的转角，因此为了使位移协调， 则挠度w 及挠度的斜率娑和譬必须连续，称为c l 连续。然而，在有限元中要 g k 砂 求形函数满足c 1 连续条件非常困难。但如果仅要求挠度w 及其斜率在节点上满 足连续性条件，而允许挠度的斜率在单元交界处不满足协调条件，那么这样的形 函数还是比较容易构造的。这种单元称为非协调元。文献【4 6 】认为，若这种单元 满足“分片试验”，其收敛性还是能够得到保证。由于协调元离散后得到的刚度 阵会比真实情况大，非协调元在一定程度上缓解了这种刚度过大的问题，因此非 协调元的精度一般比协调元要好r 。除了非协调元外，还有许多方法可以避免形 函数的c l 连续条件。例如采用混合单元，离散k i r c h h o f f 理论( d k t ) 单元和应力 杂交元等。 混合单元是由加权余量法或h e l l i n g e r - r e i s s n e r 变分原理给出的保留弯矩m 和挠度w 而消去转角口的泛函表达式建立得到的。此时我们可以发现，无论是弯 矩m 还是挠度w 的插值函数都仅要求c 0 连续。 另一种避免c l 连续条件的方法是离散k i r c h h o f f 理论( d k t ) 的薄板单元。该方 法是，首先在比较多的节点上，对挠度w 和转角口进行独立的插值，通过配点法、 配线法或子域法等离散约束方法，使挠度转角关系式秒一v = o 在一些离散区域 内得到满足，并消去一些节点自由度。这样对挠度w 和转角口的插值函数也仅需 要c 0 连续。由于单元形状和离散约束方法不一样，研究者们提出了很多不同的 d k t 单元。i r o l l s 5 0 1 提出了s e m il o o p 单元；龙驭球等郴1 根据相似的原理提出了 广义协调元的概念。 应力杂交单元是基于修正最小余能原理得到的，即使用l a g r a n g r e 乘子将单 4 中山大学硕士学位论文 元边界的平衡条件加入到最小余能原理中，且该l a g r a n g r e 乘子就是单元边界的 位移。单元的边界节点上仅有位移分量，边界上的位移通过边界节点的位移插值 得到；而应力部分不是像传统有限元那样用节点应力插值得到( 此时节点上也没 有应力分量) ，而是先假设应力的分布形式。所以关于应力的未知量并非节点应 力，而是应力分布表达式的系数。虽然应力杂交元也同时包含应力和位移，但与 混合单元中应力和位移同时存在于单元内部和交接面不同，应力杂交元中的位移 只存在于单元的交接面。 以上论述的是基于k i r c h h o f f 假设的薄板单元。对于中厚板，直法线假设不 再适用，应该代之以r e i s s n e r - m i n d l i n 假设。为了考虑k i r c h h o f f 薄板理论没有考 虑到的横向剪应力，r e i s s n e r 和m i n d l i n 假设弯曲前中面的法线在弯曲后仍然保 持为直线，但并不垂直于中面。此时原中面法线的转角是一个独立的量，并不等 于挠度的导数。因此面内位移可表示成该点的转角与厚度方向坐标z 的系数( 坐 标原点在中面上) 乘积的形式。由此可知，一阶剪切变形理论的横向剪应力是常 数，并非由三维理论得到的抛物线分布，因此在考虑剪切变形能时要考虑剪切修 正因子。应该指出，离散约束的方法不仅可以用在薄板单元，还可以用在厚板单 元。同样对挠度m 和转角p 进行独立的插值，与薄板情况不同的是，厚板时是对 】，= p v 国而不是口一v 国= 0 进行离散约束r ”。 随后w h i t n e y 等p 叫又发展了高阶剪切变形理论，即假设面内位移中除了包含 z 硕外，还有z 2 和z 3 项等。该理论是对r e i s s n e r - m i n d l i n 理论的发展。因为多了 自由度，因此计算量也增加。r e d d y p 纠发展了一种简单的三阶板理论，利用剪应 变在板的上下为零的条件，使自由度数与普通的一阶剪切变形理论一样。高阶剪 切变形理论的横向剪应力是呈抛物线分布，因此不需要考虑剪切修正因子。 此外，还有基于半解析法建立的有限条法。该方法是将一个方向做f o u r i e r 展开，其他方向用普通的多项式形函数离散得到的。该方法可以用于三维分析也 可以用于平面问题，当然也可以用于板壳问题h “。由于三角函数的正交性，f o u r i e r 展开的每一项的系数，即未知量都是非耦合的。因此该方法具有很高的计算效率。 该方法可以用于承受非对称载荷的轴对称实体r 。但由于方法本身的局限，完全 解耦的有限条法只能用于规则形状的构件。后来又发展了不完全耦合的有限条法 l 3 , f i 。该方法与半解析法不同的地方在于，解析方向并不用三角级数展开，而采用 其他级数展开或者用样条插值。而样条有限元方法可以用来分析任意形状的板问 题【5 5 】o 常用的工程结构除了板外，还有壳体结构。壳是板的衍生结构。与板相似， 壳也是其中一个方向的尺度远小于其他两个方向；但与板不同的是，板的中面是 平面，而壳的中面是曲面，因此，壳的受力形式也与板不同。板只承受垂直于中 面的力，而壳所受的力除了有中面法线方向的分量外，还有与中面相切的分量。 5 广义精细积分法及其应用 因此，与板相比，薄壳需要有三个方向的位移；而厚壳的位移模式除了有两个方 向的转角和垂直中面的挠度外，还有与中面相切的两个位移。 正因为板与壳有很多相同之处，因此，可以用平板来代替壳。如果单元划分 足够密的话，计算总会收敛。这种壳单元称为平板壳或“折板 单元。考虑到平 板壳单元与其他单元相连接时，需要使用到中面法线的转角这一自由度，称为旋 转自由度( d r i l l i n gd e g r e e ) 。因此，可以将原来每节点5 个自由度增加到6 个。 由此带来的问题是，由于没有从物理上去考虑旋转刚度，因此局部坐标上的单元 刚度阵对应旋转自由度的那一行全为0 ，即该刚度矩阵是奇异的。其中一种解决 方法是增加转动刚度系数p 。 除了平板壳外，还可以建立从三维理论蜕化而来的，中面为曲面的曲边壳单 元一。轴对称壳受轴对称荷载的情况，可将单元简化成一条线。跟前面的情况一 样，轴对称壳也有“薄与“厚，以及直线单元和曲线单元之分一一”。对于某 些受非轴对称荷载的轴对称壳可用前面提到的半解析法求解。 1 4 本文的主要工作及论文章节安排 本文主要研究了一阶非齐次常系数常微分方程组的特解精细积分法，并将该 方法应用到非线性常微分初值问题、薄板弯曲问题和虚拟激励法的求解中。 本论文的章节安排如下： 第一章是绪论，介绍了本文工作的研究背景。首先介绍了传统的求解结构动 力学方程的直接积分法及其优缺点。其次，介绍了精细积分法的发展过程及研究 现状。本章还简单介绍了板壳有限元法及半解析法的研究现状。 第二章介绍了五种传统的求解结构动力方程的直接积分法一中心差分法、 平均加速度法、线性加速度法、n e w m a r k 法和w i l s o n - 曰法，并分析了每种算法 的精度和稳定性问题。 第三章详细介绍通解的精细积分法，并分析了该方法的精度和收敛机制。本 章还介绍了精细积分法参数的自适应选择的方法。 第四章是本文最主要的工作，即一种求解特解的广义精细积分法，计算过程 充分利用了通解精细积分的中间量，达到了通解精细积分的精度，并避免了矩阵 求逆。当非齐次项为多项式、三角函数、指数函数及三者相互乘积的情况，对 d u h a m e l 积分进行分部积分处理，大大的减少了计算量并且给出了通用的求解公 式。 第五章是本文另一个主要工作，即将广义精细积分法应用于非线性常微分方 程的求解。对于非线性问题，本文的处理方法是将非线性项看作是非齐次项来处 理。先对未知函数进行预估，然后在求解区域内用分段三次多项式去逼近非线性 项，转化成广义精细积分可以处理的形式。结合预测一校正的计算策略给出一种 6 中山大学硕士学位论文 迭代算法。 第六章简单介绍了一阶常系数常微分方程组的两点边值问题的精细积分法。 第七章和第八章也是本文的工作，即将广义精细积分法应用于薄板弯曲问题 和虚拟激励法的计算中。 7 广义精细积分法及其应用 第二章结构动力方程直接积分法 结构动力方程可表示成式( 2 1 ) 的形式 膨+ q + 戤= 厂 ( 2 一1 ) 本章介绍求解结构动力方程的中心差分法、线性加速度法、w i l s o n - o 法、n e w m a r k 法，并讨论这几种直接积分算法的精度、稳定性和数值阻尼等问题。 2 1 中央差分法 用步长址将时间离散，得到各个时间点= a t ，其中i = 0 ，1 ，2 ，。若要求结 构动力学方程( 2 - 1 ) 在离散的时间点上满足，则有： 憾+ c 墨+ q = e ( 2 2 ) 其中x ，x ，鼍和z 是分别是时刻上的位移，速度，加速度和荷载向量。若 把+ 和t 一。时刻的位移向量在时刻做t a y l o r 展开得到 x + i = - x + x 出+ 兰x ，2 + 丢x 出3 + o ( x 4 ) ( 2 - 3 ) x - i 一- - x 一毫址+ l 戈, a t 2 一吉x r 3 + d ( 址4 ) ( 2 - 4 ) 把以上两式分别相加和相减后得到 定，：丢( k + 。一墨一。) + d ( 出3 ) ( 2 - 4 ) x f 2 = ( k + l 一2 鼍+ 置一1 ) + d ( 出4 ) ( 2 5 ) 将时间步长址和址2 除到右边后得到 x = 去( x + 。一k - 1 ) + d ( 出2 ) ( 2 - 6 ) 鼍= - - 专( x , - 2 x , + x 一。) + o ( x 2 ) ( 2 - 7 ) 忽略高阶小量o ( a t 2 ) 后，即有 x = 五1 ( 鼍+ 。一x 一，) ( 2 8 ) 8 中山大学硕士学位论文 , y , = 1 - - - t ( x j q - 2 x i + x i 。 这里可以看出中心差分格式是二阶精度格式， ( 2 - 9 ) 代入方程( 2 - 2 ) ，最p # 4 - ( 2 - 9 ) 即误差在出2 量级。将式( 2 8 ) 和式 古m + 西1c ) = f 一( 置一舌m ) 置一( 古肘一西1 c ) 扎( 2 - 1 0 ) 该方程可写为 衄+ 。= t ( 2 1 1 ) 其中， 府：三m + 上c t j 三越 霉= e 一( 置一矿2m ) x 一( 古m 一去c ) x i 一。e = e 一( 置一矿m ) x 一( 寿m 一击c t 这是一个显式递推格式，只要给定初始值，则往后离散点上的值可由上式递 推求得。值得一提的是，这个算法存在如何起步的问题，这里提供一个简单的方 法。由于初始条件己知，即和毫已知，可代入方程( 2 2 ) 求得藏。再将一缸时 刻( 假设初始时刻从0 开始) 的位移向量在0 时刻做t a y l o r 展开得到 舡。嗡一毫血+ 去甍出2 ( 2 一1 2 ) 下面是利用中心差分法求解结构动力方程的算法【l 】： i 初始化阶段 ( 1 ) 形成刚度阵置，质量阵肘和阻尼阵c ； ( 2 ) 由初始条件k 和x o 解出x o ； ( 3 ) 选取时间步长出，计算出积分常数= f z 、q = f 、巳= 2 c o 、 白5 必； ( 4 ) - t 篑tx _ 。鹊一毫出+ 去藏址2 ； ( 5 ) 形成有效质量阵庙= 古m + 去c ； ( 6 ) 对有效质量阵廊做三角分解：府= l d l r 2 时间步迭代阶段 9 广义精细积分法及其应用 ( 1 ) 计算t 时刻的有效荷载向量 丘= e 一( 眉一矿2m ) x 一古m 一去c ) x “z = e 一( 眉一矿m ) x 一嘻m 一击c ) x “ ( 2 ) 由方程工时墨+ 。= t ，求解o 。时刻的位移向量置+ ， ( 3 ) 若有需要，可计算t 时刻的加速度和速度向量 再对方程( 2 1 4 ) 使用中心差分法。得到递推关系式 = x - 1 1 小止降一事卜 1 0 中山大学硕士学位论文 d ”：p j ”p - 1 ( 2 1 8 ) 因为此时d 的阶数是2 ，由线性代数理论司知只能存在两个一彤卜约当块或者一个 二阶约当块。所以，当d 有两个互异特征根时，其形式为： ，= ( 素羔) 所以，此时 儿一晏 定义p ( d ) = m t = a x l , 2h i 为递推矩阵的谱半径。若要求d ”收敛，即中心差分算法稳定， 要求p ( d ) 1 。 当d 有重根时，其形式为： t ，= ( 言羞) 当存在两个一阶约当块时，口= 0 ；当存在一个二阶约当块时，口= 1 。此时，有 卜k 口列 若口= 1 ，则要求p ( d ) 4 厂+ 争2 y 虿1 ，厂一互1 综合式( 2 3 9 ) ，( 2 - 4 0 ) ，可知，n e w m a r k 法无条件稳定的条件是 厂虿1 ，丢( 7 + 丢) 2 而文献【5 9 】则分析了带阻尼的情况，给出了递推公式 褂x i + l 篓蚓 白= 丢一一五1 一岁) a p _ 2 ( 1 训肚 乞= 1 7 一( j 1 一) 胆一2 ( 1 一厂) 肛 毛= 弋1 一) 口一2 ( 1 7 弦 口= ( 壶+ 百2 4 r + 分1 j i ：玺 n ( 2 3 6 ) ( 2 3 7 ) ( 2 3 8 ) ( 2 - 3 9 ) ( 2 - 4 0 ) ( 2 - 4 1 ) ( 2 - 4 2 ) 广义精细积分法及其应用 文献【5 9 】据此得出n e w m a r k 法无条件稳定的条件是 2pyi1(2-43) 2 3w i l s o n 一0 法 w i l s o n - 口法实际上是线性加速度法的推广。线性加速度法假定加速度从时刻 t 到时刻0 。即f j + ，时刻为线性变化，而晰l s o n _ 口法则假设加速度从时刻f j 到时 亥l j t t + 0 a ty g 线性变化，其中0 1 0 ，然后再通过插值得至l j t , + ，时刻的加速度。 令f 表示时间的增量，其中0 f o a t ，于是对从t 到时刻，+ 触的时间区间， w i l s o n - 口法假定 x 竹= x + 古( 鼍+ 触一鼍) ( 2 - 4 4 ) 对f 积分上式，并利用条件定+ 。= x 和x + 。= 置得： 墨竹= 毫+ 鼍f + 互。| 西r 2 ( x i + 触一x ) ( 2 4 5 ) x + ，= 置+ 墨f + 丢k f 2 + 意( x + 锄一警) ( 2 4 6 ) 利用式( 2 4 5 ) 和( 2 4 6 ) ，将1 取为o a t 即可求出 x + 触= 万( 鼍埘_ ：奠) 一去x 一2 x ( 2 - 4 7 ) 墨+ 细= 亳( x + 触一x ) 一2 x 0 2 a t 支，( 2 - 4 8 ) 因为假定加速度是线性变化，故所用的荷载也为线性变化，即有方程 f + 锄= e + 秒( 巧+ l f ) ( 2 - 4 9 ) 要求在时亥t jt