Home

User's Manual for BECAS - DTU Orbit

1. 0 1 0 05 y 0 1 0 05 e ky 1 f amp 1 Figure 4 17 Cross section warping displacements for cylinder cross section with two isotropic materials layered C4 displacements not to scale 60 CHAPTER 4 VALIDATION 4 2 3 Three cells In the final group of numerical examples a three cell sectional geometry is considered The geometry corresponding dimensions and finite element mesh for the three cell cross section are presented in Figure 4 18 Some cross section analysis theories rely on the integration of the shear flux over each of the closed cells The aim of this example is to illustrate the ability of BECAS to correctly estimate the stiffness properties of cross sections regardless of the number of cells in the cross section 0 01 0 08 0 01 0 01 0 2375 0 01 0 2375 0 01 0 485 0 01 Figure 4 18 Geometry and finite element mesh of three cells cross section of one isotropic material T1 Three cells cross section of isotropic material T1 In this example all faces of the three cell cross section are made of isotropic material 1 The resulting non zero entries of the cross section stiffness matrix as determined by both BECAS and VABS are presented in Table 4 15 The estimated positions Table 4 15 Non zero entries of cross section stiffness matrix for three cells cross section T1 Comparison between BECAS and
2. Ou Ow L Oz i Oz AO j e allel bp Eu Rw 0 or in matrix form M 0 a Ory m Oz Assuming solutions of the type er 0 SP eee u te Y pe and introducing the previous equation in 2 20 yields M 0 C C L E R 2 EROS The solution to 2 20 can be stated in terms of a linear combination of the solutions to the eigenvalue problem SACOS 14 CHAPTER 2 THEORY MANUAL where A are the eigenvalues and and 4j the corresponding eigenvectors The eigen vectors represent the extremity modes and the corresponding eigenvalues define a diffusion length The lowest eigenvalues are the most interesting as they propagate farther into the beam Choosing the eigenvectors corresponding to the lowest eigen values it should be possible to study the effect of the loads at the extremities of the beam Nonetheless the solution to this eigenvalue problem may be cumbersome as the size of the matrices becomes larger 2 3 2 Central solutions The central solution refer to the solutions obtained at the central part of the beam corresponding to non zero stress resultants 0 4 0 For convenience let us rewrite the equilibrium equations in 2 19 so that all the terms with derivatives of the same order are grouped l Eu Rp C C7 224124 M Ru Ay LT 2 40 At this point the reader may opt to read Section 2 4 first to get an insight into the mathematical properties of the central solutions
3. Part of the input is coming from the BECAS_Constitutive_Ks function The extra input consists of e theta0 1x6 array holding the cross section generalized forces and moments i e the variable 0 ge MT defined in Section 2 2 1 The magnitude of the entries are generally determined based on the one dimensional beam finite element solution The output is e ElementStrain GlobalCS 6xn array holding the three dimensional strain components evaluated at the center of each of the elements in the cross section finite element mesh e NodalStrain_GlobalCS 6 x ngp x Ne array holding the three dimensional strain components evaluated at the each of the Gauss points ng Gauss points at each element in the cross section finite element mesh Note function is not yet validated 5 2 6 BECAS RecoverStresses The function BECAS RecoverStress is used to determine the three dimensional stress components at each element or at each Gauss point in the cross section finite element mesh The function is meant to be used after the one dimensional beam finite element solution has been determined The stresses are evaluated based on the strains determined using the BECAS RecoverStrains The resulting stress values are given in both the global and material coordinate systems The function is called as ElementStress GlobalCS ElementStress MaterialCS NodalStress GlobalCS NodalStress MaterialCS BECAS RecoverStress ElementStrain
4. 1 0 0000 0 1 0 000 0 0 1000 Tr jo o y 100 00 z010 PERETTI M TrMT c s 0 0 0 0 s co 0 00 001000 TESI ee ae ee 0 0 0 s c 0 O 00 0 01 Chapter 4 Validation In this section numerical results obtained using BECAS the BEam Cross Section analysis Software are presented BECAS is the current implementation of the method presented in the previous sections At this point the output is the cross section compliance and stiffness matrix and the positions of the shear and elastic centers The resulting entries of the cross section stiffness matrix K as well as the positions of the shear and elastic center are compared to the results from VABS the Variational Asymptotic Beam Section analysis code Yu et al 13 14 VABS has been extensively validated against different cross section analysis tools and analytical results see Yu et al 13 14 Chen et al 15 and Volovoi et al 17 and it is therefore a benchmark for validation of new cross section analysis codes The numerical experiments presented for validation have been chosen so that different material and geometrical effects are analyzed From a material properties standpoint the aim is to analyze the effect of material anisotropy and its inhomo geneous distribution over the cross section In terms of cross section geometry we look at solid thin walled open and multi cell cross sections This chapter is organized as follows First the setup for the numeri
5. 1 1 VERSION HISTORY 3 Chapter 5 User s Manual The user s manual for the MATLAB implementation of BECAS is presented here This chapter covers the practical use of BECAS as a cross section analysis tool 1 1 Version history e Version 2 0 Authors JPBL ROBI Date 09 02 2012 Change First stable version e Version 2 1 Authors JPBL ROBI Date 23 02 2012 Change The sign of the fiber and fiber plane orientation angles have been switched such that the results beam displacements and cross section stresses now match the ABAQUS results The calculation of the shear center position now neglects the bend twist coupling terms z 0 in the calculations CHAPTER 1 INTRODUCTION Chapter 2 Theory manual This chapter describes the theory underlying the implementation of the cross section analysis tool BECAS The chapter is organized as follows In the first section Section 2 2 some general definitions are introduced The beam kinematics are subsequently described The displacement of a point in the cross section is described as the sum of a rigid body motion and a warping displacement accounting for the cross section deformation A two dimensional discretization of the warping following the typical finite element approach is introduced The principle of virtual work is then invoked in the derivation of the expressions for the external and internal virtual work per unit length The equilibrium equations for the cross section are conseque
6. 1 2000E 00 0 0000E 00 Ks 44 1 0004E 03 1 0004E 03 7 2271E 10 Ks 55 1 0002E 03 1 0002E 03 6 6684E 10 K 5 66 8 5081E 04 8 5081E 04 9 5579E 10 Ks13 0 0000E 00 0 0000E 00 0 0000E 00 Ksa6 0 0000E 00 0 0000E 00 0 0000E 00 results from BECAS and VABS For laminate orientations other than 0 and 90 the shear extension and bending torsion couplings become non zero These effects are typical of laminated composite beams The influence of these couplings on the cross section warping displacements can be observed in Figures 4 7 and 4 8 Comparing with the isotropic case in Section 4 2 1 the bending deformation resulting from the curvature 7 1 is now accompanied of an out of plane distortion induced by the torsion coupling Finally the position of shear and tension center do not depend on the orientation of fibers in the laminate i e 7 ys t y 0 However note that according to the definition of the shear center see Section 2 5 2 in this case the position of the shear center is not a property of the cross section Instead because the bend twist coupling K 46 4 0 the shear center position will depend linearly on the coordinate z 48 CHAPTER 4 VALIDATION a To 1 b r 1 y 0 05 0 05 x c Ty 21 d ty 21 0 02 0 01 0 01 0 02 0 02 0 03 0 04 0 01 0 04 0 03 0 02 0 01 0 x e T 1 f T2 21 Figure 4 7 Cross section warping displacements for square cros
7. E 100 100 a 480 Eso Ey 100 100 a 120 Gyz 41 667 41 667 a 50 Grz Gay 41 667 41 667 a 60 Vyz 0 2 0 2 a 0 26 Vrz Vry 0 2 0 2 a 0 19 Table 4 2 Catalogue of cross section properties analysed for validation Ref Geometry Material S1 Solid square Isotropic 2 Solid square Isotropic 1 Isotropic 2 S3 Solid square Orthotropic Cl Cylinder Isotropic 1 C2 Half cylinder Isotropic 1 Isotropic 1 Isotropic 2 Isotropic 1 Isotropic 2 T1 Three cells Isotropic 1 T2 Three cells Isotropic 1 Orthotropic C3 Cylinder C4 Cylinder layered VABS Recall the relations between the strains and direction of forces and moments T Ki Ks1 Ksi3 Ksia Ksi5 Ks16 Ty Ty Ks22 Ks23 Ks2a Ks25 Ks26 Ty eoe Ks33 Ks3a Ks35 Ks 36 Tz M Ksa Ksas Ksa6 Ke M K 55 Ks 56 Ky M sym Ks 66 Kiz The warping displacements are presented only for a qualitative analysis of the solu tions o E Z Figure 4 1 Cross section coordinate system 4 1 SETUP 39 The accuracy of K depends on the size of the cross section finite element mesh Just like in standard finite element analysis a mesh convergence study should be performed in order to establish the minimum size of the cross section finite element mesh required to obtain realistic values of K Nonetheless since BECAS is only being compared to VABS which is also based on a finite element discretization of the cross section a mesh converge
8. M n pdA 2 2 A where the two dimensional cross product matrix n is 0 0 y n 0 0 cx y x 0 and thus x p n p where x y z is the vector position of a point in the 2 T cross section The vector of section forces 0 T7 MT can then be written as 2 2 EQUILIBRIUM EQUATIONS 7 Mo dz Figure 2 2 Cross section resultant forces for a slice dx of the beam ERN Total deformation Rigid body In plane and translation and rotation out of plane warping Figure 2 3 Schematic description of the different contributions for the deformation of the cross section 0 J Zp dA 2 3 A where the matrix Z I3 n7 and I is the i order identity matrix 1000 0 y Z 01 00 0 z 001 y x 0 2 2 2 Kinematics The displacement s sz sy se at a point of the cross section is defined as s v gt 8 where v v vy vz is the vector of displacements associated with the rigid body translation and rotation of the cross section The vector g gz 9y 9 is the vector of warping displacements associated with the cross section deformation see Figure 2 3 Assuming small displacements and rotations the rigid displacements v can be obtained as v Zr a linear combination of the the components of r DE eT The components of x z Nas Xu Xz represent the translations of the cross section reference point T while p z Vr Py Pz can be rewritten as are the cross section rotations
9. should be followed in the rotation of the material constitutive matrix 1 The array is assembled based on the engineering strains y remembering to observe the 1 2 factor 2 The strains are then rotated employing m L L 3 Having m it is possible to assemble the vector of engineering strains Ym remembering to multiply by 2 the natural shear strain components 4 The stress strain relation Om Q y is invoked to determine Om 36 CHAPTER 3 IMPLEMENTATION MANUAL 5 The array with the stress components m is then assembled to evaluate the stresses in the problem coordinate system p L76 L 6 Assemble the vector op based on the components of Gp 7 Finally each of the stress components is a function of the strain components in Yp the direction cosines in L and the entries of the constitutive matrix Qm The coefficients multiplying each of the strain components in y are the components of the constitutive matrix in the problem coordinate system Q The procedure is used to orient the fibers in a laminate or to orient the laminate plane When orienting the fiber plane the rotation matrix L will be cosa sina 0 La sina cosa 0 0 0 1 whereas for the fiber orientation cosB 0 sing Lg 0 1 0 sinG 0 cosB which are the typical two dimensional rotational matrices 3 3 Rotation and translation of constitutive matrices The translation matrix Ty is obtained from static considerations as follows
10. specify e Material properties e Orientation of the laminate plane e Orientation of the fibers laminate Based on this input the first step consists of assembling the material constitutive matrix in the material coordinate system based on a set of material properties 3 2 1 Definition In the case of orthotropic materials the stress strain relation or generalized Hookes Law is stated as 022 Ez US i Em 22 033 Eo E33 Er 33 023 _ 0 0 20 0 0 723 012 y 0 0 0 0 0 Y12 013 0 0 0 0 d o 113 011 HEC d OO ou n Q where the material properties are given in the material coordinate system see Figure 3 3c and defined as e Lj the Young modulus of material the 1 direction e E2 the Young modulus of material the 2 direction e E33 the Young modulus of material the 3 direction e G5 the shear modulus in the 12 plane e Gig the shear modulus in the 13 plane e Go3 the shear modulus in the 23 plane e 112 the Poisson s ratio in the 12 plane e 113 the Poisson s ratio in the 13 plane e 123 the Poisson s ratio in the 23 plane o the material density 34 CHAPTER 3 IMPLEMENTATION MANUAL Z gg p l b c Figure 3 3 Determination of the material constitutive matrix at each finite element of the cross section mesh a Definition of the cross section coordinate system XY Z and element coordinate system xyz b Convention adopted for the rotation of the element coordinate system xyz into the fiber p
11. z y u zi yi 2 2 30 and replacing 2 30 in 2 29 yields 9W B dul B NTSTOSN A 3 Er ST QSN u d N B QSN u dA 3 Ap ES f u N B QBN u dA J du NTSTQBNou dA A A Oz From the equation above it is possible then to identify the matrices M NTSTQSN dA A C y NT B7QSN dA A E N B QBN dA A where M and E are symmetrical while H C C7 is skew symmetrical Rewriting the total virtual work expression OW ru Ou Se H E ae 0u 52 de u The following must hold for an arbitrary variation of du 9 u Ou M H Eu 0 2 31 Oz Oz est This concludes the derivation of the equilibrium equations for the cross section considering only the warping displacements The aim now is not to solve the equation 2 4 ON THE PROPERTIES OF THE SOLUTIONS 21 above but rather discuss the properties of the solutions for this second order linear homogeneous differential equation The general solution is sought as a linear combination of u z he Inserting the above into 2 31 yields 4M AH E h 0 where h are the eigenvectors associated with the eigenvalues resulting from the solution to XM AH E 0 According to physical considerations and Saint Venant s principle we expect to have a self balanced 0 0 solution associated with the exponentially decaying modes A 0 extremity solutions and a solution presenting no decay A 0 for which the stress resultants ar
12. 04 du 2 and dr thus leading to the following set of equati OS g to the following set of equations MOR Oz C Eu Ry 2 7 RTu Ad 0 00 T 0 Oz T It is possible to further simplify the former by differentiating the first equation with respect to z Pu ron Ov OP NEUE Oz Oz Oz du OP C Bu Ry 7 1r A Se Oz Pe TTo Oz 2 3 SOLUTIONS TO EQUILIBRIUM EQUATIONS 13 and obtain the equilibrium equations for the cross section Pu q Ou Ow M C C a HL Eu Rp 0 e Ru Ap 6 2 19 Z LA ge e This concludes the derivation of the equilibrium equations of the cross section 2 3 Solutions to equilibrium equations The set of equations in 2 19 admits two types of solutions a particular integral which depends on the boundary conditions or the internal force resultants in this case and a homogeneous integral which corresponds to the eigensolutions when 0 0 The homogeneous and particular solutions will be henceforth referred to as the extremity and the central solutions cf Giavotto et al 6 respectively A more detailed discussion on this topic is presented in Section 2 4 2 3 1 Extremity solutions The extremitiy solutions owe their name to the fact that they correspond to the solutions at the extremities of the beam where the loads are applied These are the self balanced 0 0 eigensolutions of 2 19 that is the solutions to the following set Oru T Mzz 1 C C
13. 24E 02 E 03 5 00E 0 1 68E 0 1 75E 0 2 31E 04 4 17E 04 4 37E 03 25E 02 E 04 5 00E 0 1 73E 0 1 75E 0 2 30E 04 4 17E 04 4 38E 03 25E 02 E 05 5 00E 0 1 73E 0 1 75E 0 2 30E 04 4 17E 04 4 38E 03 25E 02 VABS E1 Ez2 Ks 11 Ks 22 Ks 33 Ks 44 Ks 55 K 5 66 Ks 26 Ks 35 E 00 1 00E 00 3 49E 0 3 49E 0 5 91E 04 8 34E 04 0 00E 00 0 00E 00 E 01 5 50E 0 1 28E 0 1 92E 0 2 77E 04 4 59E 04 3 93E 03 13E 02 E 02 5 05E 0 1 38E 0 1 77E 0 2 35E 04 4 21E 04 4 33E 03 24E 02 E 03 5 01E 0 1 68E 0 1 75E 0 2 31E 04 4 17E 04 4 37E 03 25E 02 E 04 5 00E 0 1 73E 0 1 75E 0 2 30E 04 4 17E 04 4 38E 03 25E 02 E 05 5 00E 0 1 73E 0 1 75E 0 2 30E 04 4 17E 04 4 38E 03 25E 02 Diff 96 E E2 Ks 11 Ks 22 K5 33 Ks 44 K55 Ks 66 Ks 26 Ks 35 E 00 0 00E 00 7 15E 04 7 15E 04 7 20E 04 1 20E 07 0 00E 00 0 00E 00 E 01 5 45E 13 7 10E 04 7 17E 04 7 19E 04 1 19E 07 7 17E 04 0 00E 00 E 02 3 96E 13 6 89E 04 7 19E 04 7 20E 04 1 19E 07 T 19E 04 0 00E4 00 E 03 3 99E 13 6 75E 04 7 19E 04 7 20E 04 1 19E 07 7 19E 04 0 00E 00 E 04 2 00E 13 6 73E 04 7 19E 04 7 20E 04 1 19E 07 7 19E 04 0 00E 00 E 05 2 00E 13 6 72E 04 7 19E 04 7 20E 04 1 19E 07 7 19E 04 8 05E 13 of the material 2 vanishes the stiffness values and positions of shear and elastic centers converge to those which would be obtained if only half the cross section was considered Note however the variation of the entry K 11 with respect to the ration E Ez plotted in Figure 4 5 The shear st
14. 3D_Utils un nth Sele ae hey RE dup n eb 71 5 3 4 BECAS_3D_Constitutive_Ks 72 5 3 5 BECAS_3D_CrossSectionProps 72 List of Symbols 3 x HnAWE m lt NE Se ee q tax E SURE Iz Ly Ley Lm Ym Nq Ne Coordinates of a point in the cross section esses 6 Cross section ALCS aaa a 6 E M kb cc E E AE E E E Co EAS E e E AENEA 6 A A 6 Material tensor is A ia 6 Tractions on cross section ssseeeeleeeeeeee eee 6 Section forces shear and axial forces 00oooooooooroooomorrrroo 6 Section forces bending moments and torque ssssssssss 6 Vector Of section ToECeS coeno I OU GU OUR 6 Two dimensional cross product matrix sssseesseseeeeeeee 6 Auxiliary matrix for evaluation of cross section forces 7 Total displacement of a point in the cross section 0 ooooooomoo 7 Rigid body displacement of a point in the cross section 7 Warping displacement of a point in the cross section o ooo ooo 7 Cross section translation and rotation ccc cece eee eee eens 7 Translation of a reference point in the cross section 0 oooo o o 7 Cross section rotation angles 0 eee eee eee eee cor 7 Two dimensional strain displacement matrix 00 eee ee eee 8 One dimensional strain displacement matrix ssssssss 8 Sixth order auxiliary matrix for definition of the se
15. 80E 02 E 04 4 94E 02 6 24E 02 2 98E 0 35E 03 9 12E 04 7 53E 03 80E 02 E 05 4 96E 02 6 24E 02 2 98E 0 35E 03 9 12E 04 7 53E 03 80E 02 VABS E E2 Ks 11 Ks 22 Ks 33 Ks 44 Ks 55 Ks 55 Ks 34 Ks 16 E 00 1 25E 01 1 25E 01 5 96E 0 2 70E 03 2 25E 03 1 94E 18 0 00E 00 E 01 3 99E 02 6 87E 02 3 28E 0 48E 03 1 08E 03 6 78E 03 62E 02 E 02 3 75E 02 6 31E 02 3 01E 0 36E 03 9 29E 04 7 45E 03 79E 02 E 03 4 74E 02 6 25E 02 2 99E 0 35E 03 9 14E 04 7 52E 03 80E 02 E 04 4 94E 02 6 24E 02 2 98E 0 35E 03 9 12E 04 7 53E 03 80E 02 E 05 4 96E 02 6 24E 02 2 98E 0 35E 03 9 12E 04 7 53E 03 80E 02 Diff 95 Ei E2 Ks 11 Ks 22 K5 33 Ks 44 Ks 55 Ks 55 Ks 34 Ks 16 E 00 7 20E 04 7 20E 04 1 68E 13 5 59E 09 7 20E 04 2 17E 00 0 00E 00 E 01 7 20E 04 7 20E 04 3 05E 13 5 33E 09 7 20E 04 7 20E 04 2 15E 11 E 02 7 20E 04 7 20E 04 3 32E 13 4 82E 09 7 20E 04 7 20E 04 2 13E 11 E 03 7 20E 04 7 20E 04 0 00E 00 4 67E 09 7 20E 04 7 20E 04 4 45E 12 E 04 7 19E 04 7 20E 04 7 05E 12 4 66E 09 7 20E 04 7 20E 04 2 21E 12 E 05 7 19E 04 7 20E 04 1 27E 11 4 68E 09 7 20E 04 7 20E 04 2 23E 12 BECAS O VABS El Half square 0 1 0 05 AA RI El E 0 0 zi 3 a a 5 10 10 10 10 10 10 Log E E Figure 4 14 Variation of Ks 11 entry of the cross section stiffness matrix with respect to the ration E1 E for the cylinder cross section C3 Results for BECAS and VABS compared with the solution
16. GlobalCS NodalStrain GlobalCS nl 2d el 2d emat matprops The output is e ElementStress_GlobalCS 6x ne array holding the three dimensional stress components in the global coordinate system evaluated at the center of each of the elements in the cross section finite element mesh 70 CHAPTER 5 USER S MANUAL e ElementStress MaterialCS 6 x ne array holding the three dimensional stress components in the local coordinate system evaluated at the center of each of the elements in the cross section finite element mesh e NodalStress_GlobalCS 6 x ngp x Ne array holding the three dimensional strain components in the global coordinate system evaluated at the each of the Gauss points where ng are the number of Gauss points at each element in the cross section finite element mesh e NodalStress_MaterialCS 6xmgP xn array holding the three dimensional strain components in the material coordinate system evaluated at the each of the Gauss points where ng are the number of Gauss points at each element in the cross section finite element mesh Note function is not yet validated 5 2 7 BECAS Becas2Hawc2 The function BECAS Becas2Hawc2 is used to generate input for RIS s HAWC2 code for the aeroelastic analysis of wind turbines It is presumed that the functions BECAS Constitutive Ks BECAS Constitutive Ms and BECAS CrossSectionProps have been previously called Extra input is required for this function namely the
17. MATRIX 35 The vectors e are unit vectors associated with the axis at each of the coordinate systems The components of the stress in the material and problem coordinate systems are Om Orem Oyym Ciym Crizm yzm Ozzm Op Saz Oyyp Ciyp Oxzp Oyzp Ozzy The 3 x 3 arrays with the stress components are given by Ortm Fxym C zxzm Cxtp Tryp Cuz Om Fytm yym yzm gt Fp Oyap Oyyp Cyzp Cztm Fzym Ozzm Ozzy Ozyp Ozzy where the refers to the array notation Then the rotations in Equation 3 2 2 can be expressed as Gm L6 L p L mL The vector of engineering strains in the material and problem coordinate systems is T Ym E Yoram Yum V22m VYyzm VLZm Your T M azp Vyyp Vezp Vyzp Vrzp Iu The natural strain components may be defined in function of the engineering strains as 1 1 1 1 Yrm 9 Yvym 9 Yvzm Yatp Zsyp 9 Yzzp Em Zym Dm 2X gt Sp Yue Vu zyz 2 zwm 2 Y2Yym Yzzm 2 zzp 2 V2Yp Vzzp where 4 are the engineering strains Since the strains are also second order tensors the relations derived for the stresses are also valid for the strains m L pL L 4L Based on the expressions presented before it is possible to establish the equations for the transformation of the material constitutive matrix The constitutive matrix in the problem coordinate system Q is obtained by transformation of the material constitutive matrix in the material coordinate system Q The following steps
18. Or Z p data f ZTP JA Oz A A Oz Work from rigid displacement dul dp em NTp 4A du f NP da Oz A A Oz Work from warping displacements 10 CHAPTER 2 THEORY MANUAL Using Equation 2 3 it is possible to further simplify the previous equation oW Or 7 00 3257 0970 or 2 202 f NT dA our f rp dA Oz A A Note that it is possible to introduce the strain parameters in the following manner Or 00 dr 00 04 br 5 04017 rTTtO dr TT 0 A O REEL O y 0 TO rro y Ot Oz Oz apr 0 Lm dr T0 61 0 Oz The expression for the external virtual work is then given as 0 go Wes 008 2 or TT 0 670 Oz Oz Op m N p dA sul NT dA A Oz which in matrix form read as dur 4T P Oz poe du oP dr 5 110 2 14 Zz m 0 z where d NTp dA f Nr a The vector P can be seen as the nodal stresses in the cross section finite element discretization as it represents the discretized stresses acting on the cross section face Internal virtual work The internal work or the work done by the elastic strain energy per unit length can be written as OW T Az f se c dA 2 15 Observing the stress strain relation in Equation 2 1 the internal work expression n 2 15 can be restated as E 2 de Qe dA Oz A 2 2 EQUILIBRIUM EQUATIONS 11 Inserting the expression for the strain displacement relation in 2 12 in the equation above yields OW du
19. The total displacement s Zr 8 2 4 8 CHAPTER 2 THEORY MANUAL 2 2 3 Strain displacement relation The strain displacement relation can be written as c Bs 88 2 5 where B and S are defined next The strain displacement relation can then be cast as Exx 0 0x 0 0 000 Eyy 0 0 0y 0 0 0 0 85 8z 2e _ 0 0y 0 0x 0 2 000 d 2 6 ene 0 0 8 0x y 100 l ey 0 0 0 0y Es 010 PRA Ezz 0 0 0 0 0 1 B S Note that this the common linear three dimensional strain displacement relation where the terms 0 0z have been set appart Inserting 2 4 into 2 5 yields Og Or e BZr SZ Bg S gt 2 7 Oz Oz It can be shown that BZ SZT where 000 0 1 0 0001 0 0 100000 0 Te 0000 0 0 0000 0 0 0000 0 0 Utilizing the relation above Equation 2 7 can be written as f e SZ Tir Bg s Oz Oz It is possible at this point to define the strain parameters Or and thus rewrite the strain displacement relation in its final form 08 e SZ Bg S 2 9 Oz T The section strain parameters Y 7 Ty Tz Kg Ky kz are representative of the strain due to the rigid displacement of two adjacent cross sections that remain undeformed The components 7 9x2 py and Ty 2 Yx represent shear strains of the cross section in the x and y directions respectively The component 4 T Oxa is the axial elongation Furthermore kz e 2 and ky Fe a
20. V Hodges D H Cross sectional design of composite rotor blades Journal of the American Helicopter Society 53 3 240 251 2008 Blasques J P Stolpe M Maximum stiffness and minimum weight optimiza tion of laminated composite beams using continuous fiber angles Structural and Multidisciplinary Optimization DOI 10 1007 s00158 010 0592 9 2010 Hansen M O L S rensen J N Voutsinas S S rensen N Madsen H A State of the art in wind turbine aerodynamics and aeroelasticity Progress in Aerospace Sciences 42 285 330 2006 Chaviaropoulos P K et al Enhancing the damping of wind turbine rotor blades the DAMPBLADE project Wind Energy 9 163 177 2006 Giavotto V Borri M Mantegazza P Ghiringhelli G Carmaschi V Maffiolu G C Mussi F Anisotropic beam theory and applications Composite Struc tures 16 1 4 403 413 1983 Borri M Merlini T A large displacement formulation for anisotropic beam analysis Meccanica 21 30 37 1986 Borri M Ghiringhelli G L Merlini T Linear analysis of naturally curved and twisted anisotropic beams Composites Engineering 2 5 7 433 456 1992 Ghiringhelli G L Mantegazza P Linear straight and untwisted anisotropic beam section properties from solid finite elements 4 12 1225 1239 1994 Ghiringhelli G L On the thermal problem for composite beams using a finite element semi discretization Composites Part B 28B 483 495 1997 Ghiringhelli G L
21. f xo i 07 C ML 2X g f 0 Y R L A YO The total virtual work expression becomes ex ITT E CR X0 007 F 0 02X C M L 2X g SOY R LT A YO For any admissible virtual displacement 00X 50 and 00 Y it is possible to obtain an expression for the cross section compliance matrix defined as XVIE CR x F amp C M L eX Y R LT A Y The corresponding stiffness matrix K can be computed as K F 2 35 This result can be used to generate beam finite element models for which the strains can be exactly described by the six strain parameters in y The material may be anisotropic inhomogeneously distributed and the reference coordinate system may be arbitrarily located The stiffness matrix K will correctly account for any geometrical or material couplings 24 CHAPTER 2 THEORY MANUAL 2 5 2 Shear center and elastic center positions The expressions for the positions of the shear and elastic center are presented next The shear center is defined as the point at which a load applied parallel to the plane of the section will produce no torsion i e k 0 Hence assume that two transverse forces Ty and Ty are applied at a point s ys at a given cross section The moments induced by the two forces are M T L z My T L z M Teys Tyxs 2 36 The aim is to find the position s Ys for which the curvature associated with the twist x 0 Thus taking into account the cross se
22. finite element mesh are presented in Figure 4 16 It is assumedin this case that E E 1000 The resulting non zero entries of the cross section stiffness matrix as 0 1 0 05 gt 0 0 05 047 0 0 1 X Figure 4 16 Geometry and finite element mesh of cylinder cross section with two isotropic materials layered C4 isotropic material 1 light and 2 dark estimated by both BECAS and VABS are presented in Table 4 14 As can be seen there is a very good agreement between the two cross section analysis tools The Table 4 14 Non zero entries of cross section stiffness matrix for cylinder cross sec tion with two isotropic materials layered C4 Comparison between BECAS and VABS BECAS VABS Rel Diff Ks11 83114E 02 8 3115E 02 7 1996E 04 K2 83114E 02 8 3115E 02 7 1996E 04 Ks33 3 9784E 01 3 9784E 01 0 00 Kosa 1 8012E 03 1 8012E 03 0 00 K s5 1 8012E 03 1 8012E 03 0 00 Kse 1 5010E 03 1 5010E 03 7 1996E 04 warping displacements resulting from different cross section strains are presented in Figure 4 17 The deformation patterns are very similar to those obtained in the case where only the isotropic material 1 is used The differences can be seen in the shear deformation which shows some local perturbations due to the lower stiffness of middle layer 42 NUMERICAL EXAMPLES 59 0 05 0 05 0 0 1
23. in Figure4 10 As expected the the shear strains are the only which induce out of plane deformation Due to the axial symmetry the torsion strain will not induce any out of plane deformation The same holds for the shear and elastic positions which coincide with the origin of the cross section coordinate system Thus according both BECAS and VABS 2 Ys x yy 0 Table 4 9 Non zero entries of cross section stiffness matrix for cylinder cross section C1 Comparison between BECAS and VABS BECAS VABS Rel Diff Ks11 1 249E 01 1 249E 01 7 200E 04 Kso9 1 249E 01 1 249E 01 7 200E 04 Ks33 5 965E 01 5 965E 01 1 675E 13 Ksas 2 697E 03 2 697E 03 5 587E 09 K ss 2 697E 03 2 697E 03 5 587E 09 Ks66 2 248E 03 2 248E 03 7 200E 04 42 NUMERICAL EXAMPLES 51 0 1 0 1 y 0 1 0 05 Figure 4 10 Cross section warping displacements for cylinder cross section C1 displacements not to scale 52 CHAPTER 4 VALIDATION Half cylinder cross section of isotropic material C2 The cylinder cross section studied in the previous chapter is divided in two to gen erate the half cylinder cross section studied here The geometry and finite element mesh are presented in Figure 4 11 This is an open cross section and the out of plane deformation is significant and should therefore be accounted for in the estimation of the stiffn
24. is possible to conclude that 92M B 02 and so 0 0 The equation above states that the resulting forces acting on the section vary linearly along the beam 2 4 2 Warping displacements In this section we assume that S g The displacement of a point in the cross section is given here as a function of the cross section deformation only Note however that this definition of the displacements entails also the representation of the rigid translation x and rotation q This is in fact the underlying motivation for the use of constraint equations as described in Section 2 3 3 to uncouple the rigid body motions and cross section deformation The stress strain relation in this case is given as og B gt85 and the expression for the variation of the total energy is T jn Aa fos dA Oz A Oz A Recalling the generalized Hooke s law o Qe and noting that p So the expression for the total virtual work of the beam cross section considering only the 20 CHAPTER 2 THEORY MANUAL warping displacements is given by 8 o1 S6 UE af Ce SOR ay a B g S298 da Oz A Oz A Oz T Fsigaas f o 7g 08 ay Pe az f c Bog dA I ots EE da A 90 eg dA f oTB g dA 2 29 A Oz A E dA f e QB g dA ag A Oz NH g B QB g dA J 08 STQB g dA A A Oz B7 QS g dA 5 S QS g dA e Expanding g using the typical finite element approach cf Equation 2 10 g x y z N
25. matrix are presented in Table 4 4 for both BECAS and VABS As can be seen there is a very good agreement between the two approaches The shear and elastic center calculated by both BECAS and VABS are exactly situated at the origo of the cross section coordinate system Thus s Ys x y 0 using both methods The warping displacements are Table 4 4 Non zero entries of cross section stiffness matrix for square cross section S1 Comparison between BECAS and VABS BECAS VABS Diff Ks11 3 4899E 01 3 4900E 01 7 1515E 04 Ks22 3 4899E 01 3 4900E 01 7 1515E 04 Ks33 1 0000E 00 1 0000E 00 0 00 Ks44 8 3384E 04 8 3384E 04 0 00 Ks55 8 3384E 04 8 3384E 04 0 00 Ks66 5 9084E 04 5 9084E 04 0 00 presented in Figure 4 3 for each of the strain components The shear strains Ty 1 and 7y 1 induces a cross section deformation which matches the analytical results presented by Timoshenko and Goodier in 26 The same holds for the warping displacements obtained when the torsional curvature k 1 The Poisson effect is visible in the results obtained for ky 1 and ky 1 with an in plane expansion and contraction of the compression and tension sides respectively 43 42 NUMERICAL EXAMPLES 0 05 0 05 0 05 0 05 b Ty 1 a Ta 1 0 04 0 04 d amp 1 c t 1 0 05 0 05 f amp 1 e amp 1 Figure 4 3 Cross section warping displacements for square cross secti
26. section BECAS_RecoverStresses Calculation of the three dimensional stresses at each point in the cross section e BECAS_TransformationMat Rotate and translate the cross section constitu tive matrices e BECAS_Becas2Hawc2 Generate output for latest version of Hawc2 e BECAS_3D_Utils Build working arrays for BECAS calculations e BECAS_3D_Constitutive_Ks Calculation of the cross section stiffness matrix based on solid finite element models of the cross section e BECAS_3D_CrossSectionProps Calculation of the cross section properties The version of BECAS based on two dimensional finite element discretizations of the cross section is presented first The use of the BECAS_3D tool is discussed last in Section 5 3 Note that the BECAS toolbox also includes a number of functions for visualiza tion of the input and output The description of this functions has been omitted as it falls beyond the scope of theis manual Part of the MATLAB code in BECAS is generated through MAPLE There are four MAPLE files which follow together with the BECAS code 65 66 CHAPTER 5 USER S MANUAL e GenerateMatricesBECAS mw used to generate all the element stiffness matrices 5 1 which are included in the MATLAB function Quad4 used within the function BECAS_Constitutive_Ks LayerRotationBECAS mw used to generate the code for the function which ro tates the material constitutive matrix Q in the MATLAB function BECAS_Utils when orientin
27. 1 Two dimensional finite element analysis 3 1 1 Q4and Q8 elements 2 eee 3 1 2 Local and global finite element matrices 3 2 Material constitutive matrix ooa 3 2 1 Definitio 4 ue eine a aa ed Fg ed 3 2 2 Rotation s ue Bote a Me deal e RR RI eae SR CR d 3 3 Rotation and translation of constitutive matrices vii CON DD c Q 13 14 16 16 17 19 22 22 24 25 26 26 viii CONTENTS 4 Validation 37 AUD Setup ac vite Hoke A a A Bee on Po 37 4 2 Numerical examples eee 40 42 Square as oi oem ee a Eon Rig EE 41 432 2 Cylindert s 74 A aee Xd dE 50 472 3 Phreeccellss ata iion EX E Mon uif ii 60 5 User s Manual 65 bl pute aea e i a sn e usos qb xoa 66 b Examples ies mA TR ek ds REA 67 5 2 List of functions and output o ooa 67 0 2 BECAS CULATA A ino o e tum YI aed E 67 9 2 2 BECAS_Constitutive_Ks 67 5 2 8 BECAS Constitutive Ms 0 0 67 5 2 4 BECAS CrossSectionProps 68 5 2 5 BECAS RecoverStrains lee 69 5 2 0 BECAS RecoverStresses 69 5 2 7 BECAS Becas2Hawc2 een 70 5 2 8 BECAS TransformMat 70 5 2 9 Examples 200 up pais te dE o e Ry te m Enea 71 5 3 The BECAS 3D implementation ls 71 DS dnp t cx e ds ts Heels Sh Ee eas AS ee MR ed 71 5 3 2 List of functions and output soo 71 5 9 9 BECAS
28. 10 Ks31 8 48E 01 8 48E 01 1 13E 10 K 35 5 64E 02 5 64E 02 3 88E 11 Ks15 2 75E 03 2 75E 03 6 11E 09 Kso6 1 49E 02 1 49E 02 2 42E 09 Kso4 3 13E 03 3 13E 03 2 16E 09 Kse4 1 83E 02 1 83E 02 5 79E 10 4 2 NUMERICAL EXAMPLES 63 Table 4 18 Shear and elastic center positions s Ys and zi yz respectively for three cells cross section T2 Comparison between BECAS and VABS BECAS VABS Diff ts 4 300E 02 4 300E 02 Ys 0 0 zx 1 456E 02 1 456E 02 Ut 0 0 SO SiS 64 CHAPTER 4 VALIDATION 0 5 88 o N 00 B8 88 Figure 4 21 Cross section warping displacements for three cells cross section with orthotropic material at 45 in the top faces and isotropic material in the vertical shear webs T2 displacements not to scale Chapter 5 User s Manual The usage of BECAS is described in this chapter The aim is to go through the functionalities of BECAS and discuss its practical usage within a structural analysis environment BECAS is implemented as a Matlab toolbox The following functions are available e BECAS_Utils Build arrays for BECAS calculations e BECAS _Constitutive_Ks Calculation of the cross section stiffness matrix e BECAS Constitutive Ms Calculation of the cross section mass matrix e BECAS _CrossSectionProps Calculation of the cross section properties e BECAS_RecoverStrains Calculation of the three dimensional strains at each point in the cross
29. 8z0 0 and 909 9269 0 n E Ran 0 O n n zm The analysis of the former sets will be done from bottom to top The last set set n corresponds to the n order derivative of the equilibrium equations in set 2 21 It is a linear homogeneous system of equations with unknowns 0 u 82 0 and 04 02 whose solutions are O u dz 0 and 0095 0209 0 These result can now be replaced in set n 1 Hence for set n 1 it is possible to see that Q 7Uu 8z U 0 and 0 U45 02 U 0 also and so on up to set 1 In set 1 the derivative of the surface stresses 00 0z 0 and thus Qu Oz 4 0 and Ow Oz 0 as well It is therefore demonstrated that the displacements u and strain parameters 4 are at most linear functions of z The displacements are obtained from the solution to the following sets E R 0 Riau A9v 80 Oz Oz Oz r 2 24 Eu Ry C C7 24 129 Ru Ay LT22 0 where from Equation 2 19 The set first set in equation 2 24 is solved first to obtain 0u 0z and Ow Oz for a given 0 It is then possible to evaluate the right hand side of the second set and thus obtain u and vy Note that the same coefficient matrix is used twice in the solution and it is therefore possible to decrease the solution time using a proper matrix factorization 16 CHAPTER 2 THEORY MANUAL 2 3 3 Constraint equations The displacement formulation as described in 2 11 is six times red
30. C O Saint Venant s principle and end effects in anisotropic elasticity Transactions of the ASME 424 430 1977 Horgan C O Knowles J K Recent developments concerning Saint Venant s principle Advances in Applied Mechanics 23 180 262 1983 Christensen O Differentialligninger of uendelige rigs kker ed Institut for Matematik Danmarks Tekniske Universitet 30 32 2009 Hogdes D H Nonlinear composite beam theory Progress in Astronautics and Aeronautics 213 2006 Bendsee M P Sigmund O Topology Optimization Theory Methods and Applications 2nd Edition Springer Verlag Berlin 2003 Reddy J N Mechanics of laminated Composite plates and shells theory and analysis 2nd Edition CRC Press 1997 Peters S T Handbook of composites 2nd Edition Chapman amp Hall London 1998 Timoshenko S Goodier J N Theory of elasticity McGraw Hill 1951
31. CROSS SECTION MASS MATRIX 25 and so Fs44Fs53 Fs 45F5 43 Fs 44Fo 55 z F2 PFsasPs 55 Fs45F5 53 Fs 1aFs 55 Ps Tt yt which are the expressions for the position of the elastic center 2 6 Cross section mass matrix The analysis of the cross section mass properties is significantly simpler than the analysis of the cross section stiffness parameters The 6 x 6 cross section mass matrix M relates the linear and angular velocities in to the inertial linear and angular momentum in y through Myy The cross section mass matrix is given with respect to the cross section reference point as cf Hodges 22 m 0 0 0 0 MYm 0 m 0 0 0 Tm M 0 0 m MYm MIm 0 BC 0 0 mm Izz ry 0 0 0 mtm Lry Iy 0 MYm MEm 0 0 0 Izz Lyy where m is the mass per unit length of the cross section The cross section moments of inertia with respect to x and y are given by Iz and Iy respectively while I is the cross section product of inertia The term Izz Ij is the polar moment of inertia associated with the torsion of the cross section The mass and moments of inertia are obtained through integration of the mass properties on the cross section finite element mesh and defined as m 0 0 1 0 0 0 Lux Izy m 0 y zy o dA 0 0 Iy alo 0 z The off diagonal terms are associated with the offset between the mass center posi tion me m Ym and the cross section reference point The position of the mass
32. Ce l N B7QSN dxdy I NTB7QSN J d dr L N s Qsz dxdy f NTS QS2 J d dn M N s qsN dxdy f N S QSN J d dy where the integration is performed using a four point Gauss quadrature cf Figure 3 2 The global matrices are subsequently assembled following typical finite element procedures Ne Ne Ne A Y Aj R S Ron E S E i i 1 i 1 Ne Ne Ne C Y C L YL M 5 M i 1 i 1l i l where ne is the number of finite elements in the cross section mesh Having obtained each of the matrices it is possible to finally solve the cross section equilibrium equations E RD X C CT L x 0 R A 0 Y LT 0 EXE I D 0 0 Ao 0 0 de 0 E RD ox 0 R A 0 i T D 0 0 Ai 0 where D is the matrix of contraint equations defined in 2 3 3 and A and A are the corresponding Lagrange multipliers The two sets above make use of the same coefficient matrix and can be solved efficiently using a proper factorization e g LU factorization The solutions are then obtained doing a forward and backward substitution The cross section compliance matrix is then readily obtained by inserting the solutions of the previous set into T X E C R X F amp C M L ox Oz Y R LT A Y 3 2 MATERIAL CONSTITUTIVE MATRIX 33 A MATLAB implementation of BECAS according to the theory presented above is described in Chapter 5 3 2 Material constitutive matrix At each element of the cross section finite element mesh the user of BECAS must
33. On the linear three dimensional behaviour of composite beams Composites Part B 28B 613 626 1997 Ghiringhelli G L Masarati P Mantegazza P Characterisation of Anisotropic Non Homogeneous Beam Sections with Embedded Piezo Electric Materials Journal of Intelligent Material Systems and Structures 8 10 842 858 1997 13 Yu W Hodges D H Volovoi V Cesnik C E S On Timoshenko like mod elling of initially curved and twisted composite beams International Journal of Solids and Structures 39 5101 5121 2002 73 74 14 15 16 gm EST 18 19 20 21 22 23 24 25 26 BIBLIOGRAPHY Yu W Volovoi V V Hodges D H Hong X Validation of the Variational Asymptotic Beam Sectional Analysis VABS AIAA Journal 40 10 2105 2113 2002 Chen H Yu W Capellaro M A critical assessment of computer tools for calculating composite wind turbine blade properties Wind Energy 13 6 497 516 2010 Jung S N Nagaraj V T Chopra L Assessment of composite rotor blade modeling techniques Journal of the American Helicopter Society 44 3 188 205 1999 Volovoi V V Hodges D H Cesnik C E S Popescu B Assessment of beam modeling methods for rotor blade applications Mathematical and Computer Modeling 33 1099 1112 2001 Horgan C O On Saint Venant s principle in plane anisotropic elasticity Jour nal of Elasticity 2 3 169 180 1972 Choi L Horgan
34. Otherwise the necessary results to retain are that the central solutions u and y are linear combinations of polynomial functions in z with n being the highest degree of such polynomials see Equation 2 32 in Section 2 4 2 Thus TL TL UT 0 and SALE 0z oz Furthermore from the equilibrium equations for a rigid cross section see Equation 2 28 in Section 2 4 1 the following must hold 2 21 0 0 ug 2 22 52 2 22 and recall from the equilibrium equations in 2 19 that 00 T T 0 2 23 A 7 2 23 We shall now utilize each of the results above to find an expression for the central solutions based on the set 2 21 Let us consider the following sets which are obtained from the evaluation of the nt order derivative of set 2 21 Eu Ry C C7 28 p29 ma 0 lt Ru Ay L 2 0 7 T Qu Oz 4 0 and Ov 0z 4 0 9 u 9p Pu du Ow _ T Es TR N C C Oz L 022 M 023 0 0 0 1 u 9 Tou AOV rou Do tbe b 022 Oz A SY 0 0 from 2 23 2 3 SOLUTIONS TO EQUILIBRIUM EQUATIONS 15 T 0 u 02 0 and 0 8z 0 3 3 4 E RE C cru pue v Q9 oz 7 0 3 O23 Oz 0 0 0 2 T Ou Age LT u 070 R Oz x 022 023 022 YY A 0 0 from 2 22 T 0 u 8z 0 and 4 82 0 T 9 n 0u 9z D 0 and 00D Jo MD 0 u ay E Re quon 9D 8z 1 Oz da ozn oO v n 1 u E RT2C Vu A 90709 T My z n 1 Ozn 1 z 0 0 T IM u
35. Oz A op T OOX p oy 5x dA A4 se Gx E 5 Bowie an pt n dA o SZ v dA A T T Z anf AE PA e A Oz Oz A OZ Oz lt A gt Se SS OT dz T OM dz M cT SZ v dA A E _ mT 00X aM roe T M d 3 X F eTM o SZdy dA Noting that ofp TT OX _ pty roo Oz Oz it is possible to further simplify in the following manner OW OTT 06 OM 35 3 x TT ip m of zov dA TT MT JT ox rr x M dy 4 mre Ttdp TT toy J oTSZ y dA A for convenience OTT OM T x oy 4 TTt p TTX 4 MT TTt g ao SZ v dA Oz Oz Oz Oz A _ A ama OT Sp T T ae ox 4 eM HTIt p e osz dA ov Oz Oz A and isolate each of the variation terms in the equation Thus for the arbitrary variation of da we get o 280 dA A where p So The equation above is the definition of the resultant forces acting on the cross section as stated in 2 2 Subsequently for the arbitrary variation of dp and 6x we obtain respectively aM T oT Oz SP Oz 0 2 4 ON THE PROPERTIES OF THE SOLUTIONS 19 The two expressions above are the equilibrium equations or the one dimensional beam equations when only the cross section rigid displacements are considered These equations are typically obtained from simple statitcs and so this result shows the consistency between the one and three dimensional approaches Furthermore based on the results above it
36. User s Manual for BECAS A cross section analysis tool for anisotropic and inhomogeneous beam sections of arbitrary geometry Jos Pedro Blasques Riso DTU National Laboratory for Sustainable Energy Technical University of Denmark Frederiksborgvej 399 P O Box 49 Building 114 DK 4000 Roskilde Denmark P3 jpbl risoe dtu dk RIS R 1785 February 23 2012 Riso DTU National Laboratory for Sustainable Energy ii Title of report User s Manual for BECAS v2 0 a cross section analysis tool for anisotropic and inhomogeneous beam sections of arbitrary geometry Author Jos Pedro Blasques Address Riso National Laboratory for Sustainable Energy Technical University of Denmark Frederiksborgvej 399 P O Box 49 Building 114 DK 4000 Roskilde Denmark E mail jpbl risoe dtu dk Copyright and ownership All rights to this User s Manual belong exclusively to Riso DTU This User s Man ual may only be accessed when the reader has a valid license from Ris DTU to use the BECAS software A license can be obtained from Jos Pedro Blasques at jpbl risoe dtu dk Disclaimer Riso DTU disclaims all responsibility for any kind of damage including loss of profit loss of capital or any caused damage or loss which might appear by use or erroneous use of the BECAS software or Documentation even though Riso DTU should have been informed of the possibilities of such damage ii lv Preface The BEam Cross section Anal
37. VABS BECAS VABS Rel Diff Ks11 7 61E 01 7 61E 01 7 20E 04 Keo 2 93E 01 2 93E 01 6 92E 04 Ks33 2 92E 00 2 92E 00 3 35E 13 Ksas 3 29E 02 3 29E 02 5 88E 08 Kss5 2 94E 01 2 94E 01 1 63E 09 Kse 3 95E 02 3 95E 02 7 20E 04 Kso6 8 26E 03 8 26E 03 6 91E 04 K 35 5 75E 02 5 75E 02 1 55E 11 of the shear and elastic center are presented in Table 4 16 As can be seen there is a very good agreement between BECAS and VABS for all results The warping displacements obtained for different values of the components of y are presented in Figure 4 19 As can be seen the shear torsion coupling term Xs 26 and the extension bending coupling 35 are non zero The effect of these couplings on the warping deformations can be seen for 7 1 and 7 1 respectively This is due to the shift in the shear and elastic center positions due to the asymmetric geometry 4 2 NUMERICAL EXAMPLES 61 Table 4 16 Shear and elastic center positions zs ys and xi yi respectively for three cells cross section T1 Comparison between BECAS and VABS BECAS VABS Diff s 2 823E 02 2 823E 02 0 Ys 0 0 0 z 1 969E 02 1 969E 02 0 Yt 0 0 0 y 0 1 05 x r 0 1 Ne x a Ta 1 b Ty 1 Figure 4 19 Cross section warping displacements for three cells cross section T1 displacements not to scale 62 CHAPTER 4 VALIDATION Three cells cross section of orthotropic top faces and isostropic webs T2 In the last
38. ation 2 24 it is important to note that the central solutions are linear and homogeneous functions of the force resultants Thus it is possible to write u X90 on E ed dd pu 2 33 x MEC rid Inserting the expressions above in 2 24 yields EX RY C C7 9X p2Y RX AY L7 15 ox OY _ RES aie ge 2 34 Note that the set of equations above can be obtained by replacing 0 Ig in 2 24 This corresponds to determining the central solutions for six different choices of the stress resultant 0 in an orderly way i e setting one of the entries to unity and the remaining to zero In fact each of the six columns of X Y aX and ax hold the corresponding displacement solution for each of the different stress resultants An homogeneous function is such that f ox of x In this specific case one important inference is that homogeneous functions do not have an independent term 2 5 CROSS SECTION PROPERTIES 23 Restating the expression for the variation of the total energy obtained from the virtual work principle oW p v w DE dA J o Se dA 0 A The previous equation can be restated as 607 F 0 ode dA A where F is the compliance matrix of the section The strain is redefined by inserting 2 33 in 2 12 e SZY0 BNX0 SN Z which then yields for the internal energy Oz 4 Oz szy BNX sw 0 dA Z The former can be stated in matrix form as in sox E C R
39. cal ex periments is described The material properties are defined and the cross section coordinate systems as well as the cross section constitutive relation are restated Moreover the general organization of the numerical experiments and the aim of each is described The results for the numerical experiments are presented next 4 1 Setup Three material types have been considered two isotropic and one orthotropic ma terial Their stiffness properties are presented in Table 4 1 Furthermore four different cross section geometries have been considered solid square cylinder half cylinder and three cells All the combinations of cross section geometry and material properties are summarized in Table 4 2 Recall from Section 2 2 1 the orientation of the coordinate system presented again here for convenience in Figure 4 1 The non zero entries of the cross section stiffness matrix K the position of the shear and elastic center and the warping displacements are calculated for each of the numerical experiments The entries of K and the position of the shear and elastic center are compared to the results from 37 38 CHAPTER 4 VALIDATION Table 4 1 Material properties for isotropic material 1 and 2 and orthotropic material scaled values for E glass according to Handbook of Composites 25 The factor a Ej Ez shall be used in the study of extremely inhomogeneous sections Material Isotropic 1 Isotropic 2 Orthotropic
40. center m is given as Ne Ne Tm gt sume S 2 e 1 e 1 Ne Ne Um ins S 2 e 1 e 1 where Eme Yme Ve and ge are the coordinates of the centroid the volume and the density of element e respectively and Nne is the number of elements in the cross section mesh 26 CHAPTER 2 THEORY MANUAL 3D finite y element mesh wee re Figure 2 4 Three dimensional finite element mesh created inside commercial finite element package Coordinate system convention and definition of the element width Az 2 7 An alternative formulation based on solid finite el ements An alternative formulation for commercial finite element codes of the theory pre sented before is briefly described here In this case the matrices necessary to solve the sets in 2 34 and 2 35 are evaluated based on the global stiffness matrix of a three dimensional mesh of the cross section using solid finite elements The most important advantage of this approach concerns the possibility of using layered solid elements in the cross section mesh As a result it is possible to decrease the num ber of elements in the cross section and simplify significantly the generation of the cross section finite element model The theory presented here has been originally described by Ghiringhelli and Mantegazza 9 2 7 1 Evaluation of cross section stiffness matrix This formulation assumes that the cross section is generated in a finite element package The cross section
41. center positions zs ys and xi yi respectively for half cylinder cross section C2 Comparison between BECAS and VABS BECAS VABS Diff 90 a 1 206E 01 1 206E 01 0 Ys 0 0 0 a 6 051E 02 6 051E 02 0 yt 0 0 0 c tT 1 d amp 1 Figure 4 12 Cross section warping displacements for half cylinder cross section C2 displacements not to scale 54 CHAPTER 4 VALIDATION Cylinder cross section of two isotropic material half C3 The effect of extreme material inhomogeneity is investigated here following the ap proach described in Section 4 2 1 for the solid square cross section The cylindrical cross section is divided in two One half is made of the isotropic material 1 while the second half of the cylinder cross section is made of the isotropic material 2 The distribution of the two materials in the cross section can be seen in Figure 4 13 The resulting non zero entries of the cross section stiffness matrix are given 0 1 0 05 0 05 04 1 0 0 1 X Figure 4 13 Geometry and finite element mesh of cylinder cross section with two materials C3 isotropic material 1 light and 2 dark in Table 4 12 with respect to the ration E1 E As can be seen there is a very good match between BECAS and VABS As the stiffness of the isotropic material 2 vanishes the entries of the c
42. ction constitutive relation y F 0 the following holds Kz Fs 611x Fs 62Ly F s6aMa Fs65My Fs66Mz 0 Inserting 2 36 into the previous equation yields Fs e1 Fs 60 L 2 Fa 66Ya Te Fs 62 Fs cal L 2 Fo 6623 Ty 0 Since the above has to be valid for any T and Ty F 62 Fs 64 L 2 TEC Fs 66 s Pio Fs 65 L 2 Fs 66 From the previous equation it can be seen that the shear center is not a property of the cross section Instead in the case where the entries Fs 4 and Fs65 associated with the bending twist coupling are not zero the position of the shear center varies linearly along the beam length The expressions for the position of the elastic center can be determined in the same manner The elastic center is defined as the point where a force applied normal to the cross section will produce no bending curvatures i Ke Ky 0 Thus assume that a load T is applied at the point x y in the cross section The moments induced by this force are Mr Tayi My T 24 2 37 We look for the positions y for which kz ky 0 From the cross section constitutive relation Kr F 43T F F 44Mz F 45 My 0 Ky F 531 sts Fs 54Mz Fs 55 My 0 Since the previous must be valid for any force T inserting 2 37 into the previous equation will result in the following set of linear equations Fs 43 Fsaayt Fs asx 0 Fs 53 Fssayt Fs 55 0 2 6
43. ction strains 8 Section strain parameters ssssssse eee nett enn n nes 8 Matrix of finite element interpolation functions sess 9 Nodal degrees of freedom in cross section finite element mesh 9 External work per unit length 0 0 cee cece eee eee eee 9 Internal work or elastic energy per unit length 0 4 10 Total virtual work per unit length 0 00 e cece eee eee eee 12 Number of nodes in the cross section finite element mesh 16 Cross section compliance matrix 0 cece cece eee en 23 Cross section stiffness matrix sess 23 Coordinates of elastic center sss ees 24 Coordinates of shear center 0 ccc cece eee ees 24 Cross section mass matrix 1 0 0 6 cece cece cee s 25 Mass per unit length ssssssessseese eect hn 25 Mass moments of inertia 00 cc cece eee n 25 Cross section mass center 26 6 eh 25 Number of nodes in the cross section finite element mesh 31 Number of elements in the cross section finite element mesh 32 ix Nd CONTENTS Number of d o f in the cross section finite element equations 67 Chapter 1 Introduction This report describes the development and implementation of the BEam Cross section Analysis Software BECAS Cross section analysis tools are commonly employed in the development of beam models for the analysis of long slend
44. d is negligible the so called central solutions 2 4 1 Rigid motions In this section we look only at the displacements v which do not strain the section that is we do not include the warping displacements Thus Xx YPz s v Zr x n Xy 22 2 25 Xz YPz Toy recalling that x x z and qv z correspond to the translation of the cross section reference point and rotations respectively The three dimensional strain components in this case are given by 1 Ovz Ove 0 ae a x y x 1 0v Ou ey 2 2 D 0 w 93 y y 1 0v Ov 2ery HA 0 Exy al 2 7 rm Que De OXx Oy E CAS T Oz Ov Ov Ox Op Bee oil zy _ OXy z Eyz E m De T 8s Ox 1 00 vz Ox Os Oy ena as tg Oz In matrix form the previous equations can be reduced to Ox Op t aA sore PHP a TMG where 0 1 0 t 1 0 0 0 0 0 It is convenient at this point to introduce the vector of strain parameters Y Tas s Ta dens Rys amp T defined in 2 8 and restated here as Ox T rt g Oz The three dimensional strain are hence SZ 2 26 Recalling the expression for the total work has been defined as T OVE e A 1 ode dA 0 2 27 Oz A OZ A 18 CHAPTER 2 THEORY MANUAL where the displacement s in 2 18 has been simply replaced here by v Replacing 2 25 and 2 26 in 2 27 yields O pPT x pTnd gol OND Ox pond sue dA Oz A Oz A T T zi Op OX ga y Dunk otszsw dA A Oz A
45. e between BECAS commercial finite element packages and HAWC2 Riso DTU s own code for the aeroelastic analysis of wind turbines Their contributions are gratefully acknowledged All feedback and suggestions for further improvements and extensions is most welcome Jos Pedro Blasques Roskilde November 2011 vi PREFACE Contents Preface 1 Introduction 1 1 Version history ars geele 2 wf ea mee De eee Roa es Theory manual 2 1 Assumptions S gorse tds en a ee eS lee BOE es Oh Res 2 2 Equilibrium equations a 2 2 1 Basic definitions 2l 2 2 2 Iinemisities s 3 6 vise ge Em de Ged eae te aw EU Ub un Goh Ble o 2 2 3 Strain displacement relation leen 2 2 4 Virtual work principle len 2 3 Solutions to equilibrium equations o 2 3 1 Extremity solutions e e e 2 3 2 Central solutions ee ee 2 3 3 Constraint equations ee 2 4 On the properties of the solutions a aooaa 2 14 Rigid motions e o o a deem ds da ds 2 4 2 Warping displacements 2 5 Cross section properties sn 2 5 1 Cross section stiffness matrix o ooa 2 5 2 Shear center and elastic center positions 2 6 Cross section mass matrix 2 a 2 7 An alternative formulation based on solid finite elements 2 7 1 Evaluation of cross section stiffness matrix Implementation manual 3
46. e non zero 0 Z 0 central solutions see 6 Owing to the structure of matrices M E and H and mostly due to the fact that H is skew symmetric for the extremity solutions the eigenvalues will be complex and come in pairs That is to each eigenvalue A a ib corresponds a second A ib The solution corresponding to the first eigenvalue decays while z increases moving away from the first end of the beam wheras the solution associated with the second eigenvalue in the pair is identical but decays as z decreases moving away from the second end of the beam This solutions can be used to study end effects and determining the diffusion length or the distance at which the effects from the external loads become negligble e g see Horgan 18 and Choi and Horgan 19 for a discussion on the diffusion length in anisotropic elasticity The central solutions for which 0 0 on the other hand do not have any exponential decay as they correspond to eigenvalues A 0 with multiplicity p These are the solutions at the central part of the beam sufficiently away from the extremities so that the effect from the external loads is negligible on the stress field The solutions for problems where A 0 and its multiplicity p gt 2 are given as see 21 for a comprehensive presentation of this topic uj z hie u j z hye haooze AES Ajz Aiz pj 1 Ajz Up j Z hpj1e 7 hpaze 9 275 Nyy where the correspo
47. ed A M S11 S22 2 Si2 S21 1 C 5 S11 S22 S12 821 1 E S11 822 Si S21 Az Note that deriving the original equations in Ghiringhelli and Mantegazza 9 yields the 1 2 factor in matrix C This factor is not present in the derviation by Ghiringhelli and Mantegazza 9 The remaining matrices can be determined as R CZa L MZa A ZLMZc where Za is defined as OS EN Ra S o l R O eu e 100 0 O Una 0 0 0 ns 0 O 1 ya Tn 0 where n is the number of nodes in the reference cross section face i e half the nodes in the whole three dimensional model Having defined all the matrices it is possible now to solve the sets 2 34 while remembering to account for the matrix of constraint equations D presented in 2 3 3 Replacing the solutions of 2 34 into 2 35 it is straight forward to evaluate the cross section stiffness matrix K The formulation presented here has been implemented in the class of Matlab functions BECAS 3D 28 CHAPTER 2 THEORY MANUAL Chapter 3 Implementation manual The theory presented in the previous sections up to the determination of the cross section stiffness matrix and shear and elastic centers is implemented in the BEam Cross section Analysis Software BECAS This section addresses the practical im plementation of the theory The MATLAB implementation of BECAS described in Appendix 5 is according to the expressions de
48. er structures These type of models can be very versatile when compared against its equivalent counterparts as they generally offer a very good compromise between accuracy and computational efficiency When suited beam models can be advantageously used in an optimal design context see e g Ganguli and Chopra 1 Li et al 2 Blasques and Stolpe 3 or in the development of complex multiphysics codes Wind turbine aeroelastic codes for example com monly rely on these types of models for the representation of most parts of the wind turbine from the tower to the blades see e g Hansen et al 4 Chaviaropoulos et al 5 In specific the development of beam models which correctly describe the behaviour of the wind turbine blades have been the focus of many investigations The estimation of the properties of these types of structures becomes more complex as the use of different combinations of advanced materials becomes a standard It is therefore paramount to develop cross section analysis tools which can correctly account for all geometrical and material effects BECAS is a general purpose cross section analysis tool specifically developed for these types of applications BECAS is able to handle a large range of arbitrary section geometries and correctly predict the effects of inhomogeneous material distribution and anisotropy Based on a def inition of the cross section geometry and material distribution BECAS is able to determine the cro
49. erial mechanical properties given with respect to the material coordinate system have been are defined as Ej the Young modulus of material the 1 direction E the Young modulus of material the 2 direction E33 the Young modulus of material the 3 direction Gj the shear modulus in the 12 plane G3 the shear modulus in the 13 plane Gs the shear modulus in the 23 plane 112 the Poisson s ratio in the 12 plane w4a the Poisson s ratio in the 13 plane 5 2 LIST OF FUNCTIONS AND OUTPUT 67 Va the Poisson s ratio in the 23 plane o the material density The rotation of the material constitutive tensor is described in Section 3 2 5 1 1 Example AII the files necessary to run all the examples presented for the Validation in Section 4 are distributed together with BECAS These files can be used as a starting point for the development of new input for BECAS 5 2 List of functions and output The different functions included in the BECAS library and corresponding output are described in this chapter 5 2 1 BECAS Utils The function BECAS Utils is used to build working arrays It is included inside all the other functions in the 2D version of BECAS It is called using utils BECAS Utils nl 2d el 2d emat matprops 5 2 2 BECAS Constitutive Ks The function BECAS Constitutive Ks is the BECAS central function and is used primarily for the evaluation of the cross section stiffness ma
50. ernative implementation of BECAS using solid finite elements is described here The main advantage of this approach concerns the possibility of using layered solid finite elements All BECAS functions associated with this implementation start by BECAS_3D 5 3 1 Input The input to the BECAS_3D group of functions is e k3d ng x 3 array with the sparse form of stiffness matrix of the cross section slice meshed with solid finite elements where ng is the number of degrees of freedom in the solid finite element model Each row in the k3d array is in the form row number column number value where row and column number correspond to the row and column positions in the global stiffness matrix e nl 3d Nn 34 X 3 array with the list of nodal positions where each row is in the form node number x coordinate y coordinate z coordinate where np 34 is the total number of nodes in the solid finite element model The node numbering needs to be in the same order as it is listed in the k3d matrix 5 3 2 List of functions and output The following functions are part of the BECAS_3D group of functions 5 3 3 BECAS_3D_Utils The function BECAS_3D_Utils is used to build working arrays It is included inside all the other functions in the 3D version of BECAS It is called using utils BECAS_3D_Utils k3d n3qd 72 CHAPTER 5 USER S MANUAL 5 3 4 BECAS_3D_Constitutive_Ks This is the main function of the BECAS_3d and is used to determine the c
51. ess properties The value of the non zero entries in the cross section stiff 0 1 0 05 gt 0 0 05 04 0 05 0 X Figure 4 11 Geometry and finite element mesh of half cylinder cross section C2 ness matrix as estimated by both BECAS and VABS are presented in Table 4 10 As can be seen there is a very good agreement between the two methods The same is observed in the estimated positions of the elastic and shear center as presented in Table 4 11 Finally the warping displacements for different values of the strain parameters Y are presented in Figure 4 12 Note that the magnitude of the warping displacements for k 1 torsional curvature are of the order of O 10 Hence it is expected that the resulting strains and stresses will also be very large also This result serves to illustrate the fact that open cross sections are weaker in torsion then closed cross sections Table 4 10 Non zero entries of cross section stiffness matrix for half cylinder cross section C2 Comparison between BECAS and VABS BECAS VABS Rel Diff Koi 4 964E 02 4 964E 02 7 195E 04 Kyo 6 244E 02 6 244E 02 7 200E 04 K ss 2 982E 01 2 982E 01 0 000E 00 K 44 1 349E 03 1 349E 03 4 672E 09 K ss 1 349E 03 1 349E 03 4 637E 09 K og 9 120E 04 9 120E 04 7 200E 04 K ss 1 805E 02 1 805E 02 8 862E 12 Ks26 7 529E 03 7 529E 03 7 200E 04 4 2 NUMERICAL EXAMPLES 53 Table 4 11 Shear and elastic
52. following have to be defined e RadialPosition coordinate of the section along the span of the blade e DutputFilename name of file to which the output is written e g OutputFilename BECAS2HAWC2 out The function is called as BECAS Becas2Hawc2 0utputFilename RadialPosition Ks Ms ShearX ShearY ElasticX ElasticY MassX MassY AlphaPrincipleAxis ElasticCenter AreaX AreaY AreaTotal Axx Ayy Axy nl 2d el 2d emat matprops 5 2 8 BECAS TransformMat The function BECAS TransformMat is used to determine the translated and rotated values of the cross section constitutive stiffness or mass matrices The function requires extra input namely e p column array specifying the coordinates of the new reference point e g p ShearX ShearY to translate the matrix to the shear center 5 3 THE BECAS 3D IMPLEMENTATION 71 e alpha angle of rotation around the z axis in degrees defined positive in the counter clockwise direction e g alpha AlphaPrincipleAxis to align with the elastic axis The function is called as M BECAS_TransformationMat M p alpha 5 2 9 Examples All the files required for the replication of the validation examples are distributed together with BECAS The user should specify the name of the corresponding folder in the file Inputdata4RunMe m The file RunMe m which calls BECAS should then be ran to obtain the results 5 3 The BECAS 3D implementation An alt
53. for half cylinder 56 CHAPTER 4 VALIDATION Table 4 13 Shear and elastic center positions z Ys and zi yz respectively for cylinder cross section C3 with respect to the ration E1 E Comparison between BECAS and VABS BECAS VABS En E Ts Ys Ts Ys 1E 00 0 000E 00 0 000E 00 0 000E 00 0 000E 00 1E 01 9 866E 02 0 000E 00 9 866E 02 1 069E 15 1E 02 1 182E 01 0 000E 00 1 182E 01 1 415E 15 1E 03 1 203E 01 0 000E 00 1 203E 01 1 292E 15 1E 04 1 206E 01 0 000E 00 1 20G6E 01 3 740E 16 1E 05 1 206E 01 0 000E 00 1 206E 01 1 308E 15 E E Tt ut Tt Yt 1E 00 0 000E 00 0 000E 00 0 000E 00 0 000E 00 1E 01 4 951E 02 0 000E 00 4 951E 02 3 279E 17 1E 02 5 931E 02 0 000E 00 5 931E 02 1 920E 17 1E 03 6 039E 02 0 000E 00 6 039E 02 3 411E 18 1E 04 6 050E 02 0 000E 00 6 050E 02 4 960E 18 1E 05 6 051E 02 0 000E 00 6 051E 02 1 255E 18 42 NUMERICAL EXAMPLES 57 0 01 0 01 0 05 0 05 0 05 Figure 4 15 Cross section warping displacements for cylinder cross section C3 with E E2 1000 displacements not to scale 58 CHAPTER 4 VALIDATION Cylinder cross section of two isotropic materials layered C4 In the final numerical example using the cylinder cross section we assume a layered structure through the thickness The material distribution geometry and cross sec tion
54. g the fibers in a layer FiberPlaneRotationBECAS mw used to generate the code for the function which rotates the material constitutive matrix Q when orienting the fiber plane in an element in the MATLAB function BECAS_Utils Input In order to run BECAS the following input is necessary e nl 2d n x 3 array with the list of nodal positions where each row is in the form node number x coordinate y coordinate where nnp is the total number of nodes The node numbering need not be in any specific order el 2d ne x 8 array with the element connectivity table where each row is in the form element number node 1 node 2 node 3 node 4 node 5 node 6 node 7 node 8 where n is the total number of elements The element numbering need not be in any specific order The value of node 5 through node 6 has to be zero for Quad4 element to be used Otherwise Quad8 is automatically chosen emat n x 4 array with element material properties assignment where each row is in the form element number material number fiber angle fiber plane angle where ne is the total number of elements The element numbering need not be in any specific order The material number corresponds to the materials assigned in the matprops array e matprops nmat x 10 array with the material properties where each row is in the form E11 E22 E33 Giz G13 G23 V12 via V23 0 where nmat is the total number of different materials considered The mat
55. iffness as estimated by both BECAS and VABS decrease past the value obtained for half the section and have a local minima at E E2 10 Chen et al in 15 present a similar result using a different cross 4 2 NUMERICAL EXAMPLES 45 Table 4 6 Shear and elastic center positions s Ys and xi yz respectively for square cross section S2 with respect to the E E ration Comparison between BECAS and VABS BECAS VABS E Ez Ts Ys Ts Ys 1E 00 0 000E 00 0 000E 00 0 000E 00 0 000E 00 1E 01 2 045E 02 0 000E 00 2 045E 02 2 306E 17 1E 02 2 450E 02 0 000E 00 2 450E 02 5 817E 17 1E 03 2 495E 02 0 000E 00 2 495E 02 1 385E 16 1E 04 2 500E 02 0 000E 00 2 500E 02 8 010E 17 1E 05 2 500E 02 0 000E 00 2 500E 02 1 683E 16 Ey E gt Xt Ut Tt Yt 1E 00 0 000E 00 0 000E 00 0 000E 00 0 000E 00 1E 01 2 045E 02 0 000E 00 0 000E 00 0 000E 00 1E 02 2 450E 02 0 000E 00 2 450E 02 2 596E 19 1E 03 2 495E 02 0 000E 00 2 495E 02 1 266E 19 1E 04 2 500E 02 0 000E 00 2 500E 02 4 118E 19 1E 05 2 500E 02 0 000E 00 2 500E 02 9 174E 19 0 4 T e BECAS O VABS 0 344 Half square 0 3 0 25 we um S 8 18 0 15F 7 a E 0 1F 0 05rF 0 1 1 1 1 10 10 10 10 10 10 Log E E Figure 4 5 Variation of K 11 entry of the cross section stiffness matrix with respect to the ration Ej E gt for the square cross section S1 Results for BECAS and VABS compared with the solutio
56. ing defined the cross section mesh and material properties the subsequent step concerns the derivation of each of the matrices in Equation 2 16 The implementation is based on four or eight node isoparametric elements The node numbering and isoparametric coordinate system are presented in Figure 3 2 The shape functions employed in the derivation of the four node isoparametric finite 29 30 CHAPTER 3 IMPLEMENTATION MANUAL Figure 3 1 Example of the two dimensional finite element mesh of a generic wind turbine section using four node isoparametric finite elements Four node element Eight node element Figure 3 2 Isoparametric coordinate system nodal positions and position of Gauss points for the four node isoparametric plane finite element element are G n M 69 i0 9 0 9 MEn 1 0v Ns E 7 1 8 14m Na 60 Ble A In the case of the eight node isoparametric finite element the shape functions are 1 1 1 M 69 08 0 9 06 N No En 5 1 8 11 5 N5 No Ns En 5 1 8 1 41 5 No Nr N 69 4 1 8 1 n 07 No N 69 50 8 0 9 Non 5048 1 07 Nr En 7 1 1 n Nm 35 1 0 15 The position of a point in the element is given by interpolation of the nodal positions as i 1 gel i l 3 1 TWO DIMENSIONAL FINITE ELEMENT ANALYSIS 31 where ny is the number of nodes in the element and zi yi zi are the nodal posi tion
57. is represented by a three dimensional slice meshed using solid finite elements The number of nodes in the faces of the elements facing the cross section plane is not restricted However along the length it is required that exactly two Gauss points only are used in the length direction The choice of the slice thickness Az should be such that the resulting solid finite elements are not too distorted as to affect the quality of the results see Figure 2 4 Thus a length of the order of the average side length of the two dimensional finite elements is recom mended Finally the preferred commercial finite element package should necessarily allow the user to access and manipulate the global finite element stiffness matrix The finite element equilibrium equations for the three dimensional model are Sw f where S is the finite element stiffness matrix w the displacement vector and f the external load vector The following decomposition of the finite element stiffness matrix S and corresponding displacement and load vector is allowed s ls S21 S82 Wa f 2 7 ANALTERNATIVE FORMULATION BASED ON SOLID FINITE ELEMENTS27 where the index 1 and 2 refer to the contributions from the nodes at z 0 and z Az respectively Based on the above submatrices Ghiringhelli and Mantegazza 9 have derived the expressions for the matrices necessary for the computation of the cross section stiffness matrix Hence after proper derivation the following are defin
58. isotropic materials in the sixth example C3 Much like in the case of the solid square cross section 2 the aim here is to validate the results of BECAS when studying thin walled cross sections with extreme ma terial inhomgeneity A similar procedure is adopted in the seventh example where a layered type of structure is assumed through the thickness of the cylinder C4 In the eight example a three cell isotropic cross section is considered T1 The aim is to validate BECAS for the analysis of multi celled thin walled cross sections The final example assumes that the three cell cross section is made of isotropic and orthotropic materials T2 The aim in this case is to validate the results from BECAS for thin walled multi celled closed cross sections with anisotropic material properties 4 2 NUMERICAL EXAMPLES Al 4 2 1 Square The dimensions of the solid square cross section are given in Table 4 3 The cross section is presented in Figure 4 2 Table 4 3 Geometrical dimensions of solid square beam Width W 0 1m Height H 0 1 m 0 05 0 005 0 0 05 Xx Figure 4 2 Geometry and finite element mesh of square cross section with one material 42 CHAPTER 4 VALIDATION Square cross section of isotropic material S1 In this case the solid square cross section is made of isotropic material 1 The resulting non zero entries of the cross section stiffness
59. l 6 ov ZT St y uT NTBT 492 NTST Q Oz A f z szw BNu SZ dA Oz Sp ZS QSZy dA A f Sp Z S QBNu dA A f Sul ZTSTQs N dA A Z f du NTBTQSZ4 dA A I u NTB QBNu dA A i dul NT BTQsn24 dA A Oz T I 5 2 NTSTQSZV dA A Oz T ef i 52 NTSTQBNu dA A Oz du A e T eb i 5 2 _NTsTQSN A Oz where the following matrices can be identified A Z SY QSZ dA A r N B QSZ dA A E NT BT QBN dA 2 16 A C a NTSTQBN dA A L f NTSTQSZ dA A M N7S QSN dA A Hence rearranging the state variables it is possible to write the internal virtual work in matrix form as TM C L du Oz Oz pon t du CT E R u 2 17 d dap LT RT A m 12 CHAPTER 2 THEORY MANUAL Total virtual work According to Equation 2 13 the total virtual work per unit length can be written as Oz Oz For a general virtual displacement ds and virtual strain de a necessary and sufficient equilibrium condition is 9 ds Tod je jos e 2 y da seo dA 2 18 Oz A Oz A oW 0 Oz Thus inserting Equation 2 14 and 2 17 into 2 18 the following relation is ob tained du 1T a Ons a C L Oe Ju C E R u sy L R A ap Internal virtual work of the beam slice 8 T S 7 y 0 External virtual work of the beam slice ort 5 170 z Equilibrium of the beam slice The Previous equation must be true for any admissible virtual displacement
60. lane coordinate system z y z c Convention adopted for the rotation of the fiber plane coordinate system z y z into the material coordinate system 123 Moreover note that IH jg Bu us where rj is the Poisson s ration that characterizes the transverse strain in the j direction when the material is stressed in the 7 direction The natural strains e are related to the engineering strains yi by yo2 22 33 33 11 11 713 2 13 12 2 12 V23 2e23 3 2 2 Rotation The next two steps concern the two rotations of the material constitutive matrix necessary to determine the element material constitutive matrix in the fiber coor dinate system The three element coordinate systems element fiber plane and fiber and respective conventions for each of the rotations are defined in Figure 3 3 The material constitutive matrix is rotated first into the fiber plane coordinate system and subsequently into the fiber coordinate system Each of the rotations is performed following the procedure described next The following approach for the rotation of the material constitutive matrix is based on Reddy 24 We consider the relationship between the stress components in a material m and a problem p coordinate systems In tensor format the stress tensor is transformed as Gis lidia mnn alo lmilni Cij m where are the direction cosines defined as lij ei m gt ej p 3 2 MATERIAL CONSTITUTIVE
61. n for half the section section but consider only one extreme value of the stiffness ratio To the authors best knowledge there is no study reported in the literature where a similar study is performed Finally the cross section warping displacements are presented for different values of w in Figure 4 6 considering E E2 1000 As the stiffness of material 2 46 CHAPTER 4 VALIDATION vanishes the extension bending and shear torsion coupling terms arise Ky26 and K 35 in Table 4 5 These coupling effects are visible in the warping displacements For 7 1 the distortion induced by the shear strain is superimposed by a torsion type of deformation pattern When 7 1 the tension strain induces a bending type of in plane deformation visible in the arched form of the cross section Figure 4 6 Cross section warping displacements for square cross section 82 where E4 E 1000 displacements not to scale 4 2 NUMERICAL EXAMPLES 47 Square cross section of orthotropic material S3 In the last example it is considered that the solid square cross section is made of layered orthotropic material The fiber plane lies parallel to the xz plane cf Figure 4 1 and rotate around the y axis The variation of the magnitude of the non zero entries in the cross section stiffness matrix with respect to the fiber orientation are presented in Table 4 7 As can be seen there is a very good agreement between the Table 4 7 N
62. nce study has not been performed Provided the same finite element mesh is used both tools will give the same results as can be seen next 40 CHAPTER 4 VALIDATION 4 2 Numerical examples All numerical experiments conducted to validate the current implementation of BE CAS are presented in this section The first results are presented for cross section S1 This first example illustrates the ability of BECAS to handle solid square cross section In the second case the same solid square cross section geometry is analyzed although in this case it is made of two different materials S2 The aim is to analyze the behaviour of BECAS when handling cross sections made of different materials with a high contrast between stiffnesses In the third case the solid square cross section is made of orthotropic material S3 The objective is to validate the effect of material anisotropy on the cross section stiffness properties estimated by BECAS In particular we look at the estimated coupling terms arising from the ma terial anisotropy In the fourth case case the cylinder cross section made of isotropic materials is analyzed C1 The aim is to analyze the behaviour of BECAS when dealing with thin walled cross sections In the fifth example only half the cylinder is modelled C2 This experiment serves to validate the behaviour of BECAS when handling open thin walled cross sections The cylinder cross section is then divided in two and made of two different
63. nding solution in this case is u z c hie Ca ha1 zh gut SEE Cp hoi zha 20 oos git The authors would like to express its gratitude to Assoc Prof Mads Peter Si irensen DTU MAT for all the help unravelling this step of the derivation 22 CHAPTER 2 THEORY MANUAL In our specific case where A 0 the solutions in this case are of the polynomial type uj z cC1h 1 co ha1 zh Toe Cp ho zhoo 4 2855 A central solution u will be any linear combination of the solutions uj for which Aj 0 If n is the maximum of all p then Oue Oz 2 32 Thus the central solutions are polynomial functions in z of at most degree n 2 5 Cross section properties We look first for the compliance matrix of a cross section of the beam That is we are interested in finding an expression for the strain energy as function of the stress resultants and moments These are in fact the central solutions derived above a particular solution depending on the applied section forces at a given cross section subject to particular boundary conditions which guarantee that the effects of the extremity solutions are negligible A procedure is presented in this section for the practical determination of the cross section stiffness matrix based on the result of the central solutions Finally expressions for the determination of the shear and elastic center are derived 2 5 1 Cross section stiffness matrix From Equ
64. ness and strength properties of anisotropic and inhomogeneous beam cross sections According to Yu et al 13 implementations of this theory have been in fact used as a benchmark for the vali dation of any new tool emerging since the early 1980 s see e g Yu et al 13 14 and Chen et al 15 Many other cross section analysis tools have been described in the litterature The reader is referred to Jung et al 16 and Volovoi et al 17 for an assessment of different cross section analysis tools Nonetheless at this stage the Variational Asymptotic Beam Section analysis commercial package VABS by Yu et al 13 is perhaps the state of the art for these type of tools VABS has been extensively validated see Yu et al 13 14 Chen et al 15 and is therefore used in this report as the benchmark for the validation of BECAS As shall be seen the cross section properties estimated by both tools are in very good agreement The theory presented in this report concerns only the determination of the cross section stiffness properties for inhomogeneous and anisotropic beam cross sections of arbitrary geometry i e the theory implemented in BECAS Most of the relevant information which is spread across the different publications namely 6 12 and which concerns the estimation of the cross section stiffness properties is compiled here The aim was to produce a self contained document which can serve as a developer s manual for the readers wishi
65. ng to use understand and further develop BECAS This report is organized as follows Chapter 2 Theory Manual All the theory leading to the evaluation of the cross section stiffness properties is presented in this chapter The assumptions un derlying the presented theory are stated first in Section 2 1 The equibilibrium equations are established next in Section 2 2 and consequently resolved in Sec tion 2 3 Some of the mathematical properties invoked in the resolution of the equilibrium equations are described in detail in Section 2 4 Finally the expressions for the cross section stiffness matrix and positions of shear and elastic centers are determined in Section 2 5 Chapter 3 Implementation Manual The details concerning the numerical im plementation of the theory are presented in this chapter A two dimensional implementation based on four node plane finite elements is presented in Sec tion 3 1 Furthermore an implementation of the method for commercial finite element codes is described next in Section 2 7 Finally in Section 3 2 the constitutive matrix is defined and some important conventions utilized in its transformation are stated Chapter 4 Validation All numerical experiments performed for the validation of VABS are presented in this Chapter The general setup for the numerical experiments is described first in Section 4 1 The validation results obtained for the different cross sections are finally presented in Section 4 2
66. ntly established The solution to the equilibrium equations a set of second order linear differential equations is discussed in Section 2 3 As shall be seen the solution is defined by a particular integral which depends on the boundary conditions or internal force resultants in this case and a general integral which resolves into an eigenvalue problem The particular integral corresponds to solutions far from the ends of the beam where the end effects are negligible the central solutions while the general integrals corresponds to the solutions at the extremities of the beam are applied extremity solutions nomencalture according to Giavotto et al 6 At this point some mathematical properties of the solutions are invoked which are only detailed later in Section 2 4 The reader may wish to avoid this section if only a general overview of the method is required The equations for the cross section stiffness matrix are presented in Section 2 5 Based on the cross section stiffness properties it is possible to compute the positions of the shear and elastic center 2 1 Assumptions The theory presented in the next sections is valid for long slender structures which present a certain level of geometric and structural continuity Thus there should not be abrupt variations of the cross section geometry and material properties along the beam length Moreover the same should be valid for the loads applied Conse quently the gradien
67. numerical experiment the top faces of the three cells cross section are laminated using the orthotropic material oriented at 45 The vertical faces are made of isotropic material The material distribution is visible in Figure 4 20 The aim is 0 1 gt 0 0 1 0 5 0 0 5 X Figure 4 20 Geometry and finite element mesh of three cells cross section with orthotropic material at 45 in the top faces dark and isotropic material in the vertical shear webs light T2 to combine in one experiment the effects of material anisotropy and inhomogeneity using a closed thin walled multi cell cross section The non zero entries of the cross section stiffness matrix as estimated by BECAS and VABS are presented in Table 4 17 As can be seen there is a very good agreement between both methods The resulting positions of the shear and elastic center are presented in Table 4 18 Like before both BECAS and VABS present a very good agreement in the estimation of the shear and elastic center Table 4 17 Non zero entries of cross section stiffness matrix for three cells cross section with orthotropic material at 45 in the top faces and isotropic material in the vertical webs T2 BECAS VABS Rel Diff Ks11 1 79E 00 1 79E 00 1 36E 10 Kyo 3 22E 01 3 22E 01 2 13E 09 K 33 437E 00 4 37E 00 2 98E 11 Ksas 5 23E 02 5 23E 02 2 16E 10 K 55 3 81E 01 3 81E 01 1 21E 10 Kse6 7 73E 02 7 73E 02 3 29E
68. on S1 dis placements not to scale 44 CHAPTER 4 VALIDATION Square cross section of two isotropic materials S2 The solid square geometry is now divided in two The section geometry and ma terial distribution are presented in Figure 4 4 The mechanical properties of the 0 05 005 0 0 05 Figure 4 4 Geometry and finite element mesh of square cross section with two materials S2 isotropic material 1 light and 2 dark isotropic material 2 are obtained by simply dividing the mechanicals properties of the isotropic material 1 see Table 4 1 by a factor a E E The variation of the value of the non zero entries of the cross section stiffness matrix with respect to the stiffness ration E1 E is presented in Table 4 5 The estimated positions of the shear and elastic centers are presented in Table 4 6 Once again there is a very good agreement between the results from BECAS and VABS As the stiffness Table 4 5 Non zero entries of cross section stiffness matrix for square cross section S2 with respect to E Ez ration Comparison between BECAS and VABS BECAS E1 Ez2 Ks 11 Ks 22 Ks 33 Ks 44 Ks 55 Ks 66 Ks 26 K5 35 E 00 1 00E 00 3 49E 0 3 49E 0 5 91E 04 8 34E 04 0 00E 00 0 00E 00 E 01 5 50E 0 1 28E 0 1 92E 0 2 77E 04 4 59E 04 3 93E 03 13E 02 E 02 5 05E 0 1 38E 0 1 77E 0 2 35E 04 4 21E 04 4 33E 03
69. on zero entries of cross section stiffness matrix for square cross section S3 with respect to fiber orientation fiber plane parallel to xz plane Comparison between BECAS and VABS BECAS VABS Rel Diff BECAS VABS Diff 0 22 5 Ks11 5 039E 01 5 039E 01 2 129E 10 7 598E 01 7 598E 01 1 627E 10 Ks 22 4 201E 01 4 201E 01 2 123E 10 4 129E 01 4 129E 01 1 632E 10 Ks 33 4 800E 00 4 800E 00 2 035E 13 3 435E 00 3 435E 00 7 482E 11 Ks 44 4 001E 03 4 001E 03 6 634E 10 2 489E 03 2 489E 03 7 052E 10 Ks 55 4 001E 03 4 001E 03 6 686E 10 2 27AE 03 2 27AE 03 9 388E 10 K 5 66 7 737E 04 7 737E 04 9 584E 10 9 499E 04 9 499E 04 9 845E 10 Ks 13 0 000E 00 0 000E4 00 0 000E 00 7 387E 01 7 387E 01 2 335E 10 Ks 46 0 000E 00 0 000E4 00 0 000E 00 4 613E 04 4 613E 04 9 910E 10 45 67 5 Ks11 8 421E 01 8 421E 01 3 985E 10 6 039E 01 6 039E 01 2 663E 10 Ks 22 4 473E 01 4 473E 01 1 616E 10 4 883E 01 4 883E 01 1 954E 10 K 5 33 1 713E 00 1 713E 00 1 057E 10 1 241E4 00 1 241E 00 4 043E 12 Ks 44 1 326E 03 1 326E 03 6 825E 10 1 032E 03 1 032E 03 6 659E 10 Ks 55 1 274E 03 1 274E 03 9 726E 10 1 030E 03 1 030E 03 6 778E 10 Ks 66 1 018E 03 1 018E 03 9 663E 10 9 171E 04 9 171E 04 9 804E 10 Ks 13 4 017E 01 4 017E 01 6 091E 10 6 317E 02 6 317E 02 4 366E 10 Ks 46 2 422E 04 2 422E 04 1 042E 09 4 786E 05 4 786E 05 1 023E 09 90 Ks 11 5 0202E 01 5 0202E 01 1 9581E 10 Ks 22 5 0406E 01 5 0406E 01 2 1603E 10 Ksg33 1 2000E 00
70. re the ay curvatures around x and y respectively and k de z is the torsion term 2 2 EQUILIBRIUM EQUATIONS 9 Finite element discretization The warping displacements g are discretized as g x y 2 m N z y u zi yi z 2 10 where N are the typical finite element shape functions and u the nodal warping displacements Note that the latter depend on the position along the beam axis although the shape functions are only defined in the plane of the section The displacement of a point in the cross section is then given as s Zr Nu 2 11 and finally introducing 2 10 in 2 9 yields SZ BNu sn 2 12 z 2 2 4 Virtual work principle The total virtual work per unit length W is given as W We Wi The first variation of the total virtual work per unit length can be written as 2 1 9 n o Oz T9 Oz 219 where W is the work done by the internal elastic forces and W the work done by the external forces acting on the cross section External virtual work Assuming that the surface and the body forces are zero the axial derivative of the work produced by the section stresses is the only contribution to the external work W Thus cf Figure 2 2 ow 0 Ss p dA jeg A OZ g The external virtual work expression can be obtained by replacing the displacement s as defined in 2 4 into the equation above ow _ f v p dg p i ie ALL dA SI dA
71. ross section stiffness matrix K It is called using Ks dX dY X Y BECAS 3D Constitutive Ks k3d n3d d The output is e Ks 6 x 6 array storing the cross section stiffness matrix Ks e dX nq x 6 array X e dY 6 x 6 array ay e X nq x 6 array X e Y 6 x 6 array Y where ng Ny x 3 is the number of degrees of freedom in cross section finite element equations 5 3 5 BECAS 3D CrossSectionProps The function BECAS 3D CrossSectionProps returns some relevant cross section properties It presumes that the cross section stiffness matrix K has been pre viously determined The function is called through ShearX ShearY ElasticX ElasticY AlphaPrincipleAxis BECAS 3D CrossSectionProps Ks The output is e ShearX the x position of the cross section shear center ShearY the ys position of the cross section shear center ElasticX the x position of the cross section elastic center e ElasticY the y position of the cross section elastic center MassX the m position of the cross section mass center e MassY the ym position of the cross section mass center MassPerUnitLength the cross section mass per unit length e AlphaPrincipleAxis the orientation of the cross section elastic Bibliography 10 11 12 Ganguli R Chopra I Aeroelastic optimization of a helicopter rotor with com posite coupling Journal of Aircraft 32 6 1326 1334 1995 Li L Volovoi V
72. ross section stiffness matrix converge to those which were obtained when only half the cross section was considered see Section 4 2 2 Note once again the behaviour of the 11 entry with respect to the ration E1 E plotted in Figure 4 14 Just like in the case of the solid square cross section in Sec tion 4 2 1 the shear stiffness seems to decrease past the half cylinder values having a minima at E E3 10 Finally the variation of the positions of the shear and elastic centers with respect to the ration E E are presented in Table 4 13 There is a very good agreement between BECAS and VABS Also note that the positions of the shear and elastic center do not explain the behavior of the entry The warping deformations for E 1000 are presented in Figure 4 15 The warping deformations in this case are also relatively large and in the same order of magnitude as the half cylinder cross section 4 2 NUMERICAL EXAMPLES 55 Table 4 12 Non zero entries of cross section stiffness matrix for cylinder cross section C3 with respect to E4 E Comparison between BECAS and VABS BECAS E E2 Ks 11 Ks 22 Ks 33 Ks 44 Ks 55 K5 66 K5 26 Ks 35 E 00 1 25E 01 1 25E 01 5 96E 0 2 70E 03 2 25E 03 1 90E 18 0 00E 00 E 01 3 99E 02 6 87E 02 3 28E 0 48E 03 1 08E 03 6 78E 03 62E 02 E 02 3 75E 02 6 31E 02 3 01E 0 36E 03 9 29E 04 7 45E 03 79E 02 E 03 4 74E 02 6 25E 02 2 99E 0 35E 03 9 14E 04 7 52E 03
73. s In matrix form for the four node element N Mb Ne Ni Nils and for the eight node element N MI Nl NI Nil Nslz Nslz Nzi Nels The integration is performed with respect to the element coordinate system although the integrals in 2 16 are defined with respect to the cross section coordinate system To account for the change of coordinates we employ the following transformation gt 2 E 3 OC Oz where Or Oy Oz dE J 2 u a n On On Ox Oy Oz OC Oi a is the Jacobian matrix In this specific case Ox Oy m J 9 o 3 1 0 0 1 Using 2 6 and finding J from 3 1 it is possible to define the strain operator B as o B B B where 1 da 0 0 Ju 0 0 0 Ja 0 0 Ja 0 iE da B J y Ju 0 B Jo3 Ji 0 d 0 0 Jg a D 0 J O 0 J U 0 J 0 0 0 0 0 0 The strain operator is then applied in the derivation of the matrix product ON ON BN B B amp t T n On where ON N aN aN aN 2 eL S Fer Ber ON On The remaining matrix products SZ and SN are obtained by simple matrix multi plication ON ON ON3 ON4 Ung Nat SM Ael 32 CHAPTER 3 IMPLEMENTATION MANUAL 3 1 2 Local and global finite element matrices The integration is performed at the element level and hence the element matrices are evaluated as A Z s qsz axdy f Z S Qsz 7 d dn n f N B Qsz dxdy f NTBTQSZ JI d dy E A NTBTQBN dxdy f NTBTQBN J d dy
74. s section 83 with fibers oriented at 45 Out of plane left and in plane right deformation displace ments not to scale 42 NUMERICAL EXAMPLES 49 a kz 1 b kz 1 0 02 0 04 0 02 0 0 02 0 04 x d ky 1 0 01 0 02 0 03 0 04 0 06 0 04 e amp 1 f hz 1 Figure 4 8 continuation Cross section warping displacements for square cross section 83 with fibers oriented at 45 Out of plane left and in plane right deformation displacements not to scale 50 CHAPTER 4 VALIDATION 4 2 2 Cylinder In this section results are presented based on the cylinder cross section The geo metrical dimensions of the cylindrical cross section are presented in Table 4 8 and the geometry and finite element mesh are presented in Figure 4 8 Table 4 8 Geometrical dimensions of cylinder section Outer radius R 0 1m Thickness t 0 01 m 0 1 0 05 gt 0 0 05 047 0 0 1 Xx Figure 4 9 Geometry and finite element mesh of cylinder cross section with one material C1 Cylinder cross section of isotropic material C1 In the first case it is assumed that the cylinder is made of isotropic material 1 The non zero stiffness entries for both BECAS and VABS are presented in Table 4 9 As can be seen there is a very good agreement between the two approaches The warping displacements are presented
75. sX the m position of the cross section mass center MassY the ym position of the cross section mass center MassPerUnitLength the cross section mass per unit length AlphaPrincipleAxis the orientation of the cross section elastic axis deter mined at the reference point in radians AlphaPrincipleAxis ElasticCenter the orientation of the cross section elastic center determined at the elastic center in radians AreaX the x coordinate of the area centroid AreaY the y coordinate of the area centroid AreaTotal total cross section area Ixx the mass moment of inertia with respect to the x axis Iyy the mass moment of inertia with respect to the y axis Ixy the mass product of inertia Axx the area moment of inertia with respect to the x axis Ayy the area moment of inertia with respect to the y axis Axy the area product of inertia 5 2 LIST OF FUNCTIONS AND OUTPUT 69 5 2 5 BECAS_RecoverStrains The function BECAS_RecoverStresses is used to determine the three dimensional strains components at the center of each element and at each Gauss point in the cross section finite element mesh The function is meant to be used after the one dimensional beam finite element solution has been determined The resulting strain values are given in the global coordinate systems The function is called as ElementStrain_GlobalCS NodalStrain_GlobalCS BECAS_RecoverStrains theta0 dX dY X Y n1_2d el_2d emat matprops
76. scribed in this section Two different approaches have been implemented The first approach is based on a two dimensional finite element mesh and is mostly attractive for readers which work on their own finite element code The second approach is based on a three dimensional mesh of the cross section using eight node solid finite elements The approach is described in detail in Ghiringhelli and Mantegazza 9 and only the most important results are presented here The approach is implemented in BECAS for illustrative purposes only Although these are standard procedures the transformation of the material constitutive matrix is a topic where in the authors opinion there is often some ambiguity concerning specific definitions for specific implementations Hence the final section of this chapter has been dedicated to the presentation of the mate rial constitutive matrix for orthotropic materials and the corresponding necessary transformations 3 1 Two dimensional finite element analysis 3 1 1 Q4 and Q8 elements The first step in the evaluation of the cross section properties is the generation of a two dimensional finite element mesh of the cross section An example of a discretized profile section using Q4 elements is presented in Figure 3 1 The material properties fiber plane orientation and fiber directions are defined at each element of the finite element mesh Thus a layer of a certain material is defined using a layer of elements Hav
77. ss section stiffness properties while accounting for all the geometri cal and material induced couplings These properties can be consequently utilized in the development of beam models to accurately predict the response of wind turbine blades with complex geometries and made of advanced materials BECAS is based on the theory originally presented by Giavotto et al 6 for the analysis of inhomogeneous anisotropic beams The theory leads to the definition of two types of solutions of which and in accordance to Saint Venant s principle the non decaying solutions are the basis for the evaluation of the cross section stiffness properties A slight modification to the theory was introduced later by Borri and Merlini 7 where the concept of intrinsic warping is introduced in the derivation of the cross section stiffness matrix Despite the modifications no difference in the results was reported The theory was subsequently extended by Borri et al 8 to account for large displacements curvature and twist Ghiringhelli and Mantegazza 2 CHAPTER 1 INTRODUCTION in 9 presented an implementation of the theory for commercial finite element codes Finally Ghiringhelli in 10 11 and Ghiringhelli et al 12 presented a formulation incorporating thermoelastic and piezo electric effects respectively The validation results presented throughout each of the previously mentioned publications highlight the robustness of the method in the analysis of the stiff
78. trix K The function is called using Ks dX dY X Y BECAS Constitutive Ks nl 2d el 2d emat matprops The output is e Ks 6 x 6 array storing the cross section stiffness matrix Kg e dX n4 x 6 array ox e dY 6 x 6 array oY e X nq x 6 array X e Y 6 x 6 array Y where ng Ny x 3 is the number of degrees of freedom in cross section finite element equations 5 2 3 BECAS Constitutive Ms The function BECAS Constitutive Ms is used for the evaluation of the cross section mass matrix M It is called as Ms BECAS Constitutive Ms nl 2d el 2d emat matprops The output is e Ms 6 x 6 array storing the cross section mass matrix Ms 68 CHAPTER 5 USER S MANUAL 5 2 4 BECAS_CrossSectionProps The function BECAS_CrossSectionProps is used to determine a series of relevant cross section properties It presumes that the cross section stiffness K has been previously determined The function is called using ShearX ShearY ElasticX ElasticY MassX MassY MassPerUnitLength AlphaPrincipleAxis AlphaPrincipleAxis_ElasticCenter AreaX AreaY AreaTotal Ixx Iyy Ixy Axx Ayy Axyl BECAS CrossSectionProps Ks nl 2d el 2d emat matprops The output is ShearX the x position of the cross section shear center ShearY the ys position of the cross section shear center ElasticX the x position of the cross section elastic center ElasticY the y position of the cross section elastic center Mas
79. ts of the resulting strains and displacement along the beam axis should also be small All the assumptions mentioned before are not imposed along the cross section coordinates in the cross section plane Finally the theory is based on the assumptions of small displacements and rotations 6 CHAPTER 2 THEORY MANUAL a E Z Figure 2 1 Cross section coordinate system 2 2 Equilibrium equations The derivation of the equilibrium equations for the beam cross section are presented in this section 2 2 1 Basic definitions The reference coordinate system for a generic cross section with area A is presented in Figure 2 1 The displacement of a point in the section s sz sy s 7 is defined with respect to the cross section coordinate system x y z The strain and stress and c are given as T 2 zz Eyy 2 xy 2 xz 2 yz Ezz ot Cra Oyy Oxy Oxz Oyz 027 The stress and strain relate through Hooke s law o Qe 2 1 where Q is the typical material constitutive matrix It is assumed that the material is linear elastic otherwise there are no restrictions regarding the level of anisotropy see Section 3 2 The ordering of the entries in e and is such that the tractions or the components of stress acting on the cross section face can be easily isolated in p Oxz Oyz Ozz The tractions p acting upon the cross section face are statically equivalent to a force T and moment M cf Figure 2 2 T paa A
80. undant Each of the redundancies corresponds to a description of each of the rigid body motions and translations by the warping displacements u It is therefore necessary to ensure that the warping displacement u is uncoupled from the rigid displacements r The following set of constraints are therefore incorporated into the solution of the sets 2 24 Nn Nn Nn gt Uam 0 X Uyn 0 X 30298 0 n 1 n 1 n 1 Nn Nn Nn 1 Zn Un YnUzn 0 gt nUxn TnUzn 0 gt Yn lan TnUyn 0 n l n l n l where n is the number of nodes in the cross section finite element mesh and Zn Yn Zn and an Uyn Uz n are the position and displacement of node n respec tively The constraints are imposed on both the displacements u and corresponding derivatives Qu Oz and can be written in matrix form as D 0 u 0 I5 I5 T 0 p i o vie n g ss Ma where I3 is the 3 x 3 identity matrix and np is obtained from replacing the nodal coordinates n Yn Zn of node n in 2 2 1 2 4 On the properties of the solutions Some results used earlier in this manuscript to justify some of the steps in the deriva tion of the cross section stiffness matrix of beams are described in this section In particular 2 28 and 2 32 are an important result in establishing the equilibrium equations for the central solutions in Section 2 3 2 The reader may skip this section altogether if only a general overview of the method is necessar
81. y The section is divided in two parts In the first part it is assumed that only the rigid body translations and rotations contribute to the displacement vector of a point in the cross section i e the cross section deformation is not accounted for The equilibrium equations are derived and it is possible to show that the force and moment resultants O vary linearly with respect to z or along the beam length In the second part the displacement of a point in the cross section is obtained through a finite element discretization of the warping displacements The equilib rium equations are derived once again The result is a second order homogeneous linear differential equation The different types of solutions are identified according to Saint Venant s principle The solutions for which the eigenvalues are different from zero will decay as z either increases or decreases These are self balanced 0 0 solutions corresponding to the modes at the extremities of the beam where the loads are applied the extremity solutions On the other hand the solutions for which the eigenvalues are zero will not present any decay and most importantly it is possible to show that these will be polynomials in z These correspond to non zero stress resultants and account for the displacements at a section of the beam 2 4 ON THE PROPERTIES OF THE SOLUTIONS 17 sufficiently far from the extremities so that the influence of the external forces in the stress fiel
82. ysis Software BECAS is a group of Matlab functions used for the analysis of the stiffness and mass properties of beam cross sections BE CAS was originally developed under the EFP 2007 Project 33033 0075 Anisotropic beam model for analysis and design of passive controlled wind turbine blades BE CAS was later updated extended and completely rewritten throughout part of the author s Ph D project Optimal Design of Laminated Composite Beams Ph D Thesis Technical University of Denmark BECAS code and user s guide is mostly developed and maintained by Jos Pedro Blasques Riso DTU National Laboratory for Sustainable Energy Technical University of Denmark Nonetheless the author is indebted to the following people which at one point or another have given or currently give invaluable support throughout the development of BECAS e Boyan Lazarov Department of Mechanical Engineering Technical University of Denmark Participated very actively in the development of the original BECAS v1 0 Among much other invaluable work Boyan was the first to suggest the constraint equations which are used in the solution of the cross section equilibrium equations Robert Bitsche Riso DTU National Laboratory for Sustainable Energy Technical University of Denmark An active member of the current BECAS development group Robert is the main bug finder and an invaluable source of good ideas and suggestions Robert is also responsible for the interfac

User's Manual for BECAS - DTU Orbit

Contents

Download Pdf Manuals

Related Search

Related Contents