考虑截面变形的自然弯扭实心截面梁几何非线性理论

考虑截面变形的自然弯扭实心截面梁几何非线性理论刘洋，单建东南大学土木工程学院，江苏南京 (210096)摘 要：截面变形理论包含拉、弯、剪、扭荷载引起截平面内外的变形。

本文采用它分析自然弯扭实心梁截面的变形。

引用Jaumann 应力应变，两Euler 角表示的有限转动矩阵和虚局部转角，通过Hamilton 原理推导自然弯扭梁的动力平衡方程和力边界条件。

关键词：截面变形；有限转动；Jaumann 应力应变；Hamilton 原理；虚局部转角1. 引言以往对自然弯扭梁的研究往往忽略截面翘曲的影响，或者只考虑扭转引起的翘曲[4][5][6]。

后来发展为用二维截面有限元法来分析截面翘曲变化，但其精度取决于单元的个数和位移模式的精度[3][8]。

E.Petrov 和M.Geradin 提出截面变形理论分析直梁的大变形[2] ,无需离散截面，形式简单。

本文采用此理论来分析自然弯扭梁的截面变形。

由于第二Kirchhoff 应力和Green-Lagrange 应变自身定义的特点和便于计算的优点，在多数关于非线性分析的文献中，它们成为最为常用的应力应变形式。

但Jaumann 应力和应变除有相同的优点外，还更符合线弹性本构关系，因此本文采用了这种应力应变。

Hodges 和Dowell 于1974提出两Euler 角表达有限转动，其优势在于能简洁表达各变量之间的关系，便于公式推导。

但两角度之间不是独立而是相继的关系，因此不便于作为位移变量，于是Perngjin F.Pai ，Ali H.Nayfeh 引入局部虚转角的概念[3]，作为位移变量。

本文沿用两Euler 角表达和局部虚转角，借助Hamilton 原理推导出的动力平衡方程和力边界条件，较为全面的反映了各种因素的耦合。

2. 截面翘曲分析利用弹性力学中三维梁端部受轴向拉（压）力、竖向集中力、弯矩、扭矩的精确位移解，组建截面变形位移。

位移的详细推导过程见[1]、[2]。

梁截面变形位移的矩阵形式为1[]⋅w =B λ （1a ）[][][]12311212233121122332122220001[]000()21()2TTTw w w e k k k e k k k vyvyzv y z vz v y z vyzγγαρργγρρρχχϕ=−−−=−−−⎡⎤⎢⎥⎢⎥⎢⎥−−−⎢⎥⎢⎥⎢⎥−−⎣⎦w =B = λλ (1b)其中v 为泊松比，12,,χχϕ为两个方向的截面剪切变形函数和截面扭转变形函数，不沿轴坐标变化。

e 为轴向变形，1γ，2γ为剪切变形，1k α−为变形后翘曲率与初始扭率的差值，22k ρ−，33k ρ−为两个方向上变形后弯曲率与初始弯曲率的差值。

上式中用了α而没有用变形后扭率1ρ是因为考虑到横向剪切变形1γ、2γ引起的剪力流对截面转动产生影响，使横截面扭转角不再等于翘曲角。

3. 自然弯扭梁的初始曲率求解首先定义一种表达图1 初始弯扭梁的几何关系3212331210[(,,)]00Ta a a a a a a a a −⎡⎤⎢⎥=−⎢⎥⎢⎥−⎣⎦Ω （2）变形前梁上任意一点的位置矢量和正交基为23123={}={,,}y z ++r r e e e e e e （3）初始曲率定义为123{}=[(,,)]{}Tk k k ′e e Ω （4）1k 为初始扭率，23,k k 为两个方向的初始弯曲率。

仿照变形前位置矢量的形式，表达出变形后梁上任意一点的位置矢量（不考虑截面变形） 12323123 {}{,,}Tu v w y z =+++++R r e e e e e e e e e = （5） 变形后的曲率定义为123{}[(,,)]{}Tρρρ′=e e Ω （6）1ρ为变形后扭转率，23,ρρ为两个方向变形后的弯曲率。

4. 有限转动和变形后曲率求解有限转动是分析梁板壳大转动的有力工具，它建立起了变形前后曲线坐标之间的关系，也是计算变形后曲率的关键，有关系表达式 {}[]{=}e R e ，其中[]R 为有限转动矩阵。

与传统的表达方式不同，本文将采用两Euler 角表示法[3]表达有限转动矩阵，如图２。

第一个欧拉角α由三个位移{,,}u v w 所引起，叫做弯曲转动角；第二个欧拉角φ表示绕1e 的转动。

因此，转动矩阵表式如下111213212223313233100[]0cos sin [()]0sin cos R R R R R R R R R φφαφφ⎡⎤⎡⎤⎢⎥⎢⎥==⎢⎥⎢⎥⎢⎥⎢⎥−⎣⎦⎣⎦R B (7) 接下来求()αB ，因为1111122133R R R =++e e e e (8a) 其中3231211112131111u vk wk v uk wk w uk vk R R R e e e′′′+−++−−+===+++，，1e =(8b)图2 两Euler 角有限旋转 图3 弯曲转动角α与位移的关系通过确定旋转轴n1111223311n n n ×==++×e e n e e e e e (9a)其中1230,n n n ===(9b)并利用向量形式的旋转公式cos (1cos )()sin ααα+−⋅×v =v v n n +n v ，得())())111213212111311121311213121311111211[()]+1+1+1++1+R R R R R R R R R R R R R R R R R α⎡⎤⎢⎥=−−⎢⎥⎢⎥−−⎣⎦B (10) 详细的推导过程见[3]。

求得变形后的曲率矩阵和曲率为 123123[][(,,)][][][(,,)][]TT T T k k k s ρρρ∂=⋅+⋅⋅∂R R R R ΩΩ (11) ()()()33''211323123132113'1321231,,i i i i i i i i i i i i i i i R R R k R R R k s s R R R k s ρρρ===∂∂≡⋅=+≡−⋅=−+∂∂∂≡⋅=+∂∑∑∑e e e e ee (12)5. Jaumann 应力和Jaumann 应变以往文献，如[4]、[5]、[6]分析非线性问题往往采用第二Kirchhoff 应力和Green-Lagrange 应变，这种应变虽然满足客观性、共轭性、对称性，但它属于能量型应变，同时所在的分析坐标为随动坐标，不相互正交。

而Jaumann 应变是几何型应变，同时所在的分析坐标为Lagrange 坐标，相互正交。

因此，Jaumann 应力和Jaumann 应变更加符合由工程应力、应变测出物质常数的本构关系。

详细分析见[7]。

本文将放弃[3]，[7]中的局部位移、局部应变的定义，直接从Jaumann 应变的定义出发，求解应变。

利用右伸长张量U ，给出Jaumann 应变的定义−J =U I 。

通过下式计算Jaumann 应变矩阵为[][][][]2T+=−F F J I （13）[]F 为变形梯度矩阵，可通过 ()()kij i k j F =⋅⋅e G g e 计算。

其中变形前的逆变基和随体协变坐标基为12323==+=+g ,g e ,g e (14a)321yk zk =−+ (14b)122331211112211333223111233122111223312221233[1()()][()] [()]()[()] [()]()e z k y k z y w z k w sy k w w w w z y w w w z y w w wy y y yρργγρρρργγρρργγργ∂′′′′==++−−−++++−−+∂′+−++−+++−+−++∂∂∂∂==+++∂∂∂∂∂=RG e e e e e e R G e e e G 3121123()w w w z z z zγ∂∂∂=+++∂∂∂∂Re e e (14c)其中1232321112233()u v w y z y z w w w γγ=++++++++++R r e e e e e e e e ，即考虑截面变形后，梁上任意一点的位置矢量。

式中21y z γγ、为从截面翘曲中分离出的常规剪切项，1121w w y z γγ=−−。

将求得Jaumann 应变各分量，用voigt 标记法标记，然后把应变矩阵写成列阵形式，得 {}111213222333{}T J J J J J J =J (15a)其中1122331213223112121122112331113131131211122()()2()()2()()J e z k y k z y w w w wJ J z k w w z y w y wJ J y k w w w z y z ρργγρρργγγρρργργγρ′′′=+−−−++++−∂′==−−++++++−∂∂′==−++++−++∂332222232333, 2, w w ww J J J J y y z z∂∂∂∂===+=∂∂∂∂ （15b ）Jaumann 应变也表示成列阵的形式{}111213222333{S S S S S S }T =S ，本构方程为{}{}[]=⋅S D J 。

6. Hamilton 原理Hamilton 原理表达为12()d 0t H t T W δΠδδ∏δτ=−+=∫。

T 为系统动能，∏为总应变能，W 为外力。

6.1动能项变分分析系统总位移为123232*********()u v w y z y z y z w w w γγ++−−+++++D =e +e +e e e e e e e e （16）为了便于计算，引入虚刚体转动角的变分来表示局部坐标向量的变分，于是有123{}[(,,)]{}T δδθδθδθ=e e Ω ，系统位移变分整理为112112132()()(){}{}[(,,)]{}000000T T T T u y z k v w y z w δδθϕχχδαδδδθδγδδθδγ++−⎧⎫⎧⎫⎡⎤⎧⎫⎪⎪⎪⎪⎪⎪⎢⎥=⋅+⋅⋅+⋅⋅⎨⎬⎨⎬⎨⎬⎢⎥⎪⎪⎪⎪⎪⎪⎢⎥⎩⎭⎩⎭⎣⎦⎩⎭D e e e Ω(17)引入{}e 的角速度，定义为 {}{}112233tωωω∂=×=++∂e e e e e ω ω （18） 角速度性质为：{}{}{}112233123112233123123123{}[(,,)](){}[(,,)]{}[(,,)][(,,)]T T T T X X X X X X X X X X X X X X X ωωωωωω++=×=++=×××=+ e e e X e e e e X +X e e ωΩωωωΩΩΩ (19)利用角速度的性质，得{}{}1231111112311123113122{,,}({}[(,,)]{}[(,,)][(,,)]{,2,2})T T u v w w w w y z uv w w y z w y z ww w ωωωωωωω=+++++++=+−+− De e e e e e e e e e ΩΩΩ (20) 利用（17），（20），得以位移变分表达的动力变分形式{}{}{}12312121032(,,,,,,)d L u v w u T A A A v A A A A A A s w θθθαγγδδθδαδδδθδγδδθδγ⎧⎫⎧⎫⎧⎫⎪⎪⎪⎪⎪⎪−++⎨⎬⎨⎬⎨⎬⎪⎪⎪⎪⎪⎪⎩⎭⎩⎭⎩⎭∫= (21a) 其中{}{}1111123111213,,{,,}{}[][]{}[][(,,)][],,[]T T u v w A A A uv w I l l l ωωωωω=+−+ I R I R R Ω{}{}123122123212223,,{,,}[][]{}[]{}[][(,,)],,TTTA A A uv w l l l θθθωωωωω=+++ R I I I Ω{}12344123010203,,{,,}[][]{}[]{}[][(,,)]{,,}T T TA A A uv w l l l αγγωωωωω=+−+ R I I I Ω (21b)3121111411511620[][(,,)]d 0(()))0AI I w y z A I k I I sym ραγγ−⎡⎤⎢⎥==−−++⎢⎥⎢⎥⎣⎦∫I Ω 213()()[]000d 000Az y A ϕχχρ++⎡⎤⎢⎥=⎢⎥⎢⎥⎣⎦∫I2112233241251262341351362011021032332301102103222[][(,,)][(,,)]d ()()()()T Aw y z w y z AI I I k I I I k I I l k l l I I l k l l I ραγγαγγαγγαγγ=+−−−−−−−−⎡⎤⎢⎥=−+++−⎢⎥⎢⎥−+++⎣⎦∫I ΩΩ2141()()[][(,,)]000d 000Az y w y z A ϕχχρ++⎡⎤⎢⎥=⎢⎥⎢⎥⎣⎦∫I Ω343536242526000I I I I I I ⎡⎤⎢⎥=⎢⎥⎢⎥−−−⎣⎦01444514620221545515620364651662[{,(),()}d ]TA l I I I l H z y A I I I l I I I αγγρϕχχαγγαγγ++⎧⎫⎧⎫⎪⎪⎪⎪=++=++⎨⎬⎨⎬⎪⎪⎪⎪++⎩⎭⎩⎭∫ 1114151162121131231415116213214151162[{,2,2}d ]2()2()TA l I I I l w w w A I I I l I I I αγγρωωωαγγωαγγ++⎧⎫⎧⎫⎪⎪⎪⎪=−=−++⎨⎬⎨⎬⎪⎪⎪⎪++⎩⎭⎩⎭∫ 2122113121233435136232425126223435136220102103224251262301[{,2,2}[(,,)]d ]2()2() 2()()2(T T A l l ww w w y z A l I I I I I I I I I l l l I I I l ρωωαγγωαγγωαγγωαγγαγγωα⎧⎫⎪⎪=−⎨⎬⎪⎪⎩⎭+++++=++−++−++−∫ Ω021032)l l γγ⎧⎫⎪⎪⎨⎬⎪⎪++⎩⎭{}111213141516222324252633343536214445465556266111()()d ()()A I I I I I I I I I I I y I I I I z yz z y A I I I symI I z I y ρϕχχϕχχ⎡⎤⎧⎫⎢⎥⎪⎪⎢⎥⎪⎪⎢⎥⎪⎪⎪⎪=++⎢⎥⎨⎬⎢⎥⎪⎪⎢⎥⎪⎪+⎢⎥⎪⎪+⎢⎥⎪⎪⎩⎭⎣⎦∫(21c)6.2内能项变分分析Jaumann 应变矩阵变分为123111231([][])([(,,)])[][][][][(,,)][]T Tδδρρρδδδδδρρρ=⋅−⋅⋅′′+⋅⋅+⋅+⋅⋅⋅J A +A +B I B C A B λ∂Ωλλλλ−Ωλ （22a ）其中333333333333[(0,,)][] , [T y z ××××××⎡⎤⎡⎤⎢⎥⎢⎥⎣⎦⎣⎦a 0I A =A =0000Ω3333335622000[]0,,[]{0,,,0,0,0}0yz z y y z ρρρρ×××⎡⎤⎡⎤⎡⎤⎢⎥==⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎢⎥−−⎣⎦I C a =I C ,C =00 123123330000[(,,)]0000,[(,,)]00TTT yz yz y z ρρρρρρ×⎡⎤∂∂∂∂⎢⎥⎢⎥⎡⎤∂∂==⎢⎥⎢⎥∂∂⎣⎦⎢⎥∂∂⎢⎥∂∂⎣⎦0Ω∂Ω(22b)定义广义力为()123123123{,,,,,}{}[][]([(,,)])[d T T AF F F M M M A ρρρ=⋅+−∫∫S A +A B ∂Ω(23a) 123123{,,,,,}{}[T AF F F M M M A =⋅⋅∫∫S I B (23b)123132[()()]d Am S S A yzαρϕρϕ∂∂=++−∂∂∫∫ (23c){}{}322111d AMM M z y S A αχχϕ=++∫∫ (23d)113212321131313121112()d ,()d ,()d AAAm S w S w A m S w S w A m S w S w A =−=−=−∫∫∫∫∫∫(23e)其中1,41,421,211,3[][][][][][][][]ϕϕχχ−=−−−B =B B B (23f),[]m n x 表示36×的矩阵中第m 行第n 列元素为x ，其它元素为零。