預應力曲桿的Cosserat動力學模型

曹登慶，宋敉淘

（哈爾濱工業大學航天學院，150001 哈爾濱）

預應力曲桿的Cosserat動力學模型

曹登慶，宋敉淘

（哈爾濱工業大學航天學院，150001 哈爾濱）

為了更方便地使用Cosserat桿模型對各種細長結構進行動力學建模，在該桿模型的框架下通過引入描述固連在桿橫截面隨其一起運動的初始變形、初始橫截面轉動以及橫截面間的初始接觸力和力矩等變量，建立了考慮預應力的細長曲桿的Cosserat動力學模型.在得到的模型的框架下，基于相應的物理假設，從數學演繹的角度推導了橋梁建設中廣泛應用的具有初始垂度的拉索以及圓拱梁的動力學方程，結果與文獻中采用其他方法得到的相應的動力學方程一致.由于Cosserat理論采用了精確的幾何構形，這里導出的細長預應力曲桿動力學模型不僅保留了所有的幾何非線性特征，而且具有很好的普適性，通過它可以容易地導出實際工程問題中細長結構的非線性動力學模型.

Cosserat桿模型;預應力曲桿;具有初始垂度的拉索;圓拱梁

18世紀的Bernoulli和Euler認為梁之所以能抗彎是梁在運動過程中軸向伸縮的結果，且梁橫截面受到的彎矩正比于梁軸線的曲率.Euler-Bernoulli梁模型用梁軸線的撓度描述梁的運動，梁橫截面的轉動用撓度對空間坐標的導數來描述.Coulomb于1776年首次將力和力矩的平衡方程用于對梁微元進行力學分析，建立了較Euler-Bernoulli梁更為精確的動力學方程，并開始考慮梁承受扭矩的情況［1］.Kirchhoff于1859年計算了桿單元的應變，從能量的角度得到了桿的動力學方程，并提出了基于 Kirchhoff動力學比擬［2－3］的 Kirchhoff彈性桿模型.該模型忽略桿軸線的伸縮變形和橫截面的剪切變形，并將桿橫截面視為剛體，以描述桿軸線的3個線位移參數以及描述各橫截面的3個姿態參數為未知變量建立桿的動力學方程.Rayleigh于1877年引入梁橫截面慣性矩的影響，考慮軸向和橫向振動對梁進行建模，并指出在高頻振動時該修正是很必要的［1］.Timoshen-ko［4－5］引入梁橫截面間剪切力的影響對梁進行建模，并指出 Euler-Bernoulli梁和 Rayleigh梁為Timoshenko梁的兩種特殊情況.更加精確的桿的動力學模型的問世隨著工程設計和分析的需要而顯得越來越必要和迫切.

近年來，基于Cosserat兄弟和其他人的早期工作，形成了關于細長桿的Cosserat介質理論.Antman［6］在他的著作中運用該理論的基本思想建立了柔性桿和薄殼的Cosserat模型.Tucker和Wang［7］和 Cao 等［8－9］進一步完善了 Cosserat桿理論，并分別應用于對鉆柱的非線性動力學研究和MEMS中復雜結構三維非線性動力學模擬.Burton和Tucker［10］回顧了Cosserat桿模型的發展以及它近年來在各個領域中的應用.國內學者在Cosserat桿模型上也做了很多重要的工作.劉延柱等［11－12］從更偏重于力學分析的角度得到了任意初始構形圓截面Cosserat桿的動力學方程，在其框架下推導了Rayleigh梁、Kirchhoff桿、以及Timoshenko梁的動力學方程，并對基于精確的Cosserat桿模型螺旋桿進行了穩定性分析.Cosserat桿模型在Kirchhoff桿模型的基礎上考慮了桿軸向的伸縮變形以及橫截面間剪切力引起的橫截面的轉動，由于所作的假設更少，Cosserat桿模型較以往的桿模型都更為精確.

鑒于Antman所建立的 Cosserat桿模型［6］幾何形式的直觀性，本文的工作將在其框架下進行.本質上Cosserat桿模型適用于對具有任何初始構形和預應力的桿進行建模，但為了討論的方便，以往文獻中給出Cosserat桿動力學方程的同時，在Cosserat桿的本構關系中實際上都只將桿的由預應力引起的初始變形用1個抽象的符號表示，或者直接只考慮無預應力直桿(Non-Prestressed Straight Rod，NPSR)的情況［7－9］，而且在 Cosserat桿的動力學方程中，他們并沒有對初始量和動力學意義上的量進行區分.如何在Cosserat桿動力學方程以及其本構關系中加入桿的初始構形和預應力項以描述具有預應力的曲桿(Prestressed Curved Rod，PCR)的動力學行為是關鍵.本文的目的便是依據文獻中NPSR的Cosserat模型的基本思想，用Cosserat桿模型對PCR進行建模，并稱得到的模型為PCR的Cosserat動力學模型.在該模型的框架下，再根據具體的假定，便可以很明確地從數學演繹的角度推導各種特殊情況下的一維細長結構的動力學模型.本文中，矢量用小寫斜體黑體字母表示，如u，v和w;張量用大寫斜體黑體字母表示，如K，I和J;重復下標為啞指標，表示采用愛因斯坦求和約定;運算符?為兩張量間的張量積運算.

1 NPSR的Cosserat模型

考慮一根細長桿，在初始時刻t0時各橫截面質心組成的曲線(質心線)上定義拉格朗日弧坐標s，s∈［0，L0］，L0為桿的初始長度.如圖1所示，假設運動過程中t時刻桿的質心線為r(s，t)，分別以r(s，t)上的各點為原點定義局部的右手笛卡爾直角坐標系 s ∈ ［0，L0］|→{d1(s，t)，d2(s，t)，d3(s，t)}，并將坐標系固定在相應位置處桿的橫截面上.將d1(s，t)垂直于相應位置處桿的橫截面，將d2(s，t)和d3(s，t)置于相應位置處桿橫截面內.忽略桿運動和變形過程中密度和橫截面面積的變化，于是Cosserat桿模型的動力學方程可以由牛頓第二定律得到:

式中:ρ(s)為桿的密度;A(s)為桿的橫截面面積;n(s，t)和m(s，t)分別為桿橫截面間接觸力和接觸力矩;f(s，t)和l(s，t)分別為外力密度和外力矩密度(可以包括重力，阻尼力，電磁力，氣動力等);I(s，t)為截面的慣性矩張量.

圖1 Cosserat桿模型

由于Cosserat桿模型考慮了剪切力引起的桿橫截面的轉動，于是r'(s，t)不再始終垂直于相應位置處的橫截面，令

式(3)和(4)分別意味著(參見文獻［6－9］)

于是，Cosserat桿的本構關系可以表示為如下形式:

式(7)中的張量 K(s，t)和 J(s，t)以及式(1b)中的張量I(s，t)的分量形式分別為

式中:E(s)為桿的彈性模量;G(s)為桿的剪切模量;dA=dξdη，η和ξ分別為坐標系基矢量d2(s，t)和d3(s，t)對應坐標軸的坐標;Ψ(ξ，η)為橫截面的翹曲函數，由圣維南原理知其對于所有橫截面都相同［13］;κ(s)為 Timoshenko剪切系數，現有不同的理論對其進行物理意義上的解釋，Shames和Dym［14］系統探討過這個問題，但不論哪種解釋，κ(s)的值都取決于橫截面的形狀.張量I(s，t)和J(s，t)可以通過選擇 d2(s，t)和 d3(s，t)的方向達到對角化，以使得它們的形式變得更為簡潔.

下面闡述式(7)中NPSR的Cosserat模型本構關系背后的力學依據.式(7a)中出現的d1(s，t)表示t0時刻桿的質心線沒有剪切變形，只有沿著桿相應處橫截面法向的變形.由于采用了桿初始構形上的弧坐標，于是初始時刻質心線對弧坐標的微分總是為d.桿初始構形上任意一微元ds便可以用ds d作為其長度和方向的度量，它隨著桿的整體運動而運動，作為微元的初始狀態始終垂直于相應位置處的橫截面，于是變為ds d1(s，t).同時，變形后的微元為v(s，t)ds.于是微元的正應變(沿d1(s，t)的方向)以及剪切應變(分別沿著d2(s，t)和 d3(s，t)的方向)分別為

在 d1(s，t)，d2(s，t)和d3(s，t)上的投影.至此，便不難得到式(7a).式(7b)表示橫截面間的接觸力矩由橫截面間相對轉角來度量，對其分析可類似于對式(7a)所作的分析那樣，通過定義隨著桿整體運動的桿橫截面的初始轉角而實現.由于在前面定義的初始時刻局部坐標系間沒有相對轉角，于是初始時刻橫截面間無相對轉角.由上面分析可見，本構關系(7)的獲得，關鍵在于定義隨著桿橫截面一起運動的各種初始量，這給PCR的本構關系的建立及分析預應力對動力學方程的影響提供了思路.

2 PCR的Cosserat模型

假如細長桿為PCR，即桿的橫截面與橫截面間接觸力n(s，t)和接觸力矩m(s，t)在t=t0時不為零(即 n0(s)=n(s，t0)≠ 0，m0(s)=m(s，t0)≠0)，則在靜平衡位置處桿在預應力和初始外力f0=f(s，t0)和初始外力矩l0=l(s，t0)作用下保持平衡，由牛頓第二定律得

其中r0(s)=r(s，t0)為靜平衡位置時桿的質心線.桿質心線上各質點的線位移分為靜位移r0(s)和動位移rd(s，t)兩部分，即

式(11)和(12)中

在初始局部坐標系下的分量.ξij(s，t)=d(s)·dj(s，t).由第1節的內容知道u(s)由d(s)的定義決定.為了保證類似于式(5)的定義有意義，應將d(s)取為s個連續的單位矢量函數.曲桿無法做到像直桿那樣將d(s)獨立于s.不考慮桿運動過程中產生裂縫的情況，d(s)取為垂直于相應位置處的橫截面能保證d(s)連續.同時，取d(s)和d(s)為任意的連續單位矢量函數都不會影響到最后的結果.因為d(s)的選取影響著u(s)的同時，也影響著 di(s，t)，進而導致 u(s，t)－(s，t)的值不依賴于d(s)的選取.將外力(外力矩)分為初始力(初始力矩)和動態力(動態力矩)兩部分，即

將式(10)－(16)代入式(1)中，并將得到的兩方程分別減去式(9)中的兩個方程得到PCR的Cosserat模型的動力學方程為

當n0、m0和r0為零時，式(17)退化為式(1)的形式.

值得注意的是，拉索可以視為桿忽略彎扭剛度以及橫截面間剪切剛度的特殊情況，于是具有初始垂度的拉索的動力學方程可以由方程(17a)來描述，而且此時張量K需要退化成如下形式:

另外需要注意的是，前面的推導是在假定曲桿的橫截面、初始質心線、橫截面間初始接觸力以及初始接觸力矩為已知的前提下進行的.實際在一般情況下，對于曲桿，特別是具有預應力的曲桿，如何選擇桿橫截面并不是一件容易的事情.同時，求解桿的預應力也非易事.不過，這并非我們所關心的主要問題.而且，對于很多特殊情況，確定桿的上述初始量并不困難，作為例子我們將在第3節中采用本節建立的PCR的Cosserat模型推導具有初始垂度的空間拉索以及圓拱梁的動力學方程的過程中可以看到.

3 兩類特殊PCR的動力學建模

3.1 具有初始垂度的空間拉索動力學建模

如圖2所示，在全局右手笛卡爾直角坐標系{o;e1，e2，e3}下，具有初始垂度的拉索運動過程中的形狀為

其中 W1(s，t)，W2(s，t)和 W3(s，t)分別為拉索各質點在e1，e2和e3方向的線位移.

圖2 具有初始垂度的空間拉索的運動示意圖

由于不考慮拉索橫截面間剪切力的影響，于是d1(s，t)沿著拉索的切線方向，即

將式(18)代入式(2)得

全局坐標系下張量K寫成矩陣形式為

將式(19)－(21)代入式(12a)，并根據式(13)得

另一方面拉索的動應變e(s，t)為

于是由式(22)和(23)得

若拉索初始張力為H，則

將式(28)進行泰勒級數展開，忽略高于一階的小量并整理得

將阻尼從外力中分離出來，即

其中c為阻尼系數.將式(18)，(29)和(30)代入式(17a)中得到拉索動力學方程為

3.2 圓拱梁的動力學建模

如圖3所示，初始時刻圓拱梁圓心為o.如果不考慮圓拱梁的預應力，則可以通過如下方式選取圓拱梁的橫截面:通過o點，取垂直于圓拱梁所在平面的平面，并讓它與圓拱梁相截，則截出的圓拱梁的截面即可選為圓拱梁的橫截面.

圖3 圓拱梁的初始構形及坐標系的定義

在圓拱梁初始構形上對質心線每個質點s建立 右 手 笛 卡 爾 直 角 坐 標 系 {o;eθ(s)，er(s)，ez(s)}，其中eθ(s)沿著質心線在s處的切線方向并且指向s增加的方向，er(s)從o指向s，ez(s)=eθ(s)×er(s).可見，ez(s)實際上不隨s而變化.圓拱梁初始時刻質心線為

圓拱梁在運動過程中質心線的形狀為

其中Wθ(s，t)和 Wr(s，t)分別為 t時刻質點 s在eθ(s)和er(s)方向的線位移.取

記圓拱梁橫截面相較于初始構形轉過的角度為ψ(s，t)，并以逆時針轉動為正，于是

將式(32)和(33)分別代入式(2)，將式(34)代入式(6)，將式(35)分別代入式(5)和(6)，并根據式(13)－(15)得

由于不考慮圓拱梁的預應力，即

二維情況下，張量 K在坐標系{o;eθ(s)，er(s)，ez(s)}中的分量形式為

張量J和 I在坐標系{o;eθ(s)，er(s)，ez(s)}中的分量形式分別為

同樣地將阻尼從外力中分離出來，并將式(33)，(36)－(39)代入式(17)中得到圓拱梁的動力學方程為

若圓拱梁為小變形，則可以用式(40)的線性化形式來描述圓拱梁的運動，此時式(41)退化為

將式(42)和(8)代入式(40)，并忽略非線性項得

若忽略阻尼和外力，式(43)和文獻［16］中用廣義Hamilton原理得到的圓拱梁的線性動力學方程一致.

4 結論

1)建立了初始局部坐標系下桿各初始量(質心線的初始變形，橫截面的初始轉角，橫截面間的初始接觸力和接觸力矩)與運動過程中它們反映在PCR的Cosserat動力學方程和本構關系式中的相應的初始變量的具體關系.

2)將初始變量耦合進Cosserat桿動力學模型及其本構關系中建立了PCR的Cosserat動力學模型及其本構關系.

3)作為特例，橋梁中廣泛使用的具有初始垂度的空間拉索以及圓拱梁的動力學方程在PCR的Cosserat動力學模型的框架下作了推導.其過程體現了Cosserat桿模型包容性強以及數學推導的透明性和一貫性，簡化了繁雜的物理和力學層面上的分析過程.PCR的Cosserat動力學模型的建立使得各種特殊的一維細長結構的動力學建模變成了一個很純粹的數學推導過程，極大地方便了使用Cosserat桿模型對各種特殊的一維細長結構的建模.

［1］ARAFAT H N.Nonlinear response of cantilever beams［D］.Blacksburg，Virginia:Virginia Polytechnic Institute and State University，1999:17－18.

［2］LOVE A E H.A treatise on the mathematical theory of elasticity［M］.New York:Dover，1944:399 －400.

［3］劉延柱.彈性細桿的非線性力學［M］.北京:清華大學出版社，2006:14－18.

［4］TIMOSHENKO S P.On the correction for shear of the differential equation for transverse vibrations of prismatic bars［J］.Philosophical Magazine，1921，41(1):744 －747.

［5］TIMOSHENKO S P.On the transverse vibrations of bars of uniform cross section［J］.Philosophical Magazine，1922，43(1):125－131.

［6］ANTMAN S S.Nonlinear problems in elasticity［M］.New York:Springer-Verlag，1995:13－359.

［7］TUCKER R W，WANG C.Torsional vibration control and Cosserat dynamics of a drill-rig assembly［J］.Meccanica，2003，38(1):143－159.

［8］CAO Dengqing，LIU Dongsheng，WANG C H－T.Nonlinear dynamic modelling for MEMS components via the Cosserat rod element approach［J］.Journal of Micromechanics and Microengineering，2005，15(6):1334 －1343.

［9］CAO Dengqing，TUCKER R W.Nonlinear dynamics of elastic rods using the Cosserat theory:modelling and simulation［J］.International Journal of Solids and Structures，2008，45(2):460 －477.

［10］BURTON D A，TUCKER R W.Practical applications of simple Cosserat methods［C］//Mechanics of Generalized Continua.In:Advances in Mechanics and Mathematics.Berlin:Springer，2010，21:87 －98.

［11］LIU Yan-zhu.On dynamics of elastic rod based on exact Cosserat model［J］.Chinese Physics B，2009，18(1):1 －8.

［12］劉延柱，薛紜.基于精確Cosserat模型的螺旋桿穩定性分析［J］.應用數學和力學，2011，32(5):570－578.

［13］TIMOSHENKO S P，GOODIER J N.Theory of elasticity［M］.New York:McGraw-Hill，1951:275 －313.

［14］SHAMES I H，DYM C L.Energy and finite element methods in structural mechanics［M］.Chaman Offset，Delhi，India:New Age International Publishers，1991:197 －203.

［15］SOFI A，MUSCOLINO G.Dynamic analysis of suspended cables carrying moving oscillators［J］.International Journal of Solids and Structures，2007，44(21):6725 －6743.

［16］JAVID F，ESMAILZADEH E，YOUNESIAN D.An investigation into the vehicle-curved bridge dynamic interaction［J］.International Journal of Automotive Engineering，2011，1(3):235－243.

Dynamic modeling of pre-stressed curved Cosserat rods

CAO Deng-qing，SONG Mi-tao

(School of Astronautics，Harbin Institute of Technology，150001 Harbin，China)

To dynamically model slender structures conveniently by Cosserat rod theory，the dynamic equations for pre-stressed curved Cosserat rods are established within the framework of this rod theory by defining variables to describe their initial configurations，their pre-stresses and their initial deformations，which move along with their cross sections.Based on their corresponding specific assumptions，the dynamic equations for a cable with initial sag and a circular arch beam，which have wide applications in bridge structures，are explored within the framework of pre-stressed curved Cosserat rods，respectively.The results are in agreement with their corresponding equations obtained by other methods in references.Since the Cosserat theory is geometrically exact，the dynamic equations for pre-stressed curved Cosserat rods derived here are not only retaining all geometric nonlinear characteristics but also of universality，which can be effectively and efficiently exploited to nonlinear dynamically model for slender structures in practical engineering situations.

Cosserat rod;pre-stressed curved rod;cable with initial sag;circular arch beam

O326;O341

A

0367－6234(2012)11－0001－06

2012－02－12.

國家自然科學基金資助項目(10772056).

曹登慶(1958—)，男，教授，博士生導師.

曹登慶，dqcao@hit.edu.cn.

(編輯 張 宏)