Dynamics of Aerospace Structures
Dynamics of Aerospace Structures ASEN 5022
Popular in Course
Mrs. Preston Lehner
verified elite notetaker
Popular in Aerospace Engineering
This 242 page Class Notes was uploaded by Laila Windler on Friday October 30, 2015. The Class Notes belongs to ASEN 5022 at University of Colorado at Boulder taught by Staff in Fall. Since its upload, it has received 55 views. For similar materials see /class/232177/asen-5022-university-of-colorado-at-boulder in Aerospace Engineering at University of Colorado at Boulder.
Reviews for Dynamics of Aerospace Structures
Report this Material
What is Karma?
Karma is the currency of StudySoup.
Date Created: 10/30/15
18 Classical and FEM Solutions of Plate Vibration Chapter 18 CLASSICAL AND FEM SOLUTIONS OF PLATE VIBRATION 18 2 18 1 INTRODUCTION Classical solutions of plate Vibrations have been a subject ofmany eminent scientists and engineers in the past A compendium by Leissa l is perhaps the most comprehensive source to date As the simply supported plates SSSSSSSS are the simplest to address we recall from the preVious chapter their solution Vibration mode shapes Wx y sin ax sin yy amanm m l2yann nl2 2 Frequency equation 193m 2 w g o 193m 2 n22 32 181 E113 D 121 112 where m n denote the number of harmonics along the X and ycoordinate directions In conformity with literature we will modify the frequency equation as 0 A w a2 MB A 18quot W 182 where a is the Xdirectional plate length Hence for a rectangular plate a b we have the following frequency equation 1quot 712 m2 112 183 We now summarize classical solutions of rectangular plate Vibration 182 ASSUMED VIBRATION MODE SHAPES FOR SIX BOUNDARY CONDITIONS The classical Rayleigh method assumes the plate de ections as the product of beam functions Wx7y YCVXX 184 each of which can be chosen depending on the boundary conditions 1821 Simply supported at x 0 and x a SS SS X x that satis es the boundary conditions 82W 0 32W W010Wya 0 185 ax can be eXpressed as Xx 2 5mm a m 234 186 It should be noted that if SS boundary conditions are imposed along y 0 and y b then x a m in the above eXpression is simply replaced by y b n 1 8 2 18 3 182 ASSUMED VIBRATION MODE SHAPES FOR SIgtlt BOUNDARY CONDITIONS 1822 Clamped at x 0 and x a C C 1 39 2 1 For even harmonics Xx cos y f E SlnIIl cosh ylf E a a sinh112 tany12tanhy12 0 m 2 4 6 187 F ddh X x 1 sin122 nh x 1 I oro armon1cs x s1ny2 a 2 sinh1239 yz a 2 tany22 tanhy22 0 m 3 57 1823 Free atx 0 and x a F F Form0 Xx 1 2x Form1Xx1 a 1 39 2 1 For even harmonics Xx cos y f cosh ylf a 2 s1nhy12 a 2 133 tany12 tanhy12 0 m 246 Si1112 2 SiI1h12 2 tany22 tanhy22 0 m 3 57 1 1 For odd harmonics Xx sin ylf Sinh 2f a 2 a 2 The above assumed de ection satis es only approximately the shear condition ie the zero shear force along the free edges 1824 Clamped at x 0 and Free atx a C F h Xx cos cosh I3x wxsinm x Sinh a a cos 13 cosh 13 a a cosygcoshy310 m123 189 1825 Clamped at x 0 and Simply Supported at x a C SS x 1 M Xxs1ny2z 2 sinh 22 tany22 tanhy22 0 m 234 x l sinh yZIZa 2 1810 1826 Free atx 0 and Simply Supported atx a F SS Form1Xx 1 f a Si1112 2 SiI1h12 2 tany22 tanhy22 0 m 234 l 1811 x l X 39 x SIMAla 2 2 x 39nh SI V22a Once again m n denote the number of Vibration nodal lines lying in the X and ydirections including the boundaries as the nodal lines except when the boundary is free We now list some of the frequency parameters Amquot from Leissa 1 for various boundary conditions for a square plate 1 8 3 Chapter 18 CLASSICAL AND FEM SOLUTIONS OF PLATE VIBRATION Table 1 Frequency Parameter Amquot mmquot a2 xpD for SS C SS C Square Plate M1 M1 412 Mz M1 Ms 28946 54743 69320 94584 102213 129086 Table 2 Frequency Parameter Amquot mmquot a2 xpD for SS C SS SS Square Plate M1 M1 412 Mz M1 Ms 23646 51674 58641 86126 100259 113217 Table 3 Frequency Parameter Amquot mmquot a2 xpD for SS C SS F Square Plate 11 03 M1 M1 Mz Mz M1 413 1269 3306 4170 6301 7240 9061 Table 4 Frequency Parameter Amquot mmquot a2 xpD for SS SS SS F Square Plate 11 03 M1 M1 Mz Mz M1 413 1168 2776 4120 5907 6186 9029 1 8 4 182 ASSUMED VIBRATION MODE SHAPES FOR SIgtlt BOUNDARY CONDITIONS Table 5 Frequency Parameter Amquot mmquot a2 xpD for SS F SS F Square Plate 11 03 411 412 413 421 422 423 98696 1613 3672 3948 4674 7075 Table 6 Frequency Parameter Amquot mmquot a2 0D for C C C C Square Plate 11 03 411 421 422 431 432 A41 3510 7290 10747 13163 16439 21035 Table 7 Frequency Parameter Amquot mmquot a2 0D for F F F F Square Plate 11 03 Fifth Mode Sixth Mode First Mode Second Mode Third Mode Fourth Mode 195961 242702 351565 636870 775896 134728 Frequency parameters for other boundary conditions can be found in Leissa 1 1 8 5 Chapter 18 CLASSICAL AND FEM SOLUTIONS OF PLATE VIBRATION 18 6 183 FINITE ELEMENT ANALYSIS OF PLATE VIBRATION Before we utilize a typical nite element software for plate Vibration analysis it is instructive to understand how the boundary conditions are treated Let us begin with the simplest case that is the geometric boundary condition For plate it is the clamped edge 8Wx7y Wxy 0 and an 0 1812 where n is the normal directional component from the clamped boundary edge X Fig 181 Plate bending element and degrees of freedom per node Referring to Figure 181 there are three discrete degrees of freedom per node w 0x 0y Since we have along the edges x 0 x a W W 0y the clamped edge boundary condition is satis ed if one chooses alongx0xa 100 and 0y0 1813 alongy0yb 100 and 0x 0 However simply supported boundary conditions must satisfy BZWOC y Wxy 0 and 0 1814 8112 This boundary conditions are approximately satis ed Via along x 0 x a w 0 and 0y is free to rotate 18 15 along y 0 y b w 0 and 0x is free to rotate 39 It is important to recognize that 82W 80 alongx0xa z y 8112 8x 2 1816 1 y 0 b 8 Wxy 80 a on g y 8112 By 1 8 6 18 7 183 FINITE ELEMENT ANALYSIS OF PLATE VIBRATION Hence the second of the simply supported boundary conditions are approximately satis ed as the element size a b becomes smaller and smaller This can be seen from the nite difference expression of 1816 m 39y N lt01 98 lt01 98 along x 0x a 8112 Bx 2a 2a 18 17 82m y 88 lt01 9f 05 912 39 alongy0yb z 8112 By 217 2b where the superscript designates the element node as shown in Figure 181 Similarly the free edge boundary conditions 82W 83 W B3 W 82Wxy8n2v 0 and MHz WM 0 1818 BS2 8113 811882 can be only approximately satis ed by leaving w 0x 0y as free variables along the edges Speci cally the comer conditions given by 2 m 0 1819 as an at the four comers x 0 y 0 x a y 0 x a y b x 0 y a are dif cult to satisfy and it is the level of violating these comer conditions that is largely responsible for somewhat slower convergence of FFFF plate higher modes Figure 182 shows the rst four modes and mode shapes of a square rectangular plate with 11 03 so that the FEMbased computations can be compared with classical solution tabulated in Table 7 The errors of FEM based modes are 75 0015 004 and 26 for the rst second third and fourth modes respectively The large errors of the rst mode is clearly due to not satisfying the zero shear and moment conditions along the four plate edges This aspect should be noted when modeling plates by the nite element method 1 8 7 Chapter 18 CLASSICAL AND FEM SOLUTIONS OF PLATE VIBRATION Lambda for mode 1 142761 Lambda for mode 2 195991 04 01 Y0 06 0 Y0 Xm 06 0 X m Lambda for mode 3 242802 Lambda for mode 4 360991 Fig 182 FFFF square plate Vibration modes and mode shapes 1 03 Reference 1 Leissa AW Vibration ofPlates NASA SP160 1969 1 8 8 ASEN 5022 Spring 2005 Dynamics of Aerospace Structures Lecture 9 15 February Vibration 0f String Bar and Shaft fxt T A WXt w X Tx dx W dxy I 8wxt 82wxt 8x 8x2 dx w T L K wLt lt gt u 9x t 3w dX Tx String in transverse Vibration Formulation Via Hamilton s principle I2 5Wn0nconsdr O L T V f1 L T2 000212 dx 1 O L V 2 Txds dx Kw2Lt O From the preceding gure we nd ds 1 wipe 01 dx m 1 w32 x tdx 2 Hence substituting 2 into 1We obtain L V2 you w xrdxKw2Lt 3 0 Transverse displacement wx t m Time derivative 2139 iii 1 m s Spatial derivtive wx 2 1 m m String tension T x N Mass per unit string length 0x kg m String length ds m End spring K Nm Since there are impotant differences between the varia tions of discrete model vs continuum model let s perform 8L below 5L5T 5V L 8T 5 0xu 2x r dx 0L 4 5V 5 you wipe r dx 0 8Kw2L a Evaluation off2 8Tdt L 8T5 0xu 2xtdx O L 2 px5u 2xtdx 5 O L 2 pxwxr5wxrdx O Evaluation off2 8Tdt cont 61 With 5T let s perform the de nite time integral of 5T 12 I2 L 5Tdt 2 0xlix t 521 X t dxdt 6 1 1 1 0 t x 0 t1 7 L 20 pxlixt 5wxt i f2 0xzIxt 8wxt dt dx f1 Evaluation off2 8Tdt concluded Since wx t1 and wx t1 are considered known every where over the spatial domain 0 f x f L we have 5wxt1 8wxt1 0 8 so that t2 5Tdt reduces to f1 12 I2 L 5Tdt 0xzIx t 5wx tdxdt 9 t1 I1 0 Evaluation off2 5 th Since the potential energy V does not involve time deriva tive of wx t it is adequate to carry out simply 8 V Viz L 5V 2 5 you wipe r dx 8Kw2L r O L 2 Tx wxx t8wxx t dx 10 O KwL t5wL t Evaluation off2 5 th cont d The rst term in the previous equation becomes L Tx wxx t8wxx t dx 0 L Tm wxltx 05w un ojff 11 L Tx wxxx t5wx t dx 0 Evaluation off2 5 th concluded L Tx wxx t8wxx t dx 0 wXxa t5wxa 10sz wXxa t5wxa txO L Tx wxxx t5wx t dx 0 12 Virtual work due to nonconservative force f x t L 5VTnoncons 2 fx9 051005 00135 13 0 Substituting 13 12 10 and 9 into 1 one obtains Hamilton s principle 1 2 5T V 5Wn0nconsdt 0 t1 1 u u Hamilton is principle for string 1 2 L f lt pltxgtwltx r I1 0 wxxxatfxat 5U39Xaltdx Tx wxx 1 Kwx 15wx txL Tx wxx t5wx txo gt611 0 In order for the above variational expression to hold rst the term inside the brace must vanish Which yields the governing equation of motion pltxgtagtltx r M wxxoc r foe rgt lt15 Similarly the two terms in the bracket must also vanish ltTltxgt wxltx r Kwltx rgtgt5wltx rgt1xL 16gt Tltxgt wxltx 05w rgt1xo 17 Boundary conditions for continuum string Two types of boundary conditions Essential or geometric boundary conditions When displacements are speci ed Natural or force boundary conditions When the boundary conditions come from force moment bal ance considerations Boundary conditions at x L From 16 we nd TL wxLt KwLtO 18 or wL t O 19 Since the right end is left to move 19 is not appropriate Hence the correct boundary conditon is given by 18 Which is a natural boundary condition Boundary conditions at x O From 17 we nd TO wx O t O 20 or w0 t O 21 Since the right end is xed 20 is not appropriate Hence the correct boundary conditon is given by 21 Which is an essential or geometric boundary conditon Boundary conditions concluded Observe that the variational equation viz Hamilton s prin ciple 14 provides the information to determine the boundary conditions depending on the end con guration or constraint condition This is a distinct advantage of the variational formulation of continuum models as opposed to Newton s approach Analogy Torsional Shaft and axial Roa String Shaft Axial bar Variable wx t m 6x t rad ux t m Stif iess Tx N GJxNm2 EAx N Hence the equations of motion for bar and shaft can be derived by appropriate parameter changes Strings in transverse vibration Equation 15 With f x t 0 offers the continuum eigen value problem for strings together the boundary condi tions 18 and 21 To this end we seek the solution of the homogeneous equation of 15 in the form of wx t Wx Ft 22 by employing the separation of variable used for the so lution of partial differential equations Substituting 22 into 15 18 and 21 we obtain d dWx ETx dx Ft W d2Ft O L x a lt x lt dW Tx dix KWxxL Ft o Wxlxo Fm 0 Applying the temporal part Ft in the form of Ft F e10 F 75 0 24 we obtain the following boundaryvalue problem set Tx a 0xWxO O ltx ltL 25 dx dx dWx Tx dx KWXxL 0 26 Wxxo 0 27 When the tension T x and mass 0x are constant along the string length equation 25 can be rearranged as d2Wx 2 2 18 W 0 dx 28 52w2pT Oltx ltL Whose general solution is known as Wx C1 sin x C2 cos x 29 Substituting the preceding solution into the x L boundary condition 26 yields T C1 cos8L C2 sin8L 30 KC1 s1n8L C2 cos8L O Similarly the x 0b0undary condition 27 yields C1 sin8 gtxlt O C2 cos8 gtxlt O 0 gt C2 O 31 With C2 O 30 reduces to T cos L i Ksin Lzo 32 Fixed end at x L When K gt 00 Viz gt 0 we have sin8L20 gt nl2 33 So that we obtain the frequency equation for the xed rightend case as k T wkzl k 12 34 L p Free end at x L When K gt O we have 2k 1 cos8LO zi jzl w 35 2L 0 So that we obtain the frequency equation for the xed rightend case as 2k 171 T 2 k1 2 36 60k I09 7 7 Note on the derivation of the string equation from elasticity theory Equilibrium 80139 a xi J 0 0 Constitutive law 37 0139 2 Ciij eke where 01 is the Cauchy stress elJ is the Eulerian strain is the body force and CUM is the constitutive constant Lagrangian straindisplacement relation 87 ale 1 Buk Buk Eil k123 J 28xj Bxl lax 8x Where xj refers to the undeformed con guration Hamilton s principle for continuum t2 6L SWmde o L T U t1 Yki uwwwod W 6 Uzgdnprzl V UampdW0J V0 Where Si is the second PiolarKirchhoff stresses V and V0 refer to deformed and undeformed volume respectively A continuum model for strings Sxx E Exx 8 Exx ux ui wi ux etc 8x T we w2dVo V0 L pXi12w2dx 00 MdAo 0 A0 U 2 SxxExx dVO V0 L 20 EAux ui w 2 dx Approximation of the string strain energy 2 2 ux ltlt wx 40 41 so that one arrives at 2 2 U EAuxwx dx EAui ux wi 025wi dx 42 22 1 2 2 EAux ux wx dx Ignoring the axial part we have A A Z A Z A L U Txw dx TxEAux 0L 43 T 0xu392 dx 0 which yields Via EulerLagrange s formalism the same governing equation of motion 15 together with the two boundary conditions 16 and 17 ASEN 5022 Spring 2005 Dynamics of Aerospace Structures Lecture 12 24 February Beams with Uncertain Boundary Conditions Consider a beam whose two ends are constrained by transverse springs kwl sz as well as rotational springs km 92 as shown in the gure below Beam with unknown boundary conditions EI mx 16021 s X 16 1 k w w w2 Fig 1 Modeling of unknown boundary conditions Before we proceed further it should be noted that the mathemati cally known boundary conditions can be realized by taking the limit values of the four spring constants as follows Simply supported k sz gt 00 and km 2 92 0 wl Cantilever beam x 0 kw1 k91 gt 00 and sz 92 0 Fixed xed kw1 k91 gt 00 and sz 92 gt 00 Freefree kwl 2 km 2 0 and sz 92 0 1 For intermediate values of the unknown springs a convenient way of identifying the boundary conditions is to invoke a variational formulation This is because the appropriate boundary conditions both natural and essential are determined as part of variational process Kinetic energy 8wxt L T O mx wxz dx wxtl T 2 Potential energy of the beam L Vb E1xwxz x Pxwxz dx 0 2 lt3 8 wx t 8wx 1 WC txx T 100 0x T where E I x is the bending rigidity and Px is the prestressed axial force Potential energy of four unknown springs Vs 2kw1w0t2 k91w0z 4 sz UJL 02 92 wLa 0 External energy due to distributed applied force fx 2 shear forces Q12 and moments M02 applied at both ends L 8W 2 fx 1 8wx 1 dx Q18w0t 0 M18w0tx Q2 8wL l 5 M28wLtx Q1 Q0 If Q2 QL 1 6 M1 M0l M2 MLl Hamilton s principle l2 8T 8Vb 8Vs8Wdt0 7 I1 Kinetic energy l2 L l2 8T dt mx wx tn8wx t dt dx 8 l1 0 l1 Variation of the internal energy of the beam Vb m E1 woe rgtxx6wltx M6 L EI wx txxx8wx 0H0 9 L EI wltxrgtxxxx P wltxrgtxx 6wltxrgtdx O Variation of the potential energy of the unknown boundary forces and moments 8VS kw1w0t8w0t k91w0tx 8w0 tx 10 sz UNLJ 5wLt l k92 wLtx 8wLtx 1 1 0 z 8Wd 5Vs Vb l28T 8 L 01pxp1 xm91 xf 0 11 z xm13 1 xmxw f 7 Z1 11 1px1 0m9 W x1 0m 921 1 0m13 Z1 11 1Px1 7m9 x1 7m3921 1 7m 13 Z1 Z1 11 1P1 0m9 15 1mm I 13 X 1 0m 13H Z1 11 1P1 7m9 Z5 1mm aquot x 1 7m 13 Z1 11 1PA49SA9 qA9 19 Z1 The governing equation of motion mx wxl zz E1 wxtxxxx fxt 0 13 The boundary conditions EI wL txxx sz wL t Q2 8wL t 0 E1 w0 txxx kw1w0t Q18w0t 0 EI wL 0m 192 wL 0x M2 8wL 0x 0 EI w0 0m k91 w0 0x M18w0tx 0 14 Free vibrations ofa beam The free vibration of beams can be modeled from 13 and 14 by setting PQ1Q2M1M2f3 l 0 15 so that 13 and 14 are simpli ed to 17100 wxl zz E1 wx txxxx 0 El wL zxxx sz wL z 8wL z 0 El w0 zxxx kw1w0t8w0t 0 16 E1 wL 0m 92 wL 0x 8wL 0x 0 El w0 txx k91w0tx3w0tx 0 Free Vibration Governing Equation and Four Natural Boundary Conditions mm W on E1 W mm 0 HM wltL rgtxx1x sz ML 0 0 HM wlt0 rgtxx1x kw1wltorgt 0 17 W ML 0 k92 wltL 0x 0 W wlt0 om k91wlt0rgtx 0 61 81 13 0 x7A4 79 701 13 0 7M4 Tm 7M4 13 0 Z x0A1 W 001 0 0A1 W xxx0A1 3 0 MOOM 30M mzm Sme m 01H pmmnsoms 911M ImIIM Mama 1 xm p zuoa suozzzpuoa Ixnpunoq zommu mod The frequency equation of vibrations of a uniform beam can be obtained by assuming the solution in form of W 01 sin Bx 02 cos Bx c3 sinh Bx c4 cosh Bx E1 WC 2 8 6100s 8x 02 sin Bx 03 cosh Bx C4 sinh Bx Wxx 82 01 sin Bx 02 cos Bx c3 sinh Bx c4 cosh Bx Wxxx 83 61 cos 8x 02 sin Bx c3 cosh Bx c4 sinh 8x 20 where 33 k61 53 cos5 Ew2 sin 5 5sin5 I8 53 sin5 l w2 cos 5 4302 sinh 5 5 cos5 kg2 cos 5 IE62 sin 5 1 I33 4391 53 cosh5 5 sinh 5 k62 cosh 5 1391 keiE1L 53 sinh 5 kw2 cosh 5 5 cosh5 k62 sinh 5 I wi kwiEI A Cl 02 C3 C4 90 21 22 Freefree beam kwl sz 2 km 2 92 0 63 0 3 0 c1 0 B 0 g 02 wooss 33mg B3cosh6 53th 03 0 L B sin5 B cos sinh B cosh J 04 U 23 1 0 1 o 1 det I 0 1 0 1 0 L cos 6 sin 8 COSh 6 Sinh 6 J sin 6 cos 6 sinh 6 cosh 6 Cantilever beam km 2 92 0kw1 gt ook91 gt oo 33 k 8 3 k w1 w1 k61 IB k61 I8 det cos 8 sin 8 cosh 8 sinh 8 Z 0 24 sin 8 cos 8 sinh 8 coshB det cos 8 sin 8 cosh 8 sinh 8 25 sin 8 cos 8 sinh 8 cosh 8 u 1 0058 coshB 0 which gives 3 1875 4694 7855 26 120 Fixed end rotational spring parameter for cantilever beam rotational spring k theta ElL Error in fundamental frequency percent Effect of xed end rotational spring on fundamental frequency of a cantilever beam Simply Supported and FixedFixed Beams w0 wL 0 0r kwl sz gt 00 0 1 0 1 01 4301 g E01 396 02 det I sin 6 cos 6 sinh 6 cosh 6 C3 sin 6 cos6 6 sinh 5 cosh j C4 402 cos 5 1302sm5 k02cosh6 4025th 27 U 6 2 sin sinh 1301 sin cosh cos sinhB l gz Boos sinh 13911392 1 cos coshB 0 Several ideal cases can be obtained from 27 Simply supported ends 1 O 92 0 gt sin8 sinhB O Clamped at x O and simply supported atx L 1391 gt oo 92 O U 28 tan 8 tanh 8 O Clamped at both ends 1391 gt 00 Egg gt 00 gt l cos cosh zo When the boundary conditions are not ideal Viz o lt 1391 lt 00 and o lt 1392 lt oo 29 one needs two or more measured frequencies to determine km and keg Theoretical Basis Same as the continuum beam moa el Kinetic energy L awx r T f0 mx wx t dx wx at T 1 Potential energy of the beam L 2 V EIx wow dx 0 2 82wx r awx t wx txx T wxatx 2T External energy due to distributed applied force fx 2 shear forces Qavz and moments Mag applied at both ends L 3W 2 fx 1 8wx 1 dx 0 Q18w0t l Q2 8wL t M1 31120 tx M2 8wLtx 3 Q1 Q0 1 Q2 QL 1 M1 M0 1 M2 ML t Potential energy of four unknown springs nggemwm fleoo 4 hawanenme Hamilton s principle f2 8T 8Vb 8Vs8WdtO 5 f1 From here on FEM modeling differs from continuum modeling Key departure in FEM modeling of beams om continuum modeling 1 Instead of carrying out the variation the 3process of Hamilton s principle5 in terms of the continuum variable wx t the energy expressions T V W etc are approximated on a completely ee beam with an arbitrary length 6 area A the bending rigidity E I and mass density 0 2 The interpolation of the transverse displacement wx 2 over the completely free beam segment 0 5 x 5 E is chosen to satisfy the homogeneous differential beam equation ie txxxx 0 if possible If not it is chosen such that as the length of the beam gets smaller and smaller it approximately satis es the homogeneous equation Approximation of wx t over the element segment 6 ltlt L wx r 600 61tx 62tx2 C3tx3 7 Observe that the above interpolation satis es the homogeneous differ ential equation of beam 6 How does one determine the coe cients co cl 02 C3 The answer comes from the nite element method T he FEM beam bending element 1 Observing that wx t and wx 0 are linearly independent from the result of variational formulation we specify them at nodes 1 and 2 in the previous gure Namely two discrete variables at node 1 wx1 t and wx1 0 two discrete variables at node 2 wx2 t and wx2 0 8 2 For convenience we introduce a local coordinate system or ele mental coordinate system x y and set x1 O and x2 E Substituting 8 into the interpolation function7 we obtain w1t wx1 0 t c0t c1t0 c2t02 c3t03 610 wx1 0 ox c1t 2c2t0 3c3t02 w2t wxz E t cot c1t 62t 2 63t 3 62t wx2 e tx c1 t 02t2 c3t3 2 9 Solving for the coef cients and backsubstituting into 7 one obtains woe r H1Sw1t H2 910 S x6 Hm mm P149920 H1lt gtlt1 3 2S3 H2 S 2S2 W H3lt gt lt3 2 2S3 Hm lt 2 56 10 Element Discretizatz on Elemental kinetic energy 5 Tel mx wx tt2dx 0 1 12 0 me wltargt da x ES Elemental potential energy E V51 EIx wx 03x dx 0 13 1 E1 f0 E 3w at gd Discretizaiion of Elemental Kinetic Energy 1 T6 1 0 me we r id 1 1 AE39 TNTN d 20 0 W qtS 14 1 cf 0 Mi NT N 61S 610 qltrgtT me 610 NIH NIH where we utilized 11 viz 39 3 woe r Z ciltrgtxi Ne qltrgt 15 10 Discretl39zatl39on of Elemental Strain Potential Energy 1 E1 V5 2 O 6 MS 0 d 1E1 2 0 qu Ngg NggQl d E3 16 HNI 1 E 1qltrgtT 0 EN N 61S qltrgt qu ke qltrgt where we utilized we Ogg New qltrgt lt17 Elemental mass matrix me Elemental stt ness matrix kel 12 66 l 12 M E el k 3 66 462 6E 262 12 6E 12 6E 18 19 Approximation of elemental external energy 3 l 1 W12 fxl8wxldx fltsrgtNltsgt6qltrgt d5 0 0 Q18w0t Q2 8wLIM18w0lx M2 8wL ox 3 Let s utilize the approximation10 wx I MS 10 wot 0x NSs 10 20 and observe that w0 I 6110 wamp I 6130 w0 0x 6120 wamp m 6140 21 so that we have 5ng mewOZ NTfx t ms Q18q1l Q2 8q3l M1 8q2l M2 8q4l 22 l 6qltrgtTFl F14 NTfltxrgt d5ltQt M1 Q2 M2gtT 0 Summary of elemental energy expressions T qltrgtlte gtT me qltrgtlte gt Vb qltrgtlte gtT ke qltrgtlte gt 23 8W aqu 1 Question How do we model a long beam with beam elements Answer 1 Partition the long beam into small elements This needs to be explained 2 Generate the elemental energy expression This is done see 23 3 Sum up the elemental energy to form the total system energy 4 Perform the variation of the total system energy to obtain the equations of motion Let s now work on items I 3 and 4 I 8 Establish the relation between the elemental and assembled global degrees of eea om lt1 11 u1 For element 1 q q1 2 2 11 I12 For element 2 qm 12 3 3 1 3 4 1 For element 3 qm 12 24 For general case gt ll 2 L g 3A Sum up the total system kinetic energy 11 Ttotal 2 Tel el1 qlt1gtTm gtqlt1gt qlt1gtTmlt qlt1gt qltngtTmlt gtqltngt 11 T m 0 0 q1 T 25 q2 0 mlt2gt 0 q2 0 0 0 qn 0 0 mm qn TWO aqt m qt 33 Sum up the total system strain energy 11 total el V Z Vb el1 qlt1gtTk gtqlt1gt qlt1gtTk gtqlt1gt qnTk Qn qlt1gt T k 0 0 mm T 26 q2 0 klt2gt 0 q2 0 0 0 qn 0 0 k0 qn Vtotal q T k Qt 3C Assemble the partitioned elements back into a beam this corresponds to substituting the partitioned elemental degrees of freedom qg by the assembled or global degrees of freedom ug via 24 c Total assembled system kinetic energy Tmquot W m qt an LTm Lug g1 LTm L ug 27 total 1 T T Similarly we have Vtotal q T k qg ugT LTk Lug ugT LTk L g 28 rm 2 u Kg ug Kg 2 LTkL 0 Total assembled system external work 11 el1 3th ft sug Lng 29 l T 8W 2 Sug fg EulerLagrange 3 equations of motionl Lagrangian L Ttotal Vtotal 30 1 T 1 T ug ug From the generic equation of motion at BL BL Qk 31 dl aqk aqk we obtain the FEM equation of motion om 3 O and 31 Mg iigKg ungg 32 20 Modeling of Dynamical Systems Continuum Models Chapter 20 MODELING OF DYNAMICAL SYSTEMS CONTINUUM MODELS 20 2 201 INTRODUCTION Let s consider the following scenario A construction company has been designing a cable that will run between two towers in a mountain After a preliminary analysis and a scale model test a technical evaluation team concluded that a third exible tower may have to be placed To this end we will needto formulate the equations of motion for the cable system complete with the appropriate boundary conditions Assume that the original two supporting towers are substantially stiffer than the new one to be placed between the two endsupports of the original cable system The design team has concluded that they should model not only the exibility but also the inertia effect of the middle tower In addition the new system is to have its fundamental frequency about 50 higher than the original system which was the reason why a third tower is needed Observe that the cable tension T of the original system could not be increased any further as that will result in an unacceptable bending on the original two towers We will now proceed the formulation steps as follows Figure 201 Cable with a middle support There are two ways of selecting the coordinate system Once is to use one coordinate system for the entire system An alternative is to employ two independent coordinate systems one for the left cable and another forthe right cable as shown in Figure 20 1 above In what follow two independent coordinate systems will be adopted 202 GOVERNING EQUATIONS OF MOTION Introduce wx t and 1110 t whose domains are de ned by Forthe left cable wx t 0 5 x 5 l 201 Forthe right cable wEt l 5 E 5 L l L2 20 2 20 3 202 GOVERNING EQUATIONS OF MOTION Using these distinct displacements we have t2 8L 8Wnommdt 0 L T V 1 Z L T 2 pxw2xtalx pcw2itali 0 Z 202 We w20 t ML u L t Mm u l t I Z L V Txw xidx TEw itali 0 Z Ko w20 t KL w2L t Km WI 0 Evaluation offZ 8Talt f2 f2 Z L 8Tdt px111xt8wxtalx pE1DEt8wJEtdJE f1 f1 0 Z M01110 t 61110 t ML111L t 8wL t mea t ML 0 dt 203 Evaluation of 8 V Z L 8V 2 Tx wxxx t8wx t alx TOE wgg t8wc t ab 0 Z T00 wxx7 03111067 txZ T00 wxx7t3wxitxo 204 TO C W205 0511101 0124 TOD 11120570511105 0122 K0 wxi t3wxitxo KL 111053 051110570124 Km 11106703111067 txZ Virtual work due to nonconservative force f x t SWHOVICOVIS 0 due to free Vibration problem of interest Substituting 205 204 and 203 into 201 one obtains the following Hamilton s principle 20 3 Chapter 20 MODELING OF DYNAMICAL SYSTEMS CONTINUUM MODELS f2 5 lt t mxmoc r Toe 1mm awe mdx f1 0 L t pltxgtwxr Tltgtewemgtawirgtdx l Tx WW I Kowxt t MoUJxt 01511106 tx0 206 Tx wxx t mex t mex t8wx tx TJE 1020 01611106 t2t TOE WOW KLwOE7t MLwft6wxttxL dt 0 Since 6wx tx and 811106 med are the same Viz wx txl wOE tlt 0 207 the fourth and fth rows in 206 must be concolidated to read Tx wxx t mex t mex t8wx tx TJE w2it8wf7txt 208 Tx WWI mextt mextt TOE w2ftt3wxttxz Tx WWI TOE11120570mexttmextt3wxttxz Substituting this into 206 leads to the following variation equation r2 z lt t mxmoc r M wave m we make r1 0 L t poawoar Tltgtegt ngarnawmndx l Tx WWI K0wxt t M01Dx7 tltSUJO tx0 Tx WOW I TOE 11120 67 I mext t Mm11gtxt tltSUJO txl TOE11120570KLwJEJML11gtJEJ5UJXJXL dt 0 209 The preceding Hamilton s equation yields two governing dilTerential equations 20 4 20 5 203 VIBRTION PROBLEM The governing equation of motion 1000111067 t TOG lumbar 0 S x S l l L2 2010 100011105 t 90111220670 5 S x S L The three natural boundary conditions 400 mm K0wx r Mowxrawx 000 0 2011 Tx wxx t TOE wE t mex t Mm111x t8wx tx 0 2012 l m Woe t KLwE t ML111E 01311106 txL 0 l 2013 Together with the above three natural boundary conditions one must augment them with the con tinuity condition at the middle support given by 207 which provide the critical fourth constraint condition The three natural boundary conditions plus the continuity condition at the middle support T WAX I KowOC I M0wx7 tx0 0 T WAX t T 11120570mex7tMm1Dx7txz 0 l L2 T Woe t KLwE t ML111EIXL 0 2014 111057 txl 1110670124 0 with T Tx TOE 203 VIBRTION PROBLEM The general form of solution from the governing two equations of motion 2010 wx t Ft Wx Wx C1sin9x C2 cos 9x 0 5 x 1110 t Ft WOE W0 2 C3 sin99E C4 cos 95c l 5 9E Ft fejw S E l L 2015 When K 0 gt 00 K L gt 00 by substituting this condition and the general solution form 2015 we nd 20 5 Chapter 20 MODELING OF DYNAMICAL SYSTEMS CONTINUUM MODELS W0 0 TWWL TW kKMW waW L4QZL WL W2l 0 Womd Wamuzo Let s focus on the second equation of the above equation set First dividing by the tension T and multiplying by l we have szml ZWWL lW k W T W0m4o Second we note that with 3982 1 szml 02 Mml2 wzp Mml2 2 Mm 2 Mml I3 9l T T pl T pl pl 2 Mm 32 Mm Mmz v ly L2 pl Third introducing Kml Km 2 T equation 2017 can be written as HWdL lwmrwmmm m womdo Therefore 2016 simpli es to W0 0 HWdL lwmrwmmm m womdo WL W2l 0 Womd Wumuzo 20 6 20 6 2016 2017 2018 2019 2020 2021 20 7 Finally using Wx C1 sin8x C2 cos8x 0 W0 2 C3 sin8x C4 cos8x l S S tlon 0 3 COS3 Km M32 Si113 0 sino coso 7 sino x 5c 203 VIBRTION PROBLEM 1 3 Si113 Km Mm32 C053 3 00573 0 sin2 2022 the preceding boundary and constraint set 2021 leads to the following characteristic equa 0 0 0 0 20023 By successive manipulations of the third and fourth columns it can be shown that the above equation simpli es to 0 l 3 0053 Km Mm32 Sin3l 3 SiI13 Km M132 0053 0 0 sin 0053 AC0 13 sin3 0 20024 This simpli cation is equivalent to taking the left and right cable coordinates independently for bothcases as 0 5 x 5 l 0 5 E 5 Setting detA 0 one nds with um Kml 2 cos8l T W 40 sin8l Km Since the desired frequency is 2157 gt M157z2 we nd the spring constant to be Kl Km 2 Tm 1026404145599745 2025 2026 2027 The fundamental mode shape is plotted in Figure 202 using the computation routines listed below 20 7 Chapter 20 MODELING OF DYNAMICAL SYSTEMS CONTINUUM MODELS 20 8 Two cables supported by a exible middle suppor l l l l l l middle spring km 10264 beta L 47124 Mode sh ape l 05 06 Beam span Figure 202 Fundamental mode of a cable with a middle support a Computations of two cable vibrations K using a general characteristic equation clear all X Km determined from the detCcable km1026404145599745 Z compute mode shape for this frequency betaell 15pi2 Z target frequency Z compute Ccable 4x4 matrix Ccable CmatrixTwoCablekm betaell X it is assumed that the minor Cbeam8x8 matrix is nonsingular Cminor Ccable24 24 bCcable241 Z mode shape coefficients assuming c4 1 coef Cminorb xLcoordzeros10 yLcoord zeros1 O xRcoordzeros10 yRcoord zeros1 0 20 8 20 9 204 A COUNTRYSIDE ENGINEER39S MODEL for i0100 span ilOO arg betaellspan sx sinarg cx cosarg WL sx coef1cx WR coef2sx coef3cx xLcoord xLcoord spaI yLcoord yLcoord WL xRcoord xRcoord spaI yRcoord yRcoord WR end xcoord000005100 ycoordyLcoord1100 yRcoord Knormalize ycoord ycoord ycoordmaxabsycoord figure1 plotxcoord ycoord xlabel Beam span ylabel Mode shape legend middle spring km num2strkm betaL num2str2betaell title Two cables supported by a flexible middle support Z end of the program Z compute the determinant of twocable with middle support function Ccable CmatrixTwoCablekm betabar Z characteristic matrix for two cable with middle flexible support Zvibration problem X the two end supports are fixed Zinput beta betaL2 km middle spring m is assumed Zto be half of the cable weight Zoutput Ccable 4x4 K K data r K c cosbetabar s sinbetabar Z Ccable O 1 O O betabarckmsbetabar 2s betabarskmcbetabar 2c betabar O O O s c s c 0 1 return 20 9 Chapter 20 MODELING OF DYNAMICAL SYSTEMS CONTINUUM MODELS 20 10 204 A COUNTRYSIDE ENGINEER S MODEL Suppose that you are in the countryside without access to advanced modeling tools and you are asked to design the same task If that individual happens to be this instructor here is an approach Step 1 Model each half of the cable by a lumped massspring system M K such that the fundamental frequency is col 27 4KM This is illustrated in Figure 203 Mm I Left cable Right cable EquValent Simple mOdEl Equivalent simple model ll Km ll K O K M M Left cable 4 Right cable O The spring distribution factor a is chosen such that when Mm Km 0 the fundamental frequency is halfof the individual frequency viz from 2pi to pi Thusa is found to be a 075 Figure 203 A simple equivalent model of a cable with a middle support Step 2 Form the left and right cables to be in parallel and attach them to the support system characterized by Mm K m Step 3 In doing so distribute the cable springs such that if the middle support is removed the fundamental frequency of the system would be half of the individual frequency This is because for cables the fundamental frequency is given by mr T col 2028 L p which states that if the cable length is doubled the fundamental frequency n 1 will be halved A simple weighted distribution of the cable springs indicates that the distribution factor should be a 34 2029 20 10 20 11 204 A COUNTRYSIDE ENGINEER39S MODEL Step 4 Form the model equation set given by K Mxl Kx1 19 fl K sz Kx2 19 f2 2030 K K mem Km 3xm 1x1x2 fm whose characteristic equation is given by K Ma2 0 det 0 K M 02 2 0 2031 K K K T T Km7 Mm 602 For computational expediency divide each of the three rows of the above characteristic matrix by M and using aquot 4KM we have a 02 0 l det o a 02 41 0 2032 12 DZ 12 Tquot Tquot KmM 7quot MmM wz Introducing the following nondimensional parameters Km Km K 2 Mm V M Kmaquot Mm V 2033 equation2032 becomes a 02 0 l det o a 02 40 5 0 2034 12 12 2 124 2 Tn Tn men l39Tn me For this characteristic equation one nds Km for the case of a 157 with um l to be Km 2 13747 2035 which when compared with the exact solution given by 2027 yields about 34 off One possible improvement for estimating the middle support spring is to distribute the spring constant K m as shown in Figure 204 20 1 1 Chapter 20 MODELING OF DYNAMICAL SYSTEMS CONTINUUM MODELS 20 12 Left cable Rig ht cable Left cable Right cable A Possible Improvement of estimating the middle support spring Km by distributing Km into three ways as shown on the right in the above figure Figure 204 Possible Improvement in modeling the support spring Whether this conjecture will be valid or not is left as an exercise to curious minds Better yet one may employ each half of the cable by a twoDOF model instead of the oneDOF model Iused herein In any case you now have a glimpse of oldfashioned modeling approaches Let me know if you would like to learn more about oldfashioned modeling practices 20 12 ASEN 5022 Spring 2004 Dynamics of Aerospace Structures Lecture 11 February 19 Cables with End Masses and End Springs Cables with End Masses and End Springs Formulation Via Hamilton s principle I2 5L l 5Wn0nconsldt 07 L T V f1 L Tzf 0xu 2xt dx 0 1 M0 w20 t ML w2L t L V 2 Tx wibc t dx Ko w20 t KL w2L t Symbol De nitions Transverse displacement wx t m Time derivative w iii 1f m s Spatial derivative wx 2 1 m m String tension T x N Mass per unit cable length 0x kg m End springs K0 KL Nm End masses M0 ML Nm Evaluation off2 5 T dt f2 5Tdt 2 f1 l2 L f f0 pltxgtwltxrgt5wltxrgtdx 2 1 M0110 t 5w0 t MLIML t 5wL t dt Evaluation 0f5 V 5V L Tx wxxx t5wx t dx 0 Tx waC 0511106 txL 3 Tx waC 0511106 tx0 K0 wx t5wx 0sz KL wx t5wx tx0 Virtual work due to nonconservative force f x t L 5anoncons A fx7 t5wx7 001x 4 Substituting 4 3 and 2 into 1 one obtains Hamilton s principle l2 5T V 5Wn0nconsdt 0 r1 1 u u Hamilton 3 principle for cable with end masses and end Springs l2 L f lt pltxgtwltx r t1 0 Tx wxxx t fx t 5wx tdx Tx wxxtKLwxf 5 MLIIJOC t5wx 0sz Tx waC t KLwOCa 1 M011x t5wx tx0 dt 0 The governing equation of motion pxi11x t T06 wxxOC 1 fx t The boundary conditions Tx 111x06 t KLwOC t MLI39DOC t5wx txL 0 Tx wxx t KLwOC t tx0 0 6 7 8 Observations 1 The governing equation of motion for the case of end masses and end springs is the same as the case Without any end conditions This means the governing equation of motion is the same for all possible boundary condi tions 2 Therefore the general form of solution wx t FtC1sin8x C2 cos 8x 9 is applicable to all the cable Vibration prob lems Boundary conditions at x L From 8 we nd TL wxL i KLwLt l ML139DLt 0 10 or wLt0 11 Boundary conditions at x 0 From 9 we nd T0 wx0 t K0w0 t M0110t 0 12 or w0 t 0 13 From 1013 one nds that the two essential boundary conditions given by 11 and 13 are special cases of the two natural boundary con ditions given by 10 and 12 since the former two are obtained from the latter two if K L gt 00 and K0 gt 00 respectively Substituting 9 into 10 and 12 together with F t 2 F610 yields the following characteristic equation B cosB l 3SiQ3 C1 0 KL uL zmn KLLL ZNOS ll l 14 3 K0 0182 J C2 0 I0 Z L a T K L K L KL L K0 2 0 15 T T M L ML ILL L M0 0 The characteristic equation of 14 is found by taking the determinant to be zero K0 uo 39om um sin 39 quot cos 39nxo m2 KL we 16 B z Sing 2 0 We now examine several special cases in the fol lowing I FixedFixedEnds kL gt 00 k0 gt 00 This case corresponds to kL gt 00 k0 gt 00 so that by dividing 16 by kL k0 yields sin zsin L 0 gt BL k71k 12 ll k T wk k123 L p 17 which was already discussed previously 2 FixedFree Ends kL gt 0 k0 gt oonusL M0 0 This case leads to the following characteristic equation from 16 cos8 sin L 0 2k 1 u 18 2k 1 T wk k123 2L 0 Which was also treated before 3 Flexible Supports IQ 2 kg 75 0 LL 2 u0 0 This case leads to the following characteristic equation from 16 k2 8 2tan8 2k3 0 k kL 2 kg 19 Other general cases can be evaluated using 17 This is left for exercises Mode shape of the first natural frequency of cable with flexible supports 1 l l l l l l l 0 l l l mode shape 0 D l l 05 7 7 04 7 7 7 frequency in Hertz 01039 spring k 10 betaL 26277 03 7 7 l l l l l l l l l 0 01 02 03 04 0 5 06 07 08 09 1 cable span Mode Shapes of Cable with End Springs Mode shape of the first natural frequency of cable with flexible supports 1 l l l l l l l o m l l mode shape 0 0391 l l 04 7 7 03 7 7 02 7 7 01 7 frequency in Hertz 01242 spring k 10000 betaL 3141 0 l l l l l l l l l 0 01 02 03 04 0 5 06 07 08 09 1 cable span Mode Shapes of Cable with End Springs ASEN 5022 Spring 2005 Dynamics of Aerospace Structures Lecture 15 March 10 Finite Element Modeling of Beam Bending Vibration Let s consider the modeling of a beam under general boundary conditions by the nite element modeling A W Beam with unknown boundary conditions k k M1 2 EI mx S S X k lt L Theoretical Basis Same as the continuum beam moa el Kinetic energy L 8 t T mx wxz dx wxtl 1 0 Potential energy of the beam L Vb E1x wxz x dx 0 2 2 8 wx t 31006 100619 T 1006 0x T External energy due to distributed applied force fx 2 shear forces Q12 and moments M02 applied at both ends L 8W 2 fx 1 8wx 1 dx 0 Q18w0Z l Q2 8wL t M18w0tx l M2 8wLtx 3 Q1 Q0 1 Q2 QL 1 M1M0l M2MLl Potential energy of four unknown springs Vs 2kw1w0t2 k91w0z 4 sz UJL 02 l 92 wLa 0 Hamilton s principle l2 8T 8Vb 8Vs8Wdt0 5 I1 From here on FEM modeling differs from continuum modeling Key departure in FEM modeling of beams from continuum modeling 1 Instead of carrying out the variation the 8pr0cess 0f Hamilton s principle5 in terms of the continuum variable wx t the energy expressions T V W etc are approximated on a completely free beam with an arbitrary length E area A the bending rigidity E I and mass density 0 2 The interpolation of the transverse displacement wx 2 over the completely free beam segment 0 f x f E is chosen to satisfy the homogeneous differential beam equation ie Elwbcanxxxx 0 if possible If not it is chosen such that as the length of the beam gets smaller and smaller it approximately satis es the homogeneous equation We now consider a oompletely free isolated small segment of a beam A completely free element of length taken out from the beam A l 4 W 1 l gtx AW an el ment isolated Z kel kez s a gt X 1 2 kwl EI mX km W lt L Approximation of wx t over the element segment t ltlt L wx z eot 61tx 62tx2 C3 tx3 7 Observe that the above interpolation satis es the homogeneous a i er ential equation of beam 6 How does one determine the eoe ieients co cl 62 63 I The answer comes from the nite element method The FEM beam bending element 1 Observing that wx t and wx tx are linearly independent from the result of variational formulation we specify them at nodes 1 and 2 in the previous gure Namely two discrete variables at node 1 wx1 t and wx1 tx 8 two discrete variables at node 2 wx2 t and wx2 tx 2 For convenience we introduce a local coordinate system or ele mental coordinate system x y and set x1 0 and x2 E Substituting 8 into the interpolation function7 we obtain w1l wx1 o z col 01ZO c2z02 c3l03 610 wx1 0 ox c1r 2c2r0 3c3r02 9 w2r wX2 E I com 010 czw2 03053 920 wx2 3 0x 610 02l2 c3 1332 Solving for the coef cients and backsubstituting into 7 one obtains W 0 H1lt5gtw1ltrgt Hm 61m 5 we Hm wzm Hm 62m 2 3 2 3 10 H1lt5gtlt1 35 2 gt H25 2 w 2935 352 253 Hm 52 53 For notational clarity we express 10 in the form wxl N 11 N H1 H2 H4 11 qltrgtTltw1rgt 61m w2l 62w l W1 91 3 A completely free beam element and its nodal degrees of freedom Element Discretization Elemental kinetic energy 3 Tel 2 mx wx Dick 0 1 12 O me wlt5rgt d5 x a Elemental potential energy K V61 EIx wx Dix dx 0 13 1 E f0 E ngm gd Discretizatl39on of Elemental Kinetic Energy Z 1 2 T6 iO 1725 we twig 1 Z A Mi qu NT N 6100 5 14 iqT 01 pew NT N d5 3910 iqu mel 610 where we utilized 11 viz 3 we r Zcimx MS 10 15gt i0 Diseretz39zatl39on of Elemental Strain Potential Energy 1 E1 V51 1 f0 7 we 0 d5 1 E1 1 f0 3 qltrgt1 N515 Nss qltrgtd5 16 1 E1 1qltrgt1 0 73 39Nss d5 qltrgt 1qltrgt1 ke1qltrgt where we utilized we 0 Neg qltrgt 17gt Elemental mass matrix me mpA Elemental sti ness matrix k 39 12 EI 65 k 3 E 12 6t 6t 4E2 6 2E2 12 6E 6E 18 6E 2E2 1 6E I 4E2 J 19 Approximation of elemental external energy 3 Z 1 swelz fxt6wxtdx fgzNg3qzzdg 0 0 Q16w0t Q2 6112LtM1 6w0tx M2 8wL 0x Let s utilize the approximation10 woe t NS qm 1 woe 0x ZNlt5g qm and observe that so that we have Z T 6ng 6qlttgtT N foe 1 ads Q1 66111 Q2 66130 M1 66120 M2 66141 0 Z 6qlttgtT 1 21 fgl NTfltxzgted5ltQ1 M1 Q2 M2gtT 0 3 20 21 22 Summary of elemental energy expressions T ma ew me web V51 ma ew 181 an 23 8W 6qltrgtlte T tel Question How do we model a long beam with beam elements Answer 1 Partition the long beam into small elements This needs to be explained 2 Generate the elemental energy expression This is done see 23 3 Sum up the elemental energy to form the total system energy 4 Perform the variation of the total system energy to obtain the equations of motion Let s now work on items I 3 and 4 A Partition a beam into many elements Partition a beam into finite elements element 1 I 2 I 3 l 4 I l global node 1 392 l3 lementl element 2 element 3 global node i 2 2 3 4 elelental node 1 2 l 2 l 2 element 1 element 2 element 3 u u u u u 1 U2 2 3 3 4 l l 2 2 3 3 ql qz ql C12 C11 C12 Partitioning and establishing the Boolean relation between the global degrees of freedom u and elemental degrees of freedom q 3 Establish the relation between the elemental and assembled global degrees of freedom 1 1 1 111 For element 1 q qugt uz 2 2 qi 112 For element 2 q 12 113 3 2 q ua For element 3 q q m u 10 0 qm 01001111 100 uz 12 0 loolol 113 I q L8 8 0 H 114 24 For general case gt qt L g 3A Sum up the total system kinetic energy 71 Ttotal Z Tel el1 qlt1gtTmlt q1gt qlt1gtTmlt q1gt qltngtTmlt gtqltngt 11 T m 0 0 11 T 19 0 mlt2gt 0 1km 0 0 0 qm L 0 0 mltngtJ TtotaZ t m it 25 107 3B Sum up the total system strain energy 11 total e V Z Vb el1 qlt1gtTklt qlt1gt qlt1gtTklt qlt1gt qngtTklt gtqltngt q1 T k 0 0 q1 T q2 0 Id 0 q2 OJ qn 0 0 L0 0 km Vtotal qu kqg 26 qn 3C Assemble the partitioned elements back into a beam this corresponds to substituting the partitioned elemental degrees of freedom qg by the assembled or global degrees of freedom ug via 24 c Total assembled system kinetic energy 27 c Total assembled system external work 8Wt0tal el1 593T fl 2 an Lng 28 Z T WW 2 8ug fg EulerLagrange 5 equations of motion I Lagrangian Ttotal Vtotal 1 1 29gt EugMg g ZugKg g From the generic equation of motion a M M 30 dt 8k aqk Qk we obtain the FEM equation of motion from 29 and 30 MgtigKgugfg 31 ASEN 5022 Spring 2005 Dynamics of Aerospace Structures Lecture 07 01 February 2005 Hamilton s Principle EulerLagrange s Equations and the Method of Lagrange s Multipliers Let s begin With d Alembert s principle N ZltFi mii z39 39 51 239 0 i1 1 N d 1 5T 5W 5 39 0 261707111 r1 N 5W 2 ZE 5r 5T 2 dmi39 1 i1 Let s integrate 1 in time I2 N d 5T l 5W m1quot8r dt0 2 if Edi 1 ltgt 1 t2 t2 N d 5T 5Wd1 m1 51 139 d1 3 1 1 1 I1 139 Noting that we have I2 N d Z ml39i l39 5I l d1 H 121 dt N Z mii z39 1796 11 4 so that if n H and I l3912 are speci ed we have 51711 51712 0 5 Hence equation 2 becomes f2 f1 5T5Wdt0 which is known as extended Hamilton 5 principle 6 Hamilton s Principle for Conservative Systems In general the work done 8 W consists of two parts 5 W 5 Wcons 5 Wnoncons 7 Where subscripts cons and noneons designate conservative and noneonservative systems respectively Hamilton s Principle for Conservative Systems For conservative systems we have have from 6 With 5VVnoncons O I2 5T 5VVcons dt O 8 f1 which is known as Hamilton s principle for conservative systems Action Integral for Conservative Systems Observe that the work done on the conservative systems can be expressed in terms of the corresponding potential energy Substituting 9 into 8 one obtains f2 5T 6Vdr0 f1 1 f2 5T 6Vdr0 t1 1 f2 t1 51 5Ldt0 LzT V 10 EulerLagrange Equation 5 of Motion Let s substitute the general work expression 7 into the extended Hamilton s principle 6 to obtain f2 5T 5Wcons 5Wn0ncons d1 O 11 t1 which With 9 together With L T V becomes f2 5L 5VTnoncons dt O 12 f1 Since 8L can be expressed in terms of the generalized coordinates as quot aL aL 5L Z 5qk 5qk 13gt k1 aqk aqk one obtains t2 t2 quot M M 5Ldrf Z aqk 5qkdr 14 n n k1 aqk aqk Integrating by part the rst term of the preceding inte gral we obtain t2 quot 8L quot 8L t 5Qk d1 5Qk 2 A k2 861k Z 8 t1 t2 quot d M Zd 5qkdt 1 1 k1 tan 15 Since 8qlt1 2 5g t2 O we have 2 quot aL Z aqkdr 1 1 k1 an t2 quot d M Z d lt gt5qkdr 1 1 k1 t an Substituting 16 into 14 then introducing the resulting expression into 12 we nally obtain t2 quot d M M Z lt gta Qkaqkdr 1 1 k 1 Qk dt aqk n 17 5VVnoncons Z Qk5Qk k1 Since 5qk are arbitrary we obtain d8L 8LQ 18 dt aqk qu k which is called EulerLagrange s equations of motion A SpringMaSSBar System for EulerLagrange 5 Equations Notice how simple the problem description becomes now Kinematics Position vectors of mass M bar at C and B Where the nonconservative force F is applied I A 2 xi L rc xi l sin6i cosej 2 lt19 rBx1Lsm61 0056 L r 2 C 2 Kinetic energy T mi A 1240 Mi C to 1092 o 1 20 M m c2 MLa ce cose ML262 Potential energy 0 er Vgravz39ty ng dl C 21 I39C L Mg 1 cos 6 Total potential energy V Vsprl39ng ngity 22 N oneonservative work 5 Wnoncons 5Vknoncons 2 FB 39 5139B Fi 5xi Lcos 6i sin 6j56 F5x L 0056 56 Qx F Q9 2 FLcose 23 Note that there are only two state variables x and 6 in the system kinetic energy potential energy and the noneonser vative work EulerLagrange 5 Equations daL M Q F dt ax a x 24 d M a L FLcos dt86 89 Q9 EulerLagrange 5 Equations d E M I m c ML6 0036 kx F ML d 25 ML39CQ g sine 2 FL 0056 22 Modeling of Shock Absorbers Chapter 22 MODELING OF SHOCK ABSORBERS 22 2 221 INTRODUCTION We have studied the socalled Den HartogOrmondroyd 2DOF shock isolation model Lecture 02 and discrete modeling of cable models Lectures l9 and 20 In this lecture we reVisit the shock isolation problem in more detail which illustrates the importance of discrete modeling in the age of FEM world 222 INVARIANT POINTS OF THE DEN HARTOG ORMONDROYD MODEL Consider a TwoDOF system as shown Figure 221 The governing equations for the system is given by M1561 K1x1k2 x1 x262 x1 x2 f1t 221 quot12 962 162 962 x1 02 x2 x1 f2t xv W N W Two DOF SpringMassDamper Suspension Model Figure 221 2 DOF Shock Isolation Model Assume that the excitation acts directly to M1 in the form f1t 2 F1 sinat f2t 0 222 The solution of 22 1 assumes x1t ImXlequot x20 2 ImXzej 223 Substituting these into 221 yields the coupled equations in the frequency domain as M1a2 K1k2icwX1 k2 F1 224 k2 J39CCZX1 mzw2 k2 J39CwX2 0 22 2 22 3 223 SHOCK ABSORBERS X1 obtained from the above equation reads X1a m2w2 k2 jcw 22 5 F2 lez K1 mzw2 k2 m260219 J39wclKl M1 m2w2l Let us parameterize the model as 91VK1M17 a2kzmzy C2 m2917Mm2M17 226 sav9h mmefwM9hccmm91 The ratio of X1 w to x is the ampli cation factor from the static to the dynamic response which can be expressed as ma Q VW 1 x 2LS2321 22 22 21 2 2 22 n w mm o gto fH Frequency Response Functions for Different Damping Ratios ofa 2DOF Suspension Model 101 7 Invariant points i E lt 103 10391 101 Frequency Hz Figure 222 Frequency Response Function at Mass 1 Figure 222 plots X1 0 vs the driving or exciting frequency a for various damping coef cients Xx 1e by vary1ng 5 while keeplng other parameters constants Note that there are two pomts marked by red dots in the gures the curves pass through for all damping ratios We will call these two points as den Haltog s invariant points and will exploit its properties in two ways one is to minimize the maximum amplitude resulting in good shock absorbers and the other is to maximize its peak amplitude so that the resulting resonators possesses a minimum loss or highQ performance 22 3 Chapter 22 MODELING OF SHOCK ABSORBERS 22 4 223 SHOCK ABSORBERS A good shock absorber should have the following characteristics a The peak amplitude should be insensitive over a wide range of excitation frequencies b The cost for implementing damping should be affordable c The auxiliary mass m2 should be as light as possible d The response time duration to reach its steady state should be short We now proceeds with the task of designing a good shock absorber 2231 Determination of Undamped Natural Frequencies the natural frequencies are determined by setting 0 0 from the denominator of 227 MfzSZ 32 132 f2l 0 1 228 s4 11Mf2s2f2 0 The two distinct roots of this equation will be designated as 312quot and 3quot 2232 Determination of Invariant Points Since the two invariant points would play a key role it is necessary to nd the two coordinates Of several approaches we utilize the fact that they are independent of system damping c To this end we rearrange the frequency response function or the magni cation factor 228 as X1a 257s2AB 1 x MHSycyfy 2cs2 CwlzyswQi A21 229 BS2f22 39 C s2 1MSZ2 D IL 232 32 1S2 122 The dampingindependency of the two invariant points implies that at the invariant points the magni cation factor satis es 1 ADBC DIG A B 2210 22 4 22 5 223 SHOCK ABSORBERS for which the ampli cation factor becomes X1w xx A Hscm E The equation for determining the invariant points 2210 yields two equations l 32 1 1132 52 f2 MfZSZ 32 152 f2l l s2 1 1132 s2 f2 2 WW s2 1W FM The rst of the above relation yields the trivial solution viz s w 521 0 a static case which ifs not of interest for the present discussion The second of 2213 can be rea1ranged as 2 ms 21 f2 w 2s2 2f2 0 The two invariant points 31 and s2 can be shown to satisfy 2 2 s1lt lt32 1 M and the corresponding ampli cation factors are given by X1a 1 xx 11 03313 2233 Determination of Ampli cation Factor at the Invariant Points The foregoing analysis now supplies the coordinates of the two invariant points as l P 951071 317 7 Q x27y2 S27 I1 1 ms In order to see succinctly the magni cation factor we arrange 2214 as 2 2 M Zz2MZ107 2 2M 2M 321M 1 so that the product of the two roots of the above equation 21 gtlt 22 is given by 22 5 11MS 2211 2212 2213 2214 2215 2216 2217 2218 Chapter 22 MODELING OF SHOCK ABSORBERS 22 6 2 2 2122 M 1 2219 M M which is independent of the frequency ratio f wz 21 It should be noted that for a shock absorber the magni cation factor should remain small for a wide range of frequencies which can be realized if we set the ampli cation factors to attain their maximum at the two invariant points and at the same time their peaks to be the same The latter e condition can be realized if the coe icient associated with zterms in 2218 vanishes l 1 1M2f20 f um2Ml 2220 521 1M The ampli cation factor X 1 cu x for this particular choice at the two invariant points becomes X 1 2 1a 2 1 2221 xx l11 MS1 2l M 2234 Determination of Absorber Damping In the previous section we imposed the condition on the ampli cation factor at the two invariant points to be the same However it is possible that the ampli cation factor may be larger at other exciting frequencies Of several possible choices we impose that the maximum ampli cation factors for the frequency ranges of interest would not exceed that at the invariant points This can be carried out as follows First we differentiate 229 with respect to s w 521 set the resulting expression to zero a X 6 160 0 2222 S xx for the two invariant points obtained from 2214 3 L1 2 Fors 31 52212 u2 Cc 81 m3 1 2223 3 L5 Fors 32 52212 zu Cc 81 m3 In practice we must have one damping value for all the frequency range To this end we take the average value to be 3 g 2224 81 W3 22 6 22 7 223 SHOCK ABSORBERS 2235 Applications First in order to utilize Matlab capability we express 225 in its transfer function form X10 quot1232 k2 CS 22 25 F2 M132K1mz32k2mz32k2SCK1M1mz32 where s is now the Laplace Transform variable The preceding equation can be arranged to read X10 22 61052 6115 612 x 1 b0s4 m3 13232 1333 b4 2226 2 002176112255217 a2602 b0 1 b121us21 b2s21uw b3 2552 Proof of two invariant points for den Hartog oscillator i 9110 0295 noo1 o invariant points Amplitude Ratio X1wx5t 52 s1 Frequency 0 Hz 101 Figure 223 Invariant Points of a Den Hartog Oscillator The ampli cation factor 2215 vs the driving frequency cu is illustrated in Fig 223 with 21 10Hz 602 95Hz u 001 Observe the two invariant points at which the ampli cation factors are di erent namely at 31 and 32 points determined by 2214 and its magnitudes given by 2217 22 7 Chapter 22 MODELING OF SHOCK ABSORBERS 22 8 Optimized Amplitude for den Hartog oscillator l 01206 g 003015 Qopt 00603 on i Amplitude Ratio X1wx5t 91110 0295 uoo1 l Frequency 0 Hz 101 Figure 224 Den Hartog suspension model optimized for a wider range of frequencies When the ampli cation factor is optimized to have its peak at the invariant points by selecting the frequency ratio according to 2220 and the damping ratio according to 2224 zeta 00603 its ampli cation attains its maximum at the invariant points This is illustrated along with two nonoptimal damping coef cients in Fig 224 Observe that the optimum ampli cation factor at the invariant points are about 135 which may not be acceptable Figure 225 illustrates the reduction of the peak ampli cation factor by about half by increasing the small mass from m2 001M1 to ve times m2 005M1 In practice a suspension system is designed to operate far smaller then the invariant points or su iciently larger then the invariant point frequencies Under this scenario the infrequent disturbances close to the invariant points can still be bounded 22 8 22 9 223 SHOCK ABSORBERS Optimized Amplitude for den Hartog oscillator on Amplitude Ratio X1w5t Frequency 0 Hz Figure 225 Den Hartog suspension model optimized for a wider range of frequencies Two DOF Vehicle Shock Isolation Model Ground input Xg t Figure 226 A Two DOF Vehicle Suspension Model The twoDOF suspension model that we have been studying so far is applicable to machinery shock isolation design For vehicles subject to road roughness conditions a modi ed model is more appropriate as shown in Fig 226 An analysis of this model for achieving an optimum shock isolation is left for an exercise 22 9 Lecture 15 09 March 2004 Discrete Modeling of Dynamical Systems 1 Introduction With a plethora of engineering analysis software packages available the dynam ics specialist can perform a variety of computerized modeling and simulation of complex systems It is not uncommon these days that a typical structural dynamics model may consist of several thousands to several millions of degrees of freedom often employing the nite element method as the dominant modeling tool As a result what used to be a challenge to the dynamics specialist several decades ago viz more elaborate models that result in a large number of degrees of freedom has become routine practices At present the most timeconsuming part of engineering analysis including dy namical systems is model development tasks While still signi cant the core of computational effort has become relatively insigni cant thanks to the advances in computational power known in the computer world as Gordon Moor s law To the engineer a new challenge has emerged how to interpret a vast array of analysis results in terms of plots tables and other graphical data representations A rst step to a meaningful interpretation of computer generated analysis results is to cul tivate an ability to reduce the complex systems to a set of simpli ed models from which one can relate the observed phenomena to the complex models This lecture is aimed at introducing to the students how to construct simple models how to interpret the simple model results and how to relate the results obtained by simple models to complex simulation models 2 OneDOF Modeling of Cable via the Finite Element Method Suppose we do not know how to obtain a sensible simple discrete model and we are asked to construct a single two or at most three degrees of freedom cable model by the nite element method A quick reference to a standard nite element text provides the following twonoded elemental mass and stiffness matrices ep 21 m 12 kel 1 where T is the tension of the cable N K is the elemental length not to be confused with the total cable length m 0 is the cable mass per unit length kg In As shown in the beam nite element modeling one can construct a single two and three degrees of freedom model as follows OneDOF Model A nite elementbased construction of simple models can be constructed in two ways I From one linear element model by constraining either of the two nodal point This approximation is shown in Figure 1a Since node 1 is xed the resulting oneDOF model can be obtained from the elemental mass and stiffness matrices given by l as L T pTc39izm 2qu N 2 which yields its frequency to be 1 a L 1217321 1 3 L 0 L 0 Comparing this with the analytical solution 0 2 we nd a frequency error of about 45 by this one degree of freedom simple model Z L one linear a one dof model element model 1 2 from one element 2 g Lg Z L2 tWO lineal 39 I b one dof model element mOdel 1 2 3 from two elements 5 g L2 g L2 fixed boundaries c onedof model from one quadratic element one quadratic element model Figure 1 Finite element modeling of a cable 2 From two linear element model by constraining the two ends This approxima tion is shown in Figure 1b First the freefree twoelement model can be written as pg 2 1 0 5111 T 1 1 0 6111 f1t KP 4 1 61392t 1 2 1 q2tf2t 4 0 12 5131 0 1 1 6131 f3t Constraining ql t and 613 t from the preceding equation which is equivalent to eliminating the rst and third rows and columns yields 405 2T 5I21 76121 f2t7 K L2 5 which yields its frequency to be 12 T 34641 T a M L L 6 p p We thus conclude that a niteelement based one degree of freedom approximation yields about 10 frequency error 3 From one quadratic element model by constraining the two ends This approx imation is shown in Figure 1c First the freefree one quadratic element model can be written as M 4 2 1 5111 2T 2 2 1 6110 f1t 2 16 25l2li 2 8 2 q2tf2t 1 2 4 6150 1 2 2 6130 f3t7 Constraining wl t and 1123 t from the preceding equation which is equivalent to eliminating the rst and third rows and columns yields 16 L l6T 3 6220 Tana f2l 8 which yields its frequency to be E a 3L 10 12316227 1 9 L p L p We thus conclude that a quaadraticelement based one degree of freedom approx imation yields about less than one percent in frequency error Remark It is worth noting that one quadratic element model captures the funda mental frequency of cable with less than 10 error whereas one and two linear element models yield about 45 and 10 in frequency error respectively Note from Figure 1 that cable modeling by linear elements in general causes the cable slope discontinuous which is a poor approximation of continuous cables On the other hand modeling by quadratic elements satis es the interelement slope con tinuity It is this interelement slope continuity that plays a key role for simple cable modeling Of course if one employs a large number of linear elements linear element models will capture the analytical frequencies This is illustrated in Figure 2 102 Cable model by linear fem method Frequency error of rst cable mode oflinear cable element Figure 2 Assumed Modebased simple models of a cable The error analysis of linear elements for cable Vibration analysis indicates that about 40 linear elements are needed to capture the rst cable mode for three digit accuracy and about 200 elements for ve digit accuracy Also note that the frequency error for models of in excess of 150 elements decreases linearly a classical convergence rate result in the nite element method 2 Modeling of Cable via Assumed Mode Approximations In a broadest sense the nite element method is a polynomialbased assumed mode approximation For example a linear element assumes 2 w l1 611l1 612l 1E EL 7x 10 where the element coordinate origin is located at the center of the element Note that the cable slope for the linear element is given by wxtx Em I 2 M which suggests that the cable slope is discontinuous constant The only exception is when the cable forms a straight line a triVial state On the other hand the quadratic element approximates the displacement by w 1 1 611l1 2612l 1 613l 11 Tqm one dof model by an assumed mode I Wxt qt sin1xL I gt X Figure 3 Assumed mode based modeling of a cable Note that the cable slope for the quadratic element wx tx x I is continuous A classical oneDOF massspring modeling is to introduce the following form wx 1 Hx61l 12 where q t is the displacement where the mass is located and H x is the assumed deformation shape Thus for oneDOF approximation of a cable with two ends xed the following displacement is assumed 71x wxl SIM 611 13 If one needs more than one displacement approximation then depending on the analyst s experience and the nature of problem to be modeled wx I may be modeled as woe t sing 6111 sm2 7xgt 612m 14 Using the assumed displacement given by 13 we have 1 2 1 77x 2 2 T 0 5mm xtgtdx f0 pxsm7gt q ltde p5 57 6120 l l V Txw xtdx Txcos717x2q2tdx 15 0 0 T712 57 6121 8 Wnanwns 0 due to free Vibration modeling Th discrete EulerLagrange s equation from the above energy expressions yields the following oneDOF equation 16 Hence we observe that the assumed displacement l 3 yields the exact fundamental frequency of a cable whose ends are xed This si expected because the shape function sinnx is the exact mode shape 3 Finite Element Modeling of Cable with Flexible lVIiddle Support With the preceding background let s revisit the task of modeling of a cable with two xed end and a middle exible support as shown below a Two linear element model of cable L2 L2 b Four linear element modeling of cable L2 L2 Figure 4 A simple equivalent model of a cable with a middle support A Single Degree of Freedom M ode If one utilizes two linear elements and assem ble the middle massspring model the resulting equation before applying boundary condition is given as see equation2 as shown in Figure 4a pg 2 1 0 911 T 1 1 0 wt 0 F 1 46um 1 921 2 1 2Km 1 q2t 0 17 0 1 2 gm 0 1 1 q3t 0 Constraining out the two support conditions Viz q1 613 O we obtain p5 T gel 6Mm 6121 E2 Km 12 Z O ng 18 K LZa Hm ZMmP a Km 2 T As in the homework problem if the required frequency is to be can 2 g with pm 2 10 we obtain K m K K m L Km 2T 725275 19 This incurs about 30 error when compared with the analytical solution A Three Degrees of Freedom Model Let us utilize four linear elements and assemble the middle massspring model as shown in Figure 4b Following a similar assembling and constraining processes taken for the case of a single degree of freedom case one arrives at the following three degrees of freedom model equation pg 4 l 0 t T 2 l 0 q1t 0 F1 46Mm 1 2t2 l 2Km 1q2t0 0 l 4 Q39 0 l 2 q3t 0 20 1 L4 um Mmltpegt Km m whose eigenvalue problem after rearrangements can be eXpressed as 2 l O 302 4 l O det l 2Km l 96 I 46Mm lO 21 O l 2 O l 4 For 6L 21571 and pm 2 2 as Mmp 272 withK L4 we nd Km KmL KmL szlcmz z Q 45 gt 4T 2T 90 22 which yields about 12 error compared with the analytical solution The important nding is that as the number of elements increased the solution converges 4 Classical Lumped MassSpring Modeling of Cable with Flexible lVIiddle Support A typical classical approach not always is to derive mass and stiffness matriX for a single degree of freedom case and try to make the most out of that simple model To this end let us utilize the single degree of freedom model derived by assumed displacement as given by 16 Two cableswith L2 l4 U2 U2 i Model parameters used t 1 I Vl1 qllslmuU I MxltllslnvrrU I Z L 2 0 0 ol m 7 Modeling a cable of length L I from two cables of U2 I T z k 2E Mm Vimpi Massless node I sz Left cable Right cable Figure 5 Constructing a cable of length L from two cables of length L2 In so doing there are three strategies distributing the mass as well as stiffness and distributing only the stiffness and distributing both the mass and the stiffness In the present description only the stiffness is redistributed The approach is general enough that one can readily extend the modeling process to other types of modeling This is illustrated in Figure5 Note that if the end purpose is simply to model the cable without considering the exible middle support one may have to redistribute the mass as well There are several considerations that one must evaluate before converging onto a sensible model In the present case it is required that the frequency of the assembled longer cable should behalf of the individual cables of length L 2 This can be determined by the following two degreeof freedom eigenvalue problem which can be constructed from Figure 5 as k wzm O l ak x1 O k wzm l akx20 23 l ak l ak 21 ak xm where xm is the displacement at the massless point for now from Figure 5 The above frequency equation gives the following characteristic equation 21 ak k wzm ak wzm 0 24 The above equation gives two frequencies k T k T wmax Z 0min m K p m e p Observe that wmax is the same as the halflength cable Hence it is the lowest frequency omn has to be half of This determines the stiffness distribution factor a to be a 26 Model parameters used Cl Cl 7 E l L2 U N H U a E M 4 c m 7 E Tn392 k 21 M MmM Kml Km Figure 6 A simple equivalent model of a cable with a middle support By adhering to the stiffness redistribution approach the exible support Km is redistributed as shown in Figure 6 The equation of motion for the system shown in Figure 6 can be derived as kbK 3kbK mxl 2 mxl 2 mlxm fl 7 3k 7 m 2kEKmx2 T5Kmlxm f2 27 3k 3k 7 mem Km 7xm x1 x2 fm whose characteristic equation is given by k Km m wz 0 T 191 det 0 k Km m 602 d gKm 0 1 gm 1 191 Km Mm 602 28 For computational expediency divide each of the three rows of the above charac teristic matrix by m and using can 2 xkm we have 1 gcmm 02 0 Km1wg det 0 1 gxmmg 02 Km1wg 0 Kmlw721 Kmlw721 Km ZHm 502 29 where Km and um are de ned as Km 2 Mm 7 Km T 7T2 Hm 0K For this characteristic equation by varying b Km was computed The result is that the best value of Km is obtained when b 0 When b O is substituted with a 1571 and pm 2 1 Km is determined to be Km 2 2 2 gt KM 10 856 31 K m T n2 T which when compared with the exact solution given by Km T 10264 yields about 5 error There are several additional aspects in modeling of simple dynamical systems which include the use of lumped mass matrix diagonal mass matrix as opposed to coupled mass matrix redistribution of kinetic energy or mass distribution tran sient response analysis reduction of largeorder models to smallerorder models incorporation of experimental observations into nite element models nonlinear phenomena among others Some of these aspects not all will be discussed in the remainder of the course 5 The Method of Lagrange s Multipliers Chapter 5 THE METHOD OF LAGRANGE S MULTIPLIERS This lecture covers 0 Begin with the Lagrangian L T V 0 Identify the constraints in the system 0 Augment the constraints to the Lagrangain by multiplying with an unknown for each constraint Apply the formalisim of the EulerLagrange equations of motion by treating the constraint functional as part of the system Lagrangian 0 Demonstrate the method via example problems 5 3 The EulerLagrange s equations of motion are obtained from the system energy that consists of the kinetic energy potential energy and virtual work due to nonconservative forces In other words while enabling the derivation process simpler the EulerLagrange formalism eliminated an important information the reaction forces on joints and constraints such as boundaries Often both the designers and analysts need to know joint force levels so that necessary joint articulations can be carried without jeopardizing the system safety A case in point is the human joints which repeadedly over loaded can casue bone deformations and arthritis The question is How do we reinlroduce the constraint forces reaction forces within the Euler Lagrange formalism This was developed by Lagrange and described in his book M canique Analytique published in 1788 This is accomplished as follows Let s consider holonomic cases viz the constraint conditions that can be eXpicitly stated in terms of position vectors and consider a twodof springmass system shown below Suppose we would like to know the reaction force between mass m1 and spring kg and the reaction force between spring k1 and the attached boundary We now proceed with the procedure that leads to the determaination of the reaction forces 5 5 A procedure for the method of Lagrange s multipliers Step I Step 2 Step 3 Step 4 Step 5 Partition the system into completely ee subsystems Identify the conditions of constraints between the completely free systems Construct the energy of each of the completely free subsystems Obtain the total energy by summing the energy of each of the completely free subsystems Append the conditions of constraints by multiplying each With an unknown coef cient A multiplier to the total system energy kinematic and potential Chapter 5 THE METHOD OF LAGRANGE S MULTIPLIERS Reference solution based on the assembled system LzT V 2 1 2 1 T Emlxl fmzxz 51 V k1x 2062 x12 SW f1lSX1 f1l8xz Hence we obtain via the EulerLagrange formalism M1551k1 k2X1 kzxz 1 10 52 M2552 kzxz k2X1 f2l 5 7 Application of the method of Lagrange s multipliers Step 1 Done in the preceding gure Step 2 Conditions of constraints Between k1 and the xed end nog x0 xg 0 53 Betweenmlandkzz n132x1 X320 Step 3 Energy of two completely free subsystems For subsystem 1 T1 m 156 V1 klm x02 54 5W1 f1l5x1 Chapter 5 THE METHOD OF LAGRANGE S MULTIPLIERS 5 8 For subsystem 2 1 2 T2 m2x2 V2 2062 x32 5395 5W2 f2l3X2 Step 4 Sum the total system energy TT1T2 VV1V2 56 5W5W15W2 5 9 Step 5 Construct the constraint functional 7395 nogkog 7713 A13 57 x0 Xg log X1 X3 13 Total Lagrangian of the partitioned system L T T V V 1 2 1 2 58 7TXgxoxix3 0g 13 Chapter 5 THE METHOD OF LAGRANGE S MULTIPLIERS Partitioned equations of motion 8x1 term 8x2 term 8X3 term 8x0 term Slog term 8M3 term 5 10 quot11551 k1x1 k1X0 k13 f1l quot12552 kzxz k2x3 f2l k k 0 2X2 2x3 13 59 k1X1k1X0 10g 0 x0 xg 0 x1 x3 0 5 11 Matrix form of the partitioned equations of motion miD2 k1 0 k1 0 1 X1 f1 0 szZ k2 k2 0 0 0 x2 f2 0 k2 kg 0 I0 1 X3 0 k1 0 0 k1 H 0 x0 0 510 0 0 0 1 0 0 Mg xg i 1 0 1 0 0 0 i mi i 0 i where D2 5 52 Symbolically the above equation can be expressed as g 511 Chapter 5 THE METHOD OF LAGRANGE S MULTIPLIERS 5 12 Hence solution of 511 provides both the displacement x and the reaction forces A Question Can one obtain the equations of motion for the assembled system 52 The answer is yes To this end we observe the following relation between the partitioned degrees of freedom and the assembled ones as x 1 l 0 0 x2010 x3 100 lxol l001i Substituting the above assembling operator L into 511 we obtain X1 x2 gt XZLXg 512 xo 5 13 LT AL LTC Xg LT f CT L 0 A xb Simple matrix multiplications show 0 0 LTC0 0 10 quot11D2 k1 k2 k2 k1 KzLTAL k2 m2D2k2 0 k1 0 k1 Therefore 513 reduces to 513 514 Chapter 5 THE METHOD OF LAGRANGE S MULTIPLIERS 5 14 m1D2k1k2 k2 k1 0 X1 f1 k2 szz kg 0 I0 x2 f2 k1 0 k1 l x0 0 515 l 0 0 1 ml liogl lxgl Finally invoking the bounday condition x0 xg0 516 we arrive at quot11D2 k1k2 k2 x1 2 fi 517 k2 szZ k2 X2 f2 Note that solution of 516 provides the reaction force at the leftend support whereas 517 does not Chapter 5 THE METHOD OF LAGRANGE S MULTIPLIERS Kinematics Kinematics of Subsystem 1 FA Xi yj VA Xi J j Kinematics of Subsystem 2 1 C XCi yCj L rA I C 3sm6i COSQj L 132 390 3sm6i cos6j 5 16 518 519 5 17 Step 1 Partitioning is done in the preceding gure Step 2 Conditions of constraints 1 Between mass m and the ground y yg j 0 2 Joint at A FA FA 2 0 U L x xC s1n61 L y yc Ecos9j0 Chapter 5 THE METHOD OF LAGRANGE S MULTIPLIERS Step 3 Construct the energy of completely free two subsystems Energy of Subsystem 1 T1 m x2 ye V1 mgy kx2 Energy of Subsystem 2 T2 Mx y y39 1C92 V2 M8 yC 5 18 521 522 5 19 Nonconservalive work 8 Wnonwns 8Wn0ncons FB 39 8FB L Fi 8xCi 8yCj cosQi sm6j86 2 L 523 F8xc E 0056 86 ng F Q9 FL 0056 Chapter 5 THE METHOD OF LAGRANGE S MULTIPLIERS Step 4 The total energy of the system T mx2 y39Z Wag ya ace V mgy kx2 Mg yC Step 5 Construct the constraint functional L IT 2 y yg og x xc sm6gtkAx L y ya 3cos6tAy 5 20 524 525 5 21 Partitioned equations of motion via the method of Lagrange s multipliers 8x mi kxAAx20 8y mjimgkogAAy20 8x0 Mica ch F 5yCZ MjiCMg AAy0 L L 86 IC63c056AAxESIn6AAyFLc056 526 Slog y yg0 L 8 sz x xC sin60 L 8AAy y yC 30089 0 Note that log AAx AAy corrrespond to the reaction forces NA XAC YAC that are obtained Via Newton s Second Law ASEN 5022 Spring 2006 Dynamics of Aerospace Structures Lecture 02 19 January 2006 Vibration 0f TwoDOF Systems 1 Equations of Motion Consider a twoDOF model as shown in Fig 1 The freebody diagrams for masses M1 and m2 are also shown in Fig 1 Note that the inertia forces for both masses are associated with minus sign for the inertia forces can be considered as resisting forces cf 2 f ma 2 O Summing the forces acting on each mass the equations of motion for the coupled twomassspringdamper system can be written as Ford1112f f11 M1551k2x2 x1029 2 X1 K1961 0 For 7772 3 Zf f2t m2552 k2X2 x1029 2 X1 0 U 1 M1551 f1t K1x1k2x2 x102 X2 X1 7772 552 k2 x2 X1 02 X2 X1 Fig 1 Two DOF Simpli ed Vertical Motion Model Two DOF SpringMass Damper Model ii Free Body Diagram for Two DOF SpringMassDamper Model As the preceding equation involves two displacements x1 and xz its general solution involves complex matrix differential algebra For design considerations however important insight can be gained by considering the special case of forcing function given by f21 F2616 f1t Fiejwl so that the solution assumes the form of x1t Xlejw x2t Xzejwt Substituting 2 and 3 into 1 we obtain 602M1X1 Fi K1X1k2 Xz X1J002 Xz Xi w2m2 X2 F2 k2 X2 X1 16002 X2 X1 In order to solve for X1 and X 2 let s rearrange the preceding equation to read w2M1 1002 K1 k2 X1 Fri162 1002 X2 w2m2 16002 k2 X2 F2 k2 160021 2 3 4 5 Solving for X1 and X 2 we obtain H2CUF1 H12CUF2 1710 1720 H12260 l 2 H1 0172 H12CUF1 H1wH2wgt H122 w 6 2 1710 602M1 1002 K1 k2 1720 602m2 16002 k2 17120 k2 1002 Although 6 appears to be very complex several simpli cations are possible to aid vibration designers This is studied below 2 What Are Vibration Modes and Mode Shape It turns out that the motions of X1 and X 2 given by 6 are not randomly independent as the solution may suggest They are interlinked by the property called mode shapes Understanding the physical properties of mode shapes is important if one is tasked to design structures subject to Vibrations As a motivation let us consider the following special 2DOF case mlzmzzl k1k226187r2020f1t0 7 For this model we apply two different excitations 02 X sin 098 gtllt 7U f2t 8 sin 098 gtllt 26187rt Figure 2 illustrates the two responses subject to the two excitations speci ed in 8 Observe that for the case off2t 02 x sin 098 gtIlt 7U both mass m1 and mass m2 are moving in phase that is they move in the same direction in time However when the system is subjected to f2t sin 098 gtIlt 26187rt mass m1 and mass m2 are moving out phase that is they move in the opposite directions in time In other words depending on the excitation frequency the motions of the two masses are drastically different To understand this strange phenomena one has to understand the roles of modes and mode shapes To this end let us recast l in a matrix form with c 0 m1 0551k1k2 k2x1f11 9 0 I712 552 k2 k2 x2 f2 1 The characteristic equation of the above coupled 2dof differential equation9 can be obtained as follows First we assume the solution of their homogeneous equations in the form 961 951 jwt Mass Positions for Mode 1 2dof system response due to harmonic input at mass 2 with omega 098omega1 2dof system response due to harmonic input at mass 2 with omega 098 omeg32 25 25 2 m2 2 m Position of 7 2 Mass 2 Position of 392 Mass 2 39U 0 E 15 7 7 15 7 k2 E 9 E 1g CL 1 7 m 7 a 1 M I I J 5 Position of Position of Mass 1 Mass 1 0 5 05 7 39 k1 quot 00 0 2 0 4 0 6 1 1 4 1 6 1 3 2 04 o Ti1me 1 1 4 1 a 2 39 39 Time I I I I Figure 2 Motions of two masses under two different excitations Second substituting this into 9 with f1 2 f2 0 yields w2n1 n32k12k2 20 11 Hence the characteristic equation is obtained by requiring that the above equation has a nontrivial solution k1 k2 k2 kikz a2 0 m 1 m 2 m 1m 2 12 Note that with m1 2 m2 1 k1 2 k2 26182 the two roots of the above characteristic equation are given by def k1 k2 602ml k2 4 k2 k2 w2m2 O gt a com 2 7T 012 2618 13 These two values are called characteristic values Whose square roots are called the natural frequencies or vibration modes of the system The corresponding eigenvectors can be computed from the second row of 1 1 k 2 2 2 w l w 14 2 ax kx kx 0 gt 2 21 22 x2 k2 k2 Note that for the two modes computed in 13 we have two different expressions x1 7T2 Forwzn 1 20618 x2 26182 2 2 15 1 2618 It Fora 2618 1 1618 x2 26182 The ratios of these eigenvector sets are plotted in Figure 3 Mode Shape for Mode 1 Mode Shape for Mode 2 Fixed end and 00 06 10 1 0 10 Modal Amplitude for Mode 1 Modal Amplitude for Mode 2 6 Mass F0 nt Figure 3 Mode Shapes of 2DOF Example Problem Observe from Figure 3 that for the case of the rst mode when mass 1 moves in the same direction with its amplitude of 0616 while mass 2 moves a unit amplitude In other words the two masses move in phase as illustrated in Figure 2 As for the second mode mass 1 moves in the opposite direction by l618 while mass 2 moves a unit amplitude which is also illustraed in Figure 2 From these observations we conclude that mode shapes indicate how the system will deform when subjected to a harmonic excitation whose frequency is close to one of the natural frequencies or vibration modes of the system Thus modeshape information is useful in designing structures subjected to harmonic excitations Examples of such systems include propellered airplanes ships motor vehicles and many other machinery equipment 3 Vibration analysis of 2DOF SpringMass Problem by MATLAB In Matlab we invoke the following routine X D eigK M where X is the eigenvector D is the eigenvalues For example the eigenproblem of the 2DOF system studied in the preceding sections can now be analyzed by MATLAB with 226187I2 26187r2 1 o K 2618n2 2618n2 M o 1 16 Upon using the following MATLAB code we nd the following results o 2 dof eigenvalue analysis stiffness matrix K 22618pi 2 2618pi 2 2618pi 2 2618pi 2 mass matrix M 1 O O 1 call eigenvalue routine X D eigK M write eigenvector and eigenvalues lambda diagD disp Eigenvalues of 2dof system num2strlambda disp disp Eigenvectors of 2dof system num2strX1 disp num2strX2 compute the frequencies freq diagD extract the two diagonal entries of D matrix disp freq sqrtfreq disp Frequencies of 2dof system num2strfreq we nd the following results Eigenvalues of 2dof system 676464 986948 Eigenvectors of 2dof system 085065 052573 052573 085065 Frequencies of 2dof system 82247 31416 Note that MATLAB prints out the highest mode rst Hence the mode shapes are given by le For the second mode With wz 2618 82247 X2 x22 415257 X11 For the rst mode With col 7T 314157 X1 m2 0850650 17 The mode shapes computed by MATLAB is normalized so that their vectorial length is unity If the modal amplitude at mass 2 xzz for the second mode shape is scaled to be unity from MATLAB computed value xzz O52573 then we would have x21 2 l618 as given by 15 The same is true for the modal amplitude for mode 1 at mass 1 and mass 2 points Finally one may express X2 as X21 O85065 For the second mode With wz 2618 822474 X2 x22 052573 18 which indicates that X2 and X2 are equally valid representation 4 Den Hartog Ormondroyd Oscillators An important application of the twodof vibration model is to study shock isolation design of a suspension system as shown Figure 1 If Figure l is rearranged as shown in Figure 4 the same model can be used as a resonator model whose generic charac teristics can be used for the design of accelerometers and more recently as lters of microelectromechanical systems MEMS In either case the governing equation 1 is applicable Let us parameetrize the model as 601K1M7 602k2m2 C2 m2w27 Mm2M1 19 and normalize the frequence response X1a given in 6 by its static displacement F2 2 OI R I l 3amp1 22 N Hm W N W Two DOF SpringMassDamper Suspension Model Figure 4 2DOF Suspension Model xst1 FlKZ 20 such that from 7 we have H X1a K1H22a 1 ml H11wH22wgt Hm 21 Dividing both the nominator and the denominator by M1 m 2 and utilizing the pa rameters introduced in 20 the frequency response function H1a at mass 1 can be expressed as 19220 H1 2 A A A2 H110 H22CU LI 1120 F111 w w2 12mw2w a mi 22 P122 60 602 1200260 60 19120 1200260 60 Amplitude Frequency Response Functions for Different Damping Ratios of a 2DOF Suspension Model Invariant points 523 1o 11 Frequency Hz Figure 5 Frequency Response Function at Mass 1 Figure 5 illustrates the frequency response function FRF for different damping ratios with w1w2lO Hz u0l 23 Note that there are two points that all FRF curves pass through two invariant points It is these invariant points that were rst discovered by den Hartog and Ormondroyd in 1928 who subsequently utilized their properties for improved design of mechanical shock attenuation In their design study by using Figure l den Hartog and Ormondroyd observed that if the magnitudes of the two invariant points are made to be the same by varying the mass ratio and the damping ratio a near optimum design goal can be achieved For the design of accelerometers using Figure 4 the objective is to maximize the amplitude of the frequency response of mass 2 in order to decrease the signal to noise ratio while realizing a at plateau of frequency interval of interest For the design of resonators using Figure 4 in addition to maximizing the amplitude of FRF at mass 2 the gap between the two invariant points should be minimized MATLAB code that produced Figure 5 islisted below XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX X hartogoscillatorm X by k c park 15 january 2002 clear X Define some useful numbers dtrpi180 Hz2rps2pi rps2Hz12pi X Conversions tofrom Hz and rads X Natural frequencies of individual masses 10Hz w110Hz2rps w2w1 mu 01 for zeta O005O8 Xzeta is the damping ratio X Define a system in transfer function form aO1 a1 2zetaw2 a2 w2 2 b01 b1 2zeta1muw2 b2 W1 2w2 2 mu14zeta 2w2 2 4muzeta 2w2 2 b3 2zetamuw2 3 w2w1 2muw2 2 4zetamuw2 3 b4 w2 2w1 2muw2 2 muw2 4 sys tf a0 a1 a2 b0 b1 b2 b3 b4 X Pick a range of frequencies to look at from 0 Hz to 100 Hz wf logspaceO515400Hz2rps X Find the frequency response Hfreqrespsyswf X This ends up being a 3 dimensional array HreshapeHsizewf X this makes it the same size as wf X Plot the amplitude as a function of frequency in Hz loglogwfrps2Hzw1 2absH grid Xlabel Frequency Hz ylabel Amplitude hold on end aXis5 30 01 30 XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX References 1 R E D Bishop Vibration Cambridge University Press 1979 2 J P Den Hartog Mechanical Vibrations Dover Pub 1984 3 Maurice Roseau Vibrations in Mechanical Systems SpringerVerlag 1984 I Theoretical Basis 0 Kinetic energy of the contimuun bar L 8ux t Thar mm W 02cm w r 0 at 0 Potential energy of the continuum bar 1 L 2 8 36 1 Vbarzj uXx9t dx9 uXx9t ax O 0 External energy due to applied force ft 3W3 ft MSG dx 1 2 3 0 Kinetic energy of the discrete mass 1g Mmafagt e 0 Potential energy of springs Kwd0 MLOF 6 o Hamilton s principle for the bar and discrete massspring system 8 Tlolal 5Vbar 5V3 5Wbar dt 2 0 Tlolal Tbar 4 Ts 6 l1 0 Variation of the foregoing Hamilton s principle 12 L f f 0xiixa 1 E1406 Mm 1 5M 10736 ll 0 Miis1 K us1 La 1 f15us1 7 EAOC x I K MD I ust SMOC 1xL EAOC uxx I5ux 1x0gt 0391 0 o The governing equations of motion pxb lx t EA umx t 2 0 Maslttgt K use ML 0 ft o The boundary conditions u0 t 0 EAx uxL I K uL I Mm 0 8 9 10 11 Step 1 same as in HWl Obtain the position vectors rA a xi L rC a x1 3003 61 Sln 6 r3 a xi LCos 9i sian 12 r A a1 L rC a1 51 r3 ai Li The Virtual displacement needed for computing the generalized forces 8rB 8x i L sin 6i cos 6j89 13 Step 2 Obtain the kinetic energy of the system L9 rC x1 7 sm6100363 14 T 2 er to ace 0 L292 L2 15 mamp2 Lxe sine T 711 2 2 Step 3 Compute the potential energy due to gravity and spring and the nonconservative Virtual work Potential energy due to spring F3 2 kArA kxi drA dxi 1 O 16 Vspring Fierng drA kx2 x Potential energy due to gravity L F gamy mgi drC dxi 5 sin 6i cos 6jd6 on 391 rref s rin Vgravzyzf FCR godrC rc O O L L 2 mgdx sinedez mgx l COSQ x 2 9 2 17 Total potential energy L V Vsprmg Vgravlly kx2 mgx 1 COS 9 18 Nonconservative generalized forces 730160 F sin 9i cos Qj 81393 8xi L sin 6i cos 6j89 u 5Wn0ncons Fign60 51393 Fsin6 8x FL 36 Qx8x Q986 u sz Fsine ngFL 19 Step 4 Obtain the equations of motion Via HamiltonEulerLagrange s formalism The generalized coordinates ql x 612 6 The Lagrangian L292 L L T V m c2 LXQ sin6T kx2 mgxl Cos6 20 The generalized forces QC F sin 6 Q9 FL Substituting the preceding expressions into d 8L 8L Q 6 a Z x 2 dt aqk aqk k 611 612 we nally obtain mL mL 2 mx st 0039kx2mg Fsm6 mLz mLic39 mgL sm62 sm9FL 3 21 22 12 Modeling for Structural Vibrations Chapter 12 MODELING FOR STRUCTURAL VIBRATIONS 12 2 121 FEM MODELING OF BEAM VIBRATIONS In modeling of cable Vibration problems by linear elements it was observed that one needs to model the cable by more than 50 elements if the fundamental frequency is to be computed within 4 digit accuracy What was not discussed therein is the accuracy of mode shapes To gain further insight into nite element modeling of Vibration problems let us consider the modeling of plane beam Vibration problems For simplicity a beam with simplesimple supports is used to model 5th modes Classical theory tells us that we have the following solution k E kth mode wk 7T L p 121 knx kth mode shape Wx sin T where E is the bending rigidity 0 is the mass per unit beam length and L is the beam span 12 2 12 5 121 FEM MODELING OF BEAM VIBRATIONS In using the nite element method for modeling of beam vibrations the question arises how many elements does one need for an accurate computation of its kth mode and mode shape While there eXists a vast amount of literature on the accuracy and convergence properties of various elements which one may employ to answer the question we will adopt a posteriori assessment approach To this end let us take a simply supported beam and concentrate on its 6th modes The frequency error vs the number of beam elements used are shown in Fig 121 As can be seen vedigit accuracy is achieved with about 45 elements If all one needs is the frequency with one percent accuracy one could use only 10 elements In Figure 122 the frequency error of a xedsimply supported beam is illustrated Although the error for this case is a little higher than that of the simply supported case the frequency error converges with the same trend Let s now focus on the mode shape accuracy 12 5 Chapter 12 MODELING FOR STRUCTURAL VIBRATIONS 12 12 Figures 12 3 5 illustrate how the mode shapes converge to the exact solution as the number of elements increases Also plotted is spline tted curves that utilize only the sample points or computed discrete mode shape points It is clear that curve tting in general enhances the discrete raw data points especially for the case of crude models Scanning over the six mode shape plots Fig 1235 vs the increasing number of elements an acceptable number of element for capturing the sixth mode shape appears to be 30 Note that there are six half sine waves in the 6th mode shape Hence for the 30element model each half sine wave is sampled by ve elements or six nodal points meaning that ten elements span one full sine wave To conclude ten elements captures the sixth mode frequency with less than one percent error whereas for an adequate mode shape capturing three time of that elements 30 elements are needed This is often referred to as threetoone rule Finally the case of tenelement model satis es the socalled Shannon s sampling criterion which state that for the minimum sampling number for a sinusoidal signal is three 122 CHARACTERIZATION OF STRUCTURAL DYNAMICS EQUATIONS Let us now study the characteristics of the secondorder damped system MiiDuKu ft D aM8K 122 where or and 8 are constants The damping matrix D 05M 8 K is called a Rayleigh damping as it is proportional both to mass and stiffness of the system 12 12 12 13 122 CHARACTERIZATION OF STRUCTURAL DYNAMICS EQUATIONS The coupled equations of motion for a linear structurel22 can be decoupled by the n x n eigenvector matrix T which relates the displacement vector u to a generalized solution vector q via 11 Tq 123 Substituting 123 into 122 and premultiplying the resulting equations by TT one obtains Iqozl8Aq I Aqfq fq TTff 124 in which TTMT I 125 TTKT A diagaf w v 126 where a is the jth undamped natural frequency component of 124 The solution of 124 for its homegeneous part can be eXpressed as q cest 127 where s is in general compleX constants Substitution of 127 into 124 with fq 0 yields 21 a1 8As Ac 0 128 Equation 128 has a nontriVial solution only if det s21 a1 8As Al 0 129 12 13 Chapter 12 MODELING FOR STRUCTURAL VIBRATIONS 12 14 from which the kth solution component can be expressed as qk akes kes 1210 where sk and k are the kth complex conjugate pairs of 129 and ak and 51k are arbitrary constants Two characterizations of 129 are possible frequency characterization viz by xing the damping parameters or and 8 and varying the undamped natural frequency a and damping characterization viz by xing the undamped natural frequency a and varying the damping parameters or and 8 Figure 126 illustrates the root sh of the characteristics equation129 when the nondimensional frequency ooh is varied where h has the dimension of time In the case of stepbystep time integration of the governing equations of motion 122 h may represent the integration time step size We will examine different parameter choices to understand the root loci drawn in Fig 126 12 14 Chapter 12 MODELING FOR STRUCTURAL VIBRATIONS 12 16 1221 Frequency Characterization Based on the Root Locus Method The characteristic equation 129 represents n roots s2ot8w2SJw20 j12n 1211 Since not only the distribution but also the maximum and minimum frequencies are not known apriori a complete range of the Rayleighdamped system can be characterized by the following expression S2oz8a2sa20 0560500 1212 In order for the subsequent characterization to remain valid for all time sampling rates or integration stepsize we normalize 1212 by the time stepsize h such that 2 ah gwhf wh2 0 0 5 ooh 5 00 U 1213 dwh 26 8c62 c620 azah 3 3 Elbe Note rst that the rigidbody motion ie c coh 0 corresponds from 1213 to s1 00 and s2 ah 0 55 0 1214 These two roots are marked as s1 and s2 in Fig 126 12 16 12 17 122 CHARACTERIZATION OF STRUCTURAL DYNAMICS EQUATIONS As the undamped natural frequency coh is increased the root loci move from the initial points s1 s2 meeet at the branch point at A as double roots or hlt1 N1 cam 0 lt1 1 Sim3 0 1215 As coh is further increased the root loci split into two conjugate paths forming a half circle whose center is located at 1 8 0 with its radius of x l 65878 For example at points B indicated in Fig 126 we have 12 13 i 1 5633 1216 As coh is further increased the two complex roots merge at h1 1 o 0 1217 SC For coh gt lsc one locus branches out toward the negative in nite axis while the other approaches the origin 0 0 This is also illustrated in Fig 126 For the special case 8 0 viz mass proportional damping the branch occurs at s 0 when a 05 2 Hence the solution components with frequency coh lt exhibit overcritically damped responses as the locus becomes a straight line ERG 0 2 h In particular if or 0 the root loci coincide with the imaginary axis whose response is characterized by purely oscillatory components 12 17 Chapter 12 MODELING FOR STRUCTURAL VIBRATIONS 12 18 From the physical viewpoint the case of massproportional damping introduces higher modal damp ing for lower frequency solution components and effective damping decreases as the frequency increases This does not however necessarily mean that the response components of the high frequency modes will decay slower than those of the lowfrequency modes within a time period As a matter of fact the decay rate is uniform for all frequency components since we have leStl W for all a 1218 The root locus of the stiffnessproportional damping case or 0 touches the origin and forms the half circle as a is increased As a is further increased the locus branches out into the negative real aXis Therefore the decay rate increases as a increases As such this representation of system damping is often used in the modeling of structural damping due to joint effects acoustic noise and internal material friction 1222 Damping Characterization The conventional characterization of damping into the equations of motion for structures is to eXpress each component of 129 Si 256 a 0 1219 so that the characteristic roots can be eXpressed as sk c i 1 525 1220 where S is termed the damping ratio The equivalent damping ratio for 1211 therefore is 1 oz Seq 560 1221 2 a 12 18 12 19 123 VARIOUS DAMPING MODELS Now we like to examine for a xed a what happens when S is varied For the undamped case viz Seq 0 the roots lie on the imaginary axis s 0 tab As 5 is increased the roots rotate toward the lefthand complex plane with the constant magnitude of c coh since 2 wk 1222 and the shift angle being 6 tan 1 l 525 The two complex roots merge on the real axis when 5 1 thus the root locus is a halfcircle with radius wh For 5 gt 1 one root branches out to the in nite negative axis while the other approaches the origin see Figure 126 This invariance property that is the magnitude of the complex roots remains the same for all damping ratios plays an important role both for active controller synthesis and computational algorithms 123 VARIOUS DAMPING MODELS Modeling of damping remains a challenge in structural dynamics Over the years various damping models have been proposed Below lists a sample of parameterized damping models 1231 MassProportional Damping Model If the damping can be made proportional to the mass such as dynamic friction cases the damping matrix can be parameterized according to D 0M 1223 as already discussed in the previous section The root loci are a special case of 122 governed by 12 19 Chapter 12 MODELING FOR STRUCTURAL VIBRATIONS 12 20 s2 a c62 0 1224 which is plotted in Fig 127 as model A with its parameter or 1 Note that the root loci form a straight vertical line when a gt at 2 1232 Viscous Damping Model One of the most widely used damping models is the Viscous damping characterization This has been adopted for modeling of coated damping layers of lubricants in rotating machines and of a plethora of unknown sources Mathematically the Viscous damping model can be expressed as D 2 2T 5wP1TT gwp diag 1af 52 mg 5N 605 1225 TTMT I TTKT A Therefore the characteristic equation of a Viscous damped system is given by g 25 7 3 32 0 1226 Notethatthetwocasesp 0 5152 SN 05 and p22 5152 SN 28 have been studied in the preceding section The root loci of a Viscous damped case when combined with a massproportional damping are obtained by p21 515225N25 gt 2a2sa a20 1227 12 20 12 21 123 VARIOUS DAMPING MODELS which is plotted in Fig 127 labeled as model C with its parameters or 10 S 0025 Note that in this model the modal damping ratio is the same 5 const for all the frequencies greater than co gt a2l 5 lt 10 12 2 1