第五章：大变形问题的有限单元法

第五章大变形问题的有限单元法1.弹性大变形问题的有限元法2.弹性分支点稳定问题有限元分析3.物质描述大变形增量问题的T.L、U.L法大变形时平衡方程和虚位移原理变形体初始和现时位形如图所示，以欧拉应力表述平衡时0,:ixijFVjjijinS:这是现时位形空间描述的平衡条件。在外荷为保守力系时VFVFiidd00SSiidd000ddVVJ经推证得上式乘以J，由于0,,:ixkXijFXVjk由复合函数求导数，平衡方程改写为平衡方程成为)(:0jkijkijkkijxXJXFxXXJV0)(0ijkkijFxXJX0:00ijjiijkkijFXFxXXJVliljijxXJ对任意j恒有0)(jkkxXJX0:ijjiNS证明321332313322212312111,,,,,,,,,,,,kjiijkxxxexxxxxxxxxJnlmijilkmnjkXxXxeexXJ21以j=1为例)(21123321nmnmkmnkXxXxXxXxexXJj0)()()()(2213231231233133223223332211XxXxXxXxXXxXxXxXxXXxXxXxXxXxXJXkk同理可验证，对任意j恒有0)(jkkxXJX上述拉格朗日和克希荷夫应力表示的平衡条件都是以初始位形作参考的物质描述。利用克希荷夫应力和拉格朗日应力的关系，可将平衡条件改为0)(:0ikikjjFXxSXV0:ikjijkNXxSSikkjijxXS0:00ijjiijkkijFXFxXXJV0:ijjiNS如果考虑到变形梯度和位移梯度间的关系，对比小变形情况，可见大变形时变形对平衡的影响，是通过变形或位移梯度表现出来的。0)]([:0ikiikkjjFXuSXVjiijjiXuXx则克希荷夫应力表达的平衡方程0)(:0ikikjjFXxSXV可改为设现时位形微小虚位移在V内单值连续、在位移边界上为零。则外力总虚功为考虑到位移边界虚位移为零和应力边界条件VSijijiiSunVuFwddeViiVuFwdeSiiSud空间描述的虚位移原理将平衡方程引入，考虑到虚位移微小，则VSijijijijSunVuxwdde利用格林公式，可得dedwwVwVijij变也即虚位移原理的虚功方程为VijijVSiiiiVSuVuFddd若虚位移用虚速度、虚应变用虚应变率代替，则VijijVSiiiiVSvVvFddd柯西应变虚功率方程为建立物质描述虚功方程，先讨论能量共轭关系。空间描述中又因为变形率张量是相对现时位形定义的柯西应变的速率ijijijijVwljkiklijXxXxSJ1jlikklijXxXxVE单位体积变形功率ijijESW因此克希荷夫应力和格林应变在能量上共轭。由于变形率和欧拉应力张量是对称的，因此jiijijijxvVwkikjijXxJ1再利用欧拉和拉格朗日应力间关系，可得ijijijijXuXvW因此拉格朗日应力和初始位形位移梯度在能量上共轭。由此即可得到物质描述的虚功方程为00000000dddVijijVSiiiiVESSuVuF00000000d)(ddVijijVSiiiiVXuSuVuF由虚位移原理可得eVeVFFBSESeeδ)(dδdTeEejTT0000VˆV式中为单元结点力矩阵，为单元等效结点荷载矩阵eFjeFEPd和PE分别为直接和等效结点荷载矩阵，R为综合等效结点荷载矩阵。上式也可写为RPPSBeVeEdTd00VˆeeSVeφNFNF0000SVddTTE按集成规则集装后可得eTd)(eVRSBU000VˆeδBEˆ1.弹性大变形问题的有限元法弹性大变形问题，需要考虑变形的非线性项和变形对平衡的影响。若以初始自然平衡状态作初始位形，则物质描述的格林应变为LijijjkikjiijijXuXuXuXuE)(21式中)(21jiijijXuXujkikLijXuXu21线性部分非线性部分为便于计算机编程，将张量转换为矩阵：格林应变矩阵和张量的分量间有如下关系T312312332211222EEEEEEE对应的克希荷夫应力矩阵和张量分量间关系为T312312332211222T312312332211222LLLLLLLT312312332211SSSSSSS引入两个算子矩阵T123312321000000000XXXXXXXXXAd,d,d,d,d,d,d,d,d,式中iXd,iXuIIIjiT332313XXXXuT123312321000000000XXXXXXXXXLuuuuuuuuuA,,,,,,,,,iX,uuiXT321uuuu再引入位移梯度向量的记号H3阶单位矩阵uXXXXXXXXXuHT333212121000000000在上述符号基础上，格林应变有位移表为uHAuAELL)21(则单元格林应变为eNueeLeBBBE~其中线性应变矩阵B和线性分析一样设单元位移场为ANBHNGGAHNABLLL2121非线性部分“应变矩阵”为LBBB~333323313232232131312311IIIIIIIIIHNGmmmXNXNXNXNXNXNXNXNXN式中G为如下9×3m的矩阵eGuHeNuT1T3T2T3T1T2T3T2T1000000000LA式中AL为如下6×9的矩阵iiuX由(AL)可见，格林应变-位移关系是非线性的。非线性部分因为eδBεδθθLLAA2为用虚位移原理建立单元刚度方程，还需有应变增量和位移增量间的关系。对线性部分321T1T3T2T3T1T2T3T2T1000000000θθθθθθθθθθθθθALθθθθθθθθθθθθθθθθθθθθθθθθθθθθθθθθθθθθθLA223T11T33T22T32T11T23T32T21T13T11T33T22T32T11T23T32T21T1uHAuAELL)21(所以综上所述，格林应变增量为eLGABEδ)(eLLLLδGAθAθAε)21(如果记)2(LBBBˆeδBEˆ，则。对弹性问题，在物质描述下本构关系(克希荷夫应力和格林应变)为pqklpqklEDSojqipmnijnlmkklpqxXxXDxXxXJDoeGeeLeBBBE~LBBB~21GABLL在上述基础上，由虚位移原理可得eVeVFFBSESeeδ)(dδdTeEejTT0000VˆV式中为单元结点力矩阵，为单元等效结点荷载矩阵eFjeFEPd和PE分别为直接和等效结点荷载矩阵，R为综合等效结点荷载矩阵。上式也可写为RPPSBeVeEdTd00VˆeeSVeφNFNF0000SVddTTE按集成规则集装后可得eTd)(eVRSBU000VˆeδBEˆ将克希荷夫应力表达式代入，可得e00T00)(Vd~ˆ)(eVeRUUKRBDBUeTd)(eVRSBU000VˆeVeBDB000TVd~ˆ式中K(U)是非对称的，为对非线性弹性问题、和都是位移的函数。Bˆ0DB~eTT)d()dd(d)(deVTUKSBSBU00Vˆˆ根据非线性方程切线刚度矩阵的定义，可得)2(LBBBˆLBBB~pqklpqklEDSoeeLeBBBE~根据本构方程，则有edddBDEDSTTˆ又因GABLL21)2(LBBBˆ式中eVeVeeSAGSB0000VVˆddddTLTT123312321000000000ddddddddddTLA所以0dB又因,所以T312312332211SSSSSSSdddddddddddTMSA311232333233121222313122111SSSSSSSSSL由此可得式中333323331323322312331312311IIIIIIIIIMSSSSSSSSSedddddG321eeVeVeeMGGSBddddTT0000VVˆeG基于上述说明，可得eTT)dd(d)(deVSBSBU00VˆˆUUKBDBMGGTeeVTe)d(dV)dˆˆ(00TT如果引入如下记号eVLTLTLLTeLBDBBDBBDBk00V)d(TTTeVeMGGk00VdTGABBLLL2LBBBˆeVTeBDBk00VdTe初应力或几何刚度矩阵线弹性刚度矩阵大位移刚度矩阵则单元切线刚度矩阵为eeeee)(LTkkkKeLeeeTkkkke“结构”切线刚度矩阵为建立了切线刚度矩阵，用牛顿法等即可求解。以上的讨论没涉及具体单元，因此具有普遍性。具体单元分析时，因形函数一般是自然坐标的函数，故需作坐标变换后代入相关公式，从而建立具体单元的切线刚度矩阵。牛顿迭代法解大变形问题的具体步骤为：0Ue牛顿迭代法解大变形问题的具体步骤为：eVLTLTLLTeLBDBBDBBDBk00V)d(TTTeVeMGGk00VdTeeeee)(LTkkkK以上是全量形式的弹性大变形分析，为了保证收敛，拟用增量迭代法。极值点失稳问题弹塑性问题2.弹性分支点稳定问题有限元分析研究弹性稳定问题时材料应力应变关系是弹性的，也即ijijkleklD。式中Dijkle对线弹性问题由弹性常数表示，对非线性弹性问题由割线或切线系数表示。稳定问题可分为两大类：分支点（第一类）稳定问题；极值点（第二类）稳定问题。分支点稳定问题存在稳定平衡状态和不稳定平衡状态的分界点，在此临界状态下，系统可在两种不同的变形形式下保持平衡。极值点（第二类）稳定问题变形形式整个过程中是不变的，但当荷载达到临界值时，变形将迅速增加使系统丧失承载能力或破坏。只讨论第一类线性弹性稳定问题。2.弹性分支点稳定问题有限元分析对于分支点稳定问题时，可设外荷是比例加载的。由于第一类稳定问