（固体力学专业论文）薄板受迫振动问题的Trefftz有限元法.pdf

中文摘要 作为一种高效有力的计算工具，杂交t r e f r z 有限元法自3 0 年前提出以来， 已经受到了愈来愈多的重视。该有限元模型的特点是在单元域内和边界上独立地 建立位移插值函数。在单元区域内部，其插值函数应事先满足的控制方程；而单 元边界上的插值函数满足单元间的位移协调条件。最后用变分法或其他方法将这 两组独立的位移函数联系起来，最终形成力一位移关系的有限元列式。其求解的 方程仍以节点位移作为未知量，而且单元公式中只含有边界积分，这样就可以构 造出任意多边形甚至曲边单元。该单元可以看作是一种特殊的对称形式的边界型 求解方法。因此，t r e f f t z 有限元法具备了传统边界元法的优点，而又避免了复杂 奇异积分方程的计算。采用控制方程完备解系作为单元域内插值函数与网线函数 是h t 有限元突出的特点。到目前为止，t r e f f t z 有限元法已经应用到诸多工程领 域。然而，目前将杂交t r e f f t z 有限元法应用于薄板强迫振动方面的研究还处于 起步阶段。 本文主要有两部分内容：对于t r e f r z 有限元方法本身的一些研究和基于该 有限元法建立了薄板受迫振动问题的有限元模型。 对于薄板受迫振动问题，本文基于t r e f f t z 有限元法基本理论，将无条件稳 定的平均加速度法与修正的运动方程结合，并运用第一类k 阶变型b e s s e l 函数及 第一类k 阶h a n k e l 函数来构造薄板受追振动问题的特解及齐次解( 完备解系) ，建 立薄板受迫振动问题的有限元模型，编制了相应的有限元程序，并通过几个数值 算例证明了该法的精度及效率。 关键词：t r e f f t z 有限元受追振动薄板 a b s t r a c t a sah i g h l ye f f i c i e n ta n dw e l le s t a b l i s h e dt o o l ，t h eh y b r i d t r e f f t z ( h t ) f e m ， i n i t i a t e da b o u tt h r e ed e c a d e sa g o ，h a sn o wb e c o m em o r ea n dm o r ep o p u l a r t h em o d e l u s e st w od i f f e r e n ts h a p ef u n c t i o n si nt h en o n - c o n f o r m i n gi n t r a - e l e m e n tf i e l da n dt h e a u x i l i a r yc o n f o r m i n gf r a m ef i e l d i n t r a e l e m e n tf i e l di sc h o s e ns oa st oap r i o r is a t i s f y t h eg o v e r n i n gd i f f e r e n t i a le q u e n t i o n s ，w h i l et h ei n t e r e l e m e n tc o n t i n u i t yi se n f o r c e d v i at h ea u x i l i a r yc o n f o r m i n gf l a m ef i e l d u s i n gav a r i a t i o n a lp r i n c i p l eo rs o m eo t h e r m e t h o d sm a k e st h e s et w of i e l d sc o n f o r m e d ，a tl a s to b t a i n e dt h er e l a t i o n s h i pb e t w e e n f o r c ea n dd i s p l a c e m e n t i nt h es t i f f n e s s e q u a t i o n o ft h em o d e l ，t h en o d a l d i s p l a c e m e n t sa r es t i l lt h eu n k n o w nv a r i a b l e s ，a n dt h ef o r m u l a t i o no n l ye a u sf o r i n t e g r a t i o na l o n gt h ee l e m e n tb o u n d a r yw h i c he n a b l e sa r b i t r a r yp o l y g o n a lo re v e n c u r v e s i d e de l e m e n t st ob eg e n e r a t e d a sar e s u l t ，i tm a yb ec o n s i d e r e da sas p e c i a l ， s y m m e t r i c ，s u b s t r u c t u r e - o r i e n t e db o u n d a r ys o l u t i o na p p r o a c ha n dt h u sp o s s e s s e st h e a d v a n t a g e so ft h ec o n v e n t i o n a lb e m i nc o n t r a s t ，h o w e v e r , h tf e ma v o i d st h e i n t r o d u c t i o no fs i n g u l a ri n t e g r a le q u a t i o n s t h ei n t r a - e l e m e n ts h a p ef u n c t i o na n dt h e e l e m e n tf l a m es h a p ef u n c t i o nf e a t u r et h eh tf e mm o d e l u pt on o w , i - i tf e mh a s b e e nt h o r o u g h l ye x p l o r e di nan u m b e ro ff i e l d si ne n g i n e e r i n g h o w e v e r , i ti ss t a r t i n g p r a s et od i s c u s st h ea p p l i c a t i o n so f t h eh tf m i t ee l e m e n ta p p r o a c ht of o r c e dv i b r a t i o n a n a l y s i sa n di t sa p p l i c a t i o nt ot h i np l a t e s t h em a i nc o n t e n to ft h i sp a p e ri st w of o l d ：s o m er e s e a r c ho nh tf e mi t s e l fa n d ， d e v e l o p m e n to ff o r c e dv i b r a t i o ns o l u t i o nt os h i np l a t ea l g o r i t h mb a s e do ni t t h i sp a p e ri sb a s e do nh tf e ma p p r o a c ht of o r c e dv i b r a t i o na n a l y s i so ft 1 1 i n p l a t e s t h ei n v e s t i g a t e da p p r o a c hc o n s i s t s i nt h ec o m b i n a t i o no ft h et - e l e m e n t s o l u t i o no fas u i t a b l ym o d i f i e de q u a t i o no fm o t i o nf o ra s e q u e n c eo fu n i f o r m l ys p a c e d d i s c r e t i z e dt i m e sw i t ha nu n c o n d i t i o n a l l ys t a b l ea v e r a g ea c c e l e r a t i o nm e t h o d a n dt h i s p a p e ra l s op r e s e n t st h em o d i f i e db e s s e lf u n c t i o nw i t ho r d e rka n dt h eh a n k e lf u n c t i o n o ft h ef i r s tk i n dw i t ho r d e rkt of o r mt h ep a r t i c u l a rs o l u t i o na n dt h eh o m o g e n e o u s s o l u t i o n ( t - c o m p l e t es y s t e mo fh o m o g e n e o u ss o l u t i o n ) a tt h ee n do ft h i sp a p e r ，t h e e f f e c t i v e n e s sa n da c c u r a c yo ft h ee l e m e n tm o d e li sa s s e s s e dt h r o u t hs o m en u m e r i c a l e x a m p l e s k e yw o r d s ：t r e f f t zf i n i t ee l e m e n t ，f o r c e dv i b r a t i o n ，t h i np l a t e 独创性声明 本人声明所呈交的学位论文是本人在导师指导下进行的研究工作和取得的 研究成果，除了文中特别加以标注和致谢之处外，论文中不包含其他人已经发表 或撰写过的研究成果，也不包含为获得：叁鲞盘堂或其他教育机构的学位或证 书而使用过的材料。与我一同工作的同志对本研究所做的任何贡献均已在论文中 作了明确的说明并表示了谢意。 学位论文作者签名：赦j 孑讦 签字日期：j o 。 年占月，) 日 学位论文版权使用授权书 本学位论文作者完全了解苤壅盘堂有关保留、使用学位论文的规是、 特授权苤盗盘堂可以将学位论文的全部或部分内容编入有关数据库进行径 索，并采用影印、缩印或扫描等复制手段保存、汇编以供查阅和借阅。同意学校 向国家有关部门或机构送交论文的复印件和磁盘。 ( 保密的学位论文在解密后适用本授权说明) 学位论文作者签名： 靛日开 导师签名： 签字同期：) 口叩年6 月f 2 日 谚歹 签字日期：z 呻年6 月c 2i 天津大学硕十学何论文 第一章绪论 1 1 引言 第一章绪论 由于现代科学技术发展，许多机械与设施必须考虑在动载荷条件下的设计计 算，而板的动力问题是现代许多工程部件设计和研究的关键【l 】。有关动力响应分 析的方法，数十年来，许多学者进行了不懈的研究与探索。迄今为止，已有大量 的理论成果。除振型叠加法之外，应用比较广泛的是建立在有限元基础上的直接 积分法，它是对所求问题在空间域进行离散，在时间域采用差分法，逐步递推求 解。直接积分的方法很多，各种方法在数学上的收敛性和稳定性不同( 2 - 3 1 。其中， 平均加速度法是无条件稳定的，即时间步长与系统的周期之比都不会由于该法本 身而造成解的无限增长。 作为计算力学中的一种数值方法，杂交t r e f l t z 有限元模型的特点是在单元内 部和单元边界分别假设两组位移插值函数【4 】，单元内部插值函数的选取必须遵循 事先满足控制方程的原则；单元边界上的插值函数( 常称为网线函数) 满足单元间 的位移协调条件。最后用变分法或其它方法将这两组独立的位移函数联系起来， 使最终求解的方程仍以结点位移参数为未知数的线性方程。 由于t r e t t t z 有限元法继承了传统有限元法和传统边界元法的长处，该单元 具有以下优点【5 】：1 推导公式中只含有对单元的边界积分，这样就可以生成任意 多边形单元甚至曲边单元。因此，t r e f f t z 有限元法可被认为是一种特殊形式的边 界型求解方法，其刚度矩阵对称、计算简单，而不像传统边界元那样要进行复杂 的计算( 如计算边界奇异积分的复杂积分规则，处理非齐次方程的特殊积分等) 嗡7 1 。2 t r e f f t z 有限元在处理各种几何、荷载奇异性或局部效应( 如角点、孔洞、 裂纹、夹杂、集中载荷等) i 口- i 题有较高的效率【6 7 】。因为只需将插值函数取成既 满足控制方程，又满足给定的边界条件即可，而这种特殊函数通常可以从文献中 找到。这样，在遇到角点、孔洞、裂纹、夹杂以及由各种载荷引起的应力集中或 各种奇异性时，不需局部网格细分，只需在插值函数上作调整，就能达到所希望 的精度。 天津大学硕十学能论文第一章绪沦 由于i ：述种种优点，到目前为止，t r e f f t z 有限元法已经得到较多学者的关注， 并取得了定的成果。然而，目前对薄板强迫振动问题的研究还很不足。鉴于此， 本文从t r e f f t z 有限元基本理论出发，结合无条件稳定的平均加速度法，建立薄 板强迫振动问题的t r e f f t z 有限元模型，分析了薄板受迫振动问题。 1 2t r e m z 有限元的发展现状 t r e f f t z 方法是g 了t r e f f t z 【8 】于1 9 2 6 年首先提出的，t r e f p t z 将这一方法应用于解决 l a p l a c e 方程i 均d i r i c h l e t t h 7 题。其后，t f e f 拖法由b i m a n 【9 - l ，r a f a r s o n 【1 2 ，1 3 1 等人的研 究中得到了推广。然而直到2 0 世纪7 0 年代才逐渐引起了学术界的重视。在经典的 t r e f f t z 法中，t r e f f t z 函数中的待定系数通过边界条件的某种加权残数方式确定， 属于间接法【1 4 】。19 8 9 年，c h e u n g ，j i n , z i e n k i e w i c z 等人【1 5 】提出tt r e f f t z 直接法。该 法是将t - 完备解函数作为权函数，对控制方程的加权残余表达式应用两次格林公 式构造出边界积分方程，以节点物理量为未知数直接求解边界节点的物理量值。 s t e i n t l 6 1 、r u o 一1 7 1 、z i e n k i e w i c z 1 8 1 等人通过有限元与边界型解法的联合使用，对 t r e f f t z 方法做了早期的数值研究。随着研究的深入，t r e f f t z 法与有限元相结合的 数值方法逐渐发展起来。 t r e f f t z 有限元是在七十年代末才提出来的一种有限元模型【l9 1 。作为计算力 学中的一种数值方法，h y b r i d - t r e f f t zo a t ) 有限元法己开始引起人们的重视并在 工程中得到越来越广泛的应用。与传统有限元模型不同的是，t r e f f t z 有限元模型 利用辅助网线位移场或面力场，在一种杂交意义上将单元域内位移场关联起来。 单元域内位移场精确满足控制微分方程，它可表达为微分方程的特解加上适当截 断的t r e f f t z 完备解系( 或完备解系) 与待定参数乘积的形式。利用定义在每个 单元边界上的独立的网线位移场，单元间的连续性就在一种近似意义上得到满 足。在单元一级上消去内部待定参数后即可得到标准的力位移关系式( 即单元 刚度方程) 。t - 完备解系的数学理论主要是由h e 玎c r a 及其合作者完成的【2 0 - 2 3 1 ，并 将其命名为c 完备解系。为纪念这种非奇异解的创立者t r e f f t z 、z i e n k i e w i c z 2 4 建议将其更名为t 完备解系。z i e l i n s k ia p 和z i e n k i w i c zo c 1 2 5 j 将t r e f f t z 完备函 数引入有限元法，求解广义协调方程。j i r o u s e k 和g u e x 2 6 1 提出杂交t r e f f t z 单元， 建立了杂交t r e f f t z 有限元法( h y b r i dt r e f f t zf e m ) 。秦庆华 4 1 的专著对于t r e f f t z 有 天津大学硕十学位沦文第一章绪论 限元有详尽的阐述，可以说t r e f f t z 有限元已经被广泛地应用于平面弹性问题 f 2 7 2 8 1 、薄板弯曲【2 9 - 3 0 1 、中厚板【3 1 】、厚板【3 2 1 、p o i s s o n 方程 3 3 , 3 4 、热传导f 3 5 】、压电 问题3 s 、自由面渗流问题【3 9 删、几何非线性问题【4 1 删、材料非线性问题【4 5 - 4 7 等各个力学领域。为了推动t r e f f t 2 方法在计算力学领域的发展和应用，自1 9 9 6 年 起每四年定期召开围际性专题讨论会，目前已经成功地召开了四次。其中第一次 于1 9 9 6 年在波兰的科拉科召开，同时纪念t r e f f t z 发表论文“g e g e n s t u c kz u m r i t z s c h e nv e r f a h r e n ”7 0 周年；为庆祝j i r o u s e k 教授潜心致力于t r e f f t z 有限元公式的 推导及其软件s a f e 的开发，1 9 9 9 年在葡萄牙的辛特拉召开第二次国际研讨会； 紧接着第三次于2 0 0 2 年在英国埃塞特召开：第四次于2 0 0 5 年在斯洛伐克z i l i n a 召 开。 t r e f f t z 单元在平面弹性问题中的应用开始于j i r o u s e k 和t e o d o r e s c u t 4 8 1 的工 作。根据网线位移场是定义在整个单元边界上及单元交边界上，他们导得了平面 弹性单元的两个变分公式。随后，j i r o u s e k 和v e n k a t e s h t 4 9 l 从m u s k h e l i s h v i l i 复变 函数法出发，发展了p - 型t r e t i t z 平面弹性单元。此外，还出现了许多其它形式 的t r e f f t z 平面弹性单元，如k o m p i s 等人【5 0 1 、s l a d e k 等人【5 1 】以及秦庆华嗍的工 作。 1 3 弹性薄板的动力响应分析 由动载荷引起的结构振动称作受迫振动或动力响应分析。所谓求动力响应就 是在初始条件下求解运动方程，确定各个时刻的状态矢量。进而利用应力和节点 位移的关系，还可以得到物体内部各点的应力随时间的变化规律。对于动力响应 的分析，可采用振型叠加法、直接积分法( 或称逐步积分法) 等方法【2 】。其中，振 型叠加法的基本思想是通过坐标变换，将一个多自由度体系的n 个耦合运动方 程，分解为n 个非耦合运动方程，问题的解为n 个非耦合运动方程解的线性组合。 由于用振型叠加法求解结构的动力响应时，须先求出结构的固有频率和振型，使 得计算步骤和程序的编制比较复杂。而直接积分法与振型叠加法不同，无需先进 行振型分析，也不对运动方程进行基底变换而是直接对运动方程进行逐步数值积 分，因而适用范围较广。其基本思想是：1 对时间离散时，不要求任何时刻都满 足运动方程，而仅要求在离散点上满足运动方程。2 在时间间隔f 内位移、速度 和加速度的变化规律及其间关系是假设的，采用不同假设得到不同的直接积分 天津大学硕十学位论文 第一章绪论 法。 直接积分法的计算过程是：假设f = 0 时刻的状态向量：位移、速度和加速度 是已知的：如将时间求解域0 t t 进行离散，即可由已知的t = 0 时刻的状态向 量计算忙o + a t 时刻的状态向量，进而计算t = f + f 时刻的状态向量，直至f = t 时刻终止，这样便可得到动力相应全过程。 j a t e i x e i r ad ef r e i t a s 于1 9 9 6 年在处理频域中的弹性动力学问题的分析过 程中，曾首次将t r e f f t z 有限元概念应用于动力学问题中【5 3 】。然而，在实际中， 此法是难于实现的。其主要原因为在计算单元刚度矩阵时，修正的运动方程中仍 存在未定频率参数d ，这样将会导致在构造t r e f f t z 函数时依然存在此未知参数， 从而使杂交矩阵h 的求逆受到阻碍。鉴于此，本文将在第三章中详细介绍，基 于在离散的时间t ；( i - 0 ，1 ；2 ) 下运用直接积分法中的平均加速度法与t r e t f t z 有限元基本理论结合的方式来修正运动方程，从而避免上述问题的出现。 1 4 本文的主要工作 本文首先对平面弹性力学问题，利用m u s l d a e l i s h v f l i 复变函数构造了t r e f f t z 函数( 单元域内) ，建立平面t r e f f t z 有限元，并通过几个典型的数值算例证明了 该单元的正确性，为下面工作打下基础。 然后，建立薄板强迫振动问题的t r e i t t z 有限元。在该单元中，本文采用将 时间域进行离散，运用无条件稳定的平均加速度法对薄板横向振动微分方程进行 修正，同时利用第一类k 阶交型b e s s e l 函数和第一类k 阶h a n k e l 函数构造单元 域内位移场插值函数中的特解及齐次解，并根据t r e f f t z 有限元基本方法建立相 关矩阵，形成有限元模型，并通过几个的算例说明了该单元的精度及效率。 最后，本文对t r e f f t z 有限元的建立进行了讨论和说明。 天津大学硕十学位论文 第：章t r e f f t z 有限元法 2 1 位移插值 第二章t r e f f t z 有限元法 与传统有限元法不同，t r e f f t z 有限元法分别在单元域内和边界上独立地建立 两个位移场1 4 , 1 9 。在一给定单元p 的区域q 。内位移场为： 而 u 。= f i 。+ n j c j = f i 。+ n 。c 。 ( 在域q 。内) ( 2 1 ) j = l 式中，n 。为位移向量，c 。为待定参数向量，而为向量c 。的维数。讧。和n 。分别为 问题控制方程 m ( u ) 一6 = 0 ( 在域q 内)( 2 2 ) 的特解和齐次解，即： m 倾) = i ，m ( n ，) = 0 ( 在域q 。内)( 2 3 ) 式中，m 和i 分别为控制方程的线性算子和外力矢量。 另一方面，在单元边界上，第二组独立位移，常称为网线位移5 4 1 ，可设为 n 。= n 。d 。( 在r 。上) ( 2 - 4 ) 式中，矗为单元边晃l 上的位移分布，r 为插值矩阵，它的构造遵循相邻单元间 位移协调的原则，d 为通常的节点自由度矢量，上标“”的使用是为了与域内位 移区分开来。下节将讨论如何通过变分法在某种弱形式下【4 1 将d 和c 关联起来， 最终建立以节点位移参数为基本变量的单元刚度方程。 2 2 变分原理的构造 用于推导t r e f r z 有限元公式的泛函可对传统的余能泛函加以修正而得。为此， 对于线弹性问题，考虑如下的平衡方程、本构方程、几何方程以及给定位移边界 条件和给定外力边界条件等基本方程： 天津大学颀十学待论文第：章t r e f f t z 有限元泫 平衡方程 l 仃十i = 0( 在q 内) ( 2 5 ) 本构方程( 应力应变关系) o = d ( 在q 内) ( 2 6 ) 几何方程( 应变位移关系) = l t u( 在q 内) ( 2 - 7 ) 给定位移边界条件 u = 面 ( 在l 上)( 2 - 8 ) 给定外力边界条件 t = a o = i ( 在c 上) ( 2 - 9 ) 式中q 代表弹性体所占据的区域，l 、e 分别代表给定位移边界和给定外力边 界，并有f = lur r ，r 代表弹性体全部边界，o ，1 1 1 和t 分别表示应力、应 变、位移和边界力，1 6 f ，i 表示给定的位移和边界力，百表示单位体力，l 、d 和 a 分别是微分算子、弹性矩阵和方向余弦矩阵【5 5 】。对于平面应力问题，有 l = t ，o = 专 11 ， v1 oo o o 1 1 ， 2 h ，0 a5 i ( ；刀y 对于平面应变问题，e 、v 分别用e ( 1 - - v 2 ) 、v ( 1 一v ) 替代。这里e 、y 分别 表示弹性模量和泊松比，、捍y 表示边界上任一点处的外法线方向余弦。 对于给定了单元边界l 上的表面力t 的典型单元e 则有余能泛函 群= 脖；c 。d f l 一弦t 。d r ( 2 - 1 1 ) q 。 相应的最小余能原理铘。= o 的使用条件是：在单元域g 内应力满足平衡 方程( 2 5 ) 和在单元边界l 上满足给定外力条件( 2 - 9 ) ，其中c = d ，万代表变分 符号。 进一步地，在t r e f f t z 有限元法中方程( 2 5 ) 一( 2 9 ) 还需要增加下列单元交界面 间的连续性条件【2 7 1 ，即 位移相容条件： u 。= u 曲( 在l = 吒r 、l ) 上)( 2 1 2 ) a一砂a一西 o a 一砂 a 一良 o 天津大学硕十学位论文第：章t r e f f t z 有限元法 表面力平衡条件： t 。= t 曲 ( 在r f ，= ( i 乙厂、c 。) l ) ( 2 1 3 ) 其中l ，l 分别为相邻单元各自的边界，而r e 为其交边界，l = ln l ， l = ln r , ，并且有r ，= r 。u r 。，u r 。，。 与系统总余能n 。= e n 。对应的最小余能原理的使用条件除了( 2 - 5 ) 和( 2 9 ) 外，还要求相邻单元间满趸面力平衡条件( 2 1 3 ) 。为了放松单元边界c ，上的平衡 约束t 。= t ，定义单元边界位移矗。为拉格朗日乘子，则修正余能泛函n 。 2 6 1 由驻值原理，可得 。= n 。一( t 。一己) t 吨d r j t 汽d f ( 2 - 1 4 ) t t ti d o t l 。= 万n 删= o 0 t i 。= o j 6l 肇州。a | d q 一8s t , g d r 一6 沁i 一弧j d r - 8 i t f i ，d r q 。 lr ，。r 7 = 啦而，) d q 一乎瓦) d r 一弘一现 d r 一裂砘羁) d r = 。 ( 2 1 5 ) ( 2 - 1 6 ) 由对称佳和夫糸式 勺= ，= 去g u + g j , i ) = ， 式( 2 1 6 ) 第一项中的被积函数可写成 毛i ! ：自o 。+ c 。p 圆口0 = c 粥口声o 。= 扛i 芦。t 将上式代入( 2 1 6 ) ，再应用平衡方程，+ 瓦= o ，则 f f z d i 。p ia q = s u , a t , d r - l p ；协d f 2 = i u f i t ，d f p 1 将式( 2 - 1 7 ) 代k ( 2 - 1 6 ) 中，得到 m 所。= j 甜，8 t ，d r e 以d r 一龀一磅观+ 西，以 d r 一( f f 观+ 霉以) d r 1 w 蚺dr - 、i e i - t l = j ( g i 一巧溉d r 一心一t , - 7 ) a u ；d f + f 融，一厅i ) a t j t i 万u i ) f = 0 r -r ，。r r ( 2 - 1 8 ) 从而可以得到式( 2 8 ) 、( 2 9 ) 、( 2 - 1 2 ) 、( 2 1 3 ) 天津大学硕士学位论文 第：章t r e f f t z 有限元法 上面的证明过程表明，余能泛函h 。，的驻值条件满足边界条件和单元问连续 性条件，因此可以用来推导t r e m z 有限元公式。 2 3 单元刚度方程的形成【5 5 】 利用域内位移场( 2 - 1 ) 和关系式。= l u 。和o 。= d l u 。，很容易得到相应的应 变场t 和应力场o 。： 式中 式中 ( 2 - 2 0 ) 而 6 。= d 。= 云。+ t j c j = 云。+ t c 。；( 2 - 2 1 ) j = l t = l 豇。，b 。= l n 。( 2 - 2 2 ) 舌，= d l f i 。，t = d l n 。( 2 - 2 3 ) 一旦求得应力的特解和齐次解，单元边界力t 。= a o 。就可以写成 t 。= t + z q j 巳= t + q 。c 。 ( 2 2 4 ) = j 乞= a o 。，q 。= a t( 2 2 5 ) 下面利用修正余能泛函兀一来推导t r e f f t z 单元公式。为此，考虑式( 2 5 ) 、 ( 2 6 ) 、( 2 - 7 ) 并规定在l 上：矗。一u e ，以此可把兀。修改为 1 1 m ”= = 三豇( 】二。i l e ) 1d l t u , d f 2 - ，： e t d rd i c t r 乏，t 面cd i - ，e t 面e d 1 1 c 2 - 2 6 ， 。q ， r r dr d 对上式右边第一项应用高斯散度定理，得到 ( l t u 。) 1d l t u ，d r 2 = t 沁d r 一u j l d l t u ，d q ( 2 - 2 7 ) f cb+ p 一 = 勺 b 扁州 + c - = f 天津大学硕十学何论文第j ：章t r e f f t z 有限元法 则泛函( 2 2 6 ) 可变成如下形式 n 。= 丢肛。d r + 告腮u 。d q f t 疆d r _ ( t 。一i 。) 1 矗，d r - 胯。( 1 f ( 2 - 2 8 ) l。q r kb 再将( 2 1 ) 、( 2 - 4 ) 和( 2 2 4 ) ；j l , , k 上面的泛函式中，有 。 n 。= c j h ，c 。一c j g 。d 。+ c j 氍+ d ；h ，+ 不包含c ，或d 。的项 ( 2 2 9 ) 利用泛函驻值条件硼o c j = 0 ( 这实际上是在单元一级上使用最小余能原 理) 得 c 。= h ：1 ( g 。d ，一h 。) ( 2 3 0 ) 代回至( 2 2 9 ) 后，再利用泛函驻值条件饥懈a d j = 0 ，便得到下面的单元刚度方 程： 式中 k 。d 。= e k 。= g t 。h - 。1 g 。 p = g j h 。h 。- g 。 ( 2 - 3 1 ) ( 2 3 2 ) ( 2 - 3 3 ) 分别称为单元刚度矩阵和单元等效节点载荷列阵。这里，辅助矩阵h 。、g 。、h 。 和g 。的表达式分别为 d f = n ：t q ，d f ( 2 3 4 ) 。= l q ：n 。= l n 。q ， ( 2 。) r tl g 。= 眨丸d r ( 2 - 3 5 ) r f h 。= 吾( n hq 拽) n 吾n j 瓦m 。l。n ，f 一 ， g 。= l n ：t ed f i n - t 。- - i 。d f ( 2 - 3 6 ) ( 2 - 3 7 ) 其中对称正定矩阵h 。称为单元形变矩阵，由此知单元刚度矩阵k 。也为对称正定 矩阵。当体力忽略不计时( 即e = 0 一面，= t = 0 ) ，单元应变能u 。就等于表面 力t 。= q 。c ，作用在相应位移面，= 衬。d 。上的外力5 9 1 ，即： 天津大学硕十学位论文 第二二章t r e f f t z 有限元泫 其中 u r2 抄叫r 毛廿a r - 吉c 陬t = 扣k a弘3 8 ， h 。= f q ；n 。d r = j n j q 。d f = i 1j ( n 。r q 。+ q j n 。) d r ( 2 - 3 9 ) l 。l 2 4 二维弹性问题 2 4 1 单元域内场位移插值 将某一特定单元e 域内的非协调位移场可写为： u 。= ：) = 芝) + 喜n ，c ，= 讧。+ n 。c 。 c 2 4 。， 式中c 。为待定参数矢量。 若体力- f 为常数，那么位移特解豇。可取成 睁阱半圈 若瓦为b ，y ) 的一般函数，则特解面。可通过如下积分来表示吲 讧。= 耋：) = ：妻乏：塞乏 d q 式中甜；为平面应力问题的格林函数5 川： “；( ) = 石l t - v 【_ 一( 3 1 ，) 8 vi n r m + ( 1 + v ) 嘞，嘞。一 而”；( 嘞) 代表由在q 点沿方向的单位力引起在场点尸沿i 方向的位移，距离 函数嘞定义为： 嘞= 厄i 再瓦而( 2 删 ” 习 筇 4 4 4 q q p 天津大学硕十学何论文 第：章t r e f f t z 有限元泫 源点q 、场点p 和距离函数嘞的定义也可参见图2 - 1 。 场点p r e o 源点q o y 图2 - 1 距离函数示意图 其中x ，y 表示笛卡尔局部坐标。 平面弹性问题t r e 妇娩函数n 可借助于m u s 姓e l i s h v i l i 复变函数法来导得【2 6 】。 其措施是将所有场变量用两个待定的解析函数痧g ) 和沙( z ) 来表示【5 7 】： e g ，+ f 材y ) = ( 3 一v ) 痧g ) 一( 1 + v ) k 瓦) + 矿g ) j ( 2 - 4 5 ) 吒+ = ( z ) + g ) 一z 痧g ) 一歹( z ) ( 2 - 4 6 a ) 一= ( z ) + 声( z ) + z ( z ) + 歹g ) ( 2 - 4 6 b ； 式中z = x + y ，i = 行，矽和杪是在所考虑区域内的解析函数，变量上方的圆 点和横杠分别表示对z 求导数和对z 取共轭。如果依次将函数痧( z ) 和| ；f ，( z ) 取成： ( z ) = 拓，y ( z ) = o 痧( z ) = ，g ) = 0 ( z ) = o ， 少g ) = z ( z ) - - o ， f ，g ) = z = 1 ，2 ，) ( 2 - 4 则齐次解n ( = 1 ，2 ，) 的完备解系可由如下序列的实部和虚部构成： 铡，= 【r h e z z 。 ，z l 。= ( 3 - v ) 山( 1 川之 ，= 髓：) ，2 ：。= 3 - v 灿+ v ) 矿1 ( 2 4 8 j ( 2 - 4 9 ： 天津大学硕+ 学待论文 第：章t r e f f t z 有限元法 磊= ( 1 + y ) ( 2 5 0 ) z 0 = 一( 1 + 1 ，) z ( 2 5 1 ) 对于每个n ，相应的齐次应力解r i ：，可由( 2 - 2 3 ) 第二式给出 恐j r e 恐+ r e i m l+：=盖re耋si3。，5。=旅主一l ，i ：，+ z = 一t ，= 肥卜1 【i n l j = 圈 墨i = 影q 迸一步地，由公式( 2 2 5 ) 第二式便可由e 求得边界力齐次解q ，如下： ( 2 - 5 2 ) ( 2 - 5 3 ) ( 2 5 4 ) ( 2 5 5 ) q ，= n j ( r e r u + r e s l k ) + n y i m s l , 墨。，墨。的含义同上 ( 2 5 6 ) q 川= 1 b n ( ( r r e e r 恐2 , 。- r e s 2 k 。) ) - 一n y i r e s 2 , + r e s zn x i r n s 2 。) r 。，是。的含义同上( 2 5 7 ) u 川2 1 ，2 ，( r e 恐女 i ) 一lf 屹屯刈u 函义h 工 u 。 f q 心= l 玛。，的含义同上 ( 2 5 8 ) r 。，墨。的含义同上( 2 5 9 ) 需要注意的是，当尼= 0 时，式( 2 5 0 ) 、( 2 - 5 1 ) 沿着坐标轴将产生常量平动，即： 1 2 乙磊 m 眦眦 呱 位贴h =，i弋 = 一 3 删 掣掣肪 h o扣 汝玎 = = 七 哺r 墨 、ll：-l， s s r r 根氓帆 。 。 m置墨l r r = 弓 1 以 存彬 d d 口肛 甜 r，【 e e r r ，r，【 = + i 吸崛山以蒜 尽b 、l，j 豇鼬 h h 以 一 一 女 & 鼬耻 珞凡 一 ，j、【 = m q 天漳大学硕十学位论文 第一章t r e f f t z 有限元法 e u ，= 0 ，e u ，= 1 + ，和e u ，= 一( 1 + v ) ，e u ，= 0 。当k = l 时，式( 2 - 4 8 ) f ：生常量转 动为w = 恤，砂一彻，缸) = _ 4 ，其中e u ，= 4 y ，e u ，= 4 x 。从k = 1 和式( 2 - 4 9 ) 开始，随着多项式阶数k 的增加，序歹u ( 2 - 4 8 ) - ( 2 5 1 ) 能生成独立的位移模式( 当 k = l 时有三项独立模式；k 1 时有四项独立模式) 。因此，在从序列( 2 4 8 ) - ( 2 5 1 ) 中选用元素作为域内插值函数n ，时，为避免造成h ，阵奇异，须舍弃对应于刚体 位移模式( 即零应变模式) 的n ，项，相应的t ，和q ，项也应舍弃。 此外，针对序y u ( 2 4 8 ) 一( 2 - 5 1 ) 、( 2 5 2 ) 一( 2 5 5 ) 以及( 2 5 6 ) ( 2 5 9 ) ，我们作下面几 点说明： 1 在t r e t t t z 有限元公式中应用的杂交技术表明，生成的齐次位移解n ；将产 生非零的应力t ，这是保证前面定义的单元形变矩阵h ，或单元刚度矩阵k 。非奇 异的必要条件。为此，需舍弃齐次解矩阵n 。中的三项刚体运动模式。 2 t r e f r z 项数而对单元性能影响较大。一般来说，疡过大，势必增加h 。的 阶数，从而降低计算效率；反之，疡过小会造成虚能量模态1 4 1 。n ，n ，n 南的 最小数目建议取为 ( 2 6 0 ) 这也是保证单元刚度矩阵k 。满秩的必要( 但有时非充分) 条件，其中n d o f 是 单元节点自由度数。为了保证充分性，通常采用适当增加t r e f f - t z 项数的办法使 得确 虹。对于具体的单元类型，廊的最佳值需通过数值实验求得。 3 为了维持单元在坐标轴旋转情况下的几何不变性，采用一个简单的“截断 法则”，它以( 2 - 4 9 ) 或( 2 5 1 ) 为截断点。对于k = 1 只产生三个独立函数，因此历通 常是一个奇数。 4 。为了防止在对h ，求逆时数值溢出或出现病态，通常采用t r e f f t z 无量纲坐 标善和7 7 来代替笛卡尔局部坐标x 和y ： x = n 童，y = n q - 6 1 ) 式中 铲三n 窆i = l 拇可 ( 2 6 2 ) x = x 一五= x 一专莓置 。2 6 3 ， o 夕2y 一艺= y 一言喜y ， 、。 天津大学硕士学位论文第：章t r e f f t z 宵限元泫 这里刀为单元节点数，嘭为单元各节点到其形心的平均距离，称为单元特征长度， 以、艺为单元形心的笛卡尔总体坐标。 2 4 2 单元网线位移插值 如图2 - 2 所示，t r e f f t z p 一扩展单元是这样一种单元，它在单元角点有固定的 自由度，即通常的有限元节点参数。而在单元边中节点可按需要增加适量的虚 拟自由度数。一般地，对于平面弹性三角形或四边形单元，在任一边上通过两 端点的插值只能得到协调的线性变化位移场。若要得到高阶变化位移场，可在单 元边中节点引入额外的虚拟自由度。 域内位移场 而 u 。= 藉。+ n ，c j = f i 。+ n 。c 。娩) j = t y 2 - = 一1 f = 0 f = + 1 ) m n b n c l n c 2 卜 一 r l 7 一、 e t c 固定自由度( 圾，l d 虚拟自由度( a l ，b l ，a 2 ，b 2 ) 图2 - 2 普通t r e f f t z p 寸展元 圭( 1 一f ) ( 1 + f ) ( 1 一f 2 ) f ( 1 一f 2 ) 以图2 - 2 中a c - b 边为例，平面弹性问题p 扩展元的网线位移可设为【2 6 】： “u x 一- - - i 每h + 畏毒酝+ 皂铲衬c _ a c j ? 覆y ： 蓐b + 苻毒和+ m 护衬c b c i ( 2 6 4 ) 其中占是为了保证相邻单元间高阶位移模式的唯一定义而引入的一个系数， 它的数值取决于单元边a c