FEM Subsystem Replacement Techniques for Strength Problems in Variable Geometry Trusses

02 / 19 / 2007 Macareno 1 FEM SUBSYSTEM REPLACEMENT TECHNIQUES FOR STRENGTH PROBLEMS IN VARIABLE GEOMETRY TRUSSES. Luis M. Macareno*, Josu Agirrebeitia, Carlos Angulo, Rafael Avilés Department of Mechanical Engineering, University of the Basque Country Alameda Urquijo s/n, 48013 Bilbao, Spain * Corresponding author: Luis M. Macareno Department of Mechanical Engineering, University of the Basque Country UPV/EHU Escuela Técnica Superior de Ingeniería, ETSI Alameda de Urquijo s/n, 48013 Bilbao SPAIN Tel: +34 946017399, Fax: +34 946014215 E-mail: [email protected] Number of words: 5267 Number of figures: 15 Tables: 2 Abstract This work presents the application of a procedure for replacing FEM subsystems with a high number of dofs (possibility of including mobile internal parts) using Equivalent Parametric Macroelements (EPMs) with a much reduced number of elements to This is the accepted manuscript of the article that appeared in final form in Finite Elements in Analysis and Design 44(6/7) : 346-357 (2008),which has been published in final form at https://doi.org/10.1016/j.finel.2007.12.003. © 2007 Elsevier under CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/) FEM subsystem replacement techniques for strength problems in variable geometry trusses 02 / 19 / 2007 Macareno 2 decrease significantly the analysis time with an acceptable error. This procedure is applied to the replacement of VGT mobile joints. The equivalence criterion proposed is based on elastic strain energies absorbed by both bodies. Said replacement involves resolution of a redundant non-linear equation system, iteratively focused via linearization and subsequent resolution via the Least Squares Method. The search for initial approximation is supported by Genetic Algorithm techniques. Key words: Variable Geometry Truss (VGT), Finite Element, Equivalent Parametric Macroelement (EPM), Energy Method, Optimization. 1. Introduction. Variable geometry structures are those capable of modifying their geometry to adapt to different loads and working conditions. This is possible because some of the elements comprising them can vary their length. These elements are called actuators. Another characteristic making them interesting is their high stiffness to weight ratio, which has contributed to the application of variable geometry structures in the spatial research field. The study of these structures dates back to the 1980s [1]. A specific type within this group is based on spatial truss type structures known as Variable Geometry Truss (VGT). Its most common application is as manipulators [2, 3 and 4]. These trusses are formed via repetition of the main module whose topology can be highly varied [5, 6].The most widely known VGT is the Double Octahedral, made up of two main octahedral modules [7, 8, 9]. FEM subsystem replacement techniques for strength problems in variable geometry trusses 02 / 19 / 2007 Macareno 3 The MBAD 1 (Multi Body Analysis and Design) workgroup of the University of the Basque Country has developed a five-module VGT prototype, where the geometry of the main module is also established upon the octahedral shape. The real structure is shown in Fig. 1. Each module is a parallel kinematic mechanism in itself with actuators set on the horizontal planes. These planes are joined among them using fixed length bars called longerons. The joint between these bars and the actuators is done by special joints. These joints have also been developed by the workgroup and are patented. In Fig. 2 there is an exploded view of the joint where it can be appreciated the different elements comprising the same. Fig. 1. Five-Module VGT Prototype Fig. 2. Exploded view of the joint. 1 http://www.ehu.es/mbad Special piece for bearings Joint body Double ring Simple ring FEM subsystem replacement techniques for strength problems in variable geometry trusses 02 / 19 / 2007 Macareno 4 To perform the analysis on the structure a finite element model was created in MSC/Nastran, preceded by another detailed cinematic model created using MSC/Adams. Fig. 3 shows both models indicating the different elements comprising the five-module VGT. Fig. 3. VGT kinematic and FEM model In these structures it is necessary to perform the analysis on various configurations. A program has been developed in PCL (Patran Command Language) which automatically creates a FEM model in each position to enable diverse studies. Once the FEM model of the variable geometry structure has been created, its behavior is studied under different load cases and in different positions. As this is a detailed model, both creation and analysis time are excessive, so the need to reduce the model arises. There are different static model reduction techniques, where the most widely known are reduction via static Longeron Joint Batten Actuator Base FEM subsystem replacement techniques for strength problems in variable geometry trusses 02 / 19 / 2007 Macareno 5 condensation [10, 11], iterative techniques [12, 13], macroelements [14, 15] and substructuring techniques [16]. This paper presents the application of a static reduction technique for FEM models via replacement of certain submodels by equivalent parametric macroelements (EPMs). These macroelements will have a much lower number of elements than the submodels they replace, and consequently computation cost is drastically reduced. The technique is based on elastic strain energy equivalence between the submodel and its corresponding macroelement. This criterion has led to good results in other applications [17, 18]. The model to be analyzed is the variable geometry structure shown in Fig. 1, and the submodels to be replaced are the intermediate joints of the structure. These submodels consist of 3D tetrahedral and hexahedral type elements, likewise rigid and contact elements. Altogether there are approximately 2000 elements per submodel. Fig. 4 shows the FEM model of one of the joints. Fig. 4. Joint FEM Model FEM subsystem replacement techniques for strength problems in variable geometry trusses 02 / 19 / 2007 Macareno 6 The following section defines the mathematical development of the macroelement in detail. Once its function has been understood, the methodology of the model reduction technique is explained; defining all the variables likewise the equations arising from energy equivalence. Below are the results obtained via application of this technique for the joint of the five-module VGT. Finally, a brief commentary on these results and the conclusions are shown. 2. Mathematical definition of the joint macroelement. The macroelement as understood in this paper is a series of elements, which join together forming a new body. If this set of elements is intended to replace another model, it might satisfy certain restraints: firstly it must be cinematically equivalent to the replaced model and secondly equivalent from a structural viewpoint. The choice of the type of elements is free as far as the connections with the rest of the model and the internal mobility are preserved. In this case, the defined macroelement to replace the finite element submodel of the joint is formed by eight 3D beam elements with circular section and nine nodes, as shown in Fig. 5. Each of these nodes has six dofs, three translations (ux, uy, uz) and three rotations (θx, θy, θz). The beams are arranged so nodes 1, 3, 6 and 7 correspond to the central points of the joint spherical bearings and nodes 8 and 9 to the link points of the rings with actuator bars. It can be appreciated node 9 really correspond to the point midway between the two link points of the double ring. Nodes 2, 4 and 5 are internal and have not physical correspondence with any of the joint. Node 4 is located on the central point belonging to the rotation axis of the joint rings. In Fig. 6 it can be appreciated the correspondence of the macroelement nodes with the joint submodel. FEM subsystem replacement techniques for strength problems in variable geometry trusses 02 / 19 / 2007 Macareno 7 Fig. 5. Macroelement Fig. 6. Submodel and macroelement correspondence. 2.1 Stiffness matrix. The macroelement stiffness matrix is obtained via the sum of expanded matrices of each element, once expressed into the global system. The axes of this global system are XYZ shown in Fig. 5. As the macroelement consists of 9 nodes each with 6 dofs, the final matrix is 54 x 54. The element matrices are 12 x 12, since 3D beam type elements were used. As said before, the macroelement must be cinematically equivalent to the submodel. The rings round the joint body can rotate in relation to the central axis. These rotation FEM subsystem replacement techniques for strength problems in variable geometry trusses 02 / 19 / 2007 Macareno 8 angles shown in Fig. 5 have been defined as θ1 and θ2. The rotation condition must be included in the macroelement. This has been done releasing the rotation dof regarding the central axis on the nodes corresponding to elements 7 and 8. These elements would correspond to the rings in the submodel. The fact of releasing dofs on a node is equal to including the null force transmission condition in the same direction. The node undergoing this situation is number 4, where four elements concur numbered 3, 4, 7 and 8. Elements 7 and 8, as mentioned before, would correspond to the rings and elements 3 and 4 would represent the central axis of the real joint. Therefore, contribution of the moment in direction Z on node 4 of elements 7 and 8 must be cancelled. Fig. 7 shows released dofs, likewise the local systems of elements 7 and 8. Anyway on this node, effort continues being transmitted in that direction, although only between elements 4 and 5. Introducing this condition means the stiffness matrix of elements 7 and 8 suffer certain modifications directly influencing the macroelement stiffness matrix. Below the mathematical development is shown in detail. Fig. 7. Released dofs FEM subsystem replacement techniques for strength problems in variable geometry trusses 02 / 19 / 2007 Macareno 9 Moment according to the released dof is not transmitted, so it is null. In element 7 it corresponds to Mz on the first node of the element. Thus, using its equilibrium equation (1) and equaling to zero it can be obtained the released dof which now becomes dependent on the rest. 22 2 11 2 1 2646 0z z y z z z y z zL IE u L IE L IE u L IE M           (1) 2 )( 2 32 121 z yyz uu L      (2) Eq. (2) is replaced in all the matrix rows where term θz1 appears and the new matrix is obtained. This matrix has zeros in all the row and column corresponding to the released dof. This development is analogous for element 8. All relevant matrices are shown on the Appendix. Once the previous steps have been executed to achieve cinematic equivalence, the macroelement global stiffness matrix is obtained via coupling of all the matrices of each element in global coordinates. At this point, the macroelement has been completely defined, which depends on the physical properties (E, υ, G) and dimensions of the elements comprising the same, such as the section or radii. Any of these variables may be selected as a parameter. 2.2 Macroelement stiffness matrix condensation. FEM subsystem replacement techniques for strength problems in variable geometry trusses 02 / 19 / 2007 Macareno 16 Applying Eq. (20) for each displacement and each parameter, the sensitivity matrix of dimensions n x p can be obtained.                                                p nnn p p A V A V A V A V A V A V A V A V A V A V ... ............ ... ... 21 2 2 2 1 2 1 2 1 1 1 (21) Taking into account all the above, the following linearized equation system is formulated for k iteration, where approximate elastic strain energy is equaled to the submodel energy for each displacement case i. ni AA A V VU p j k j k j k j k i k ii ...,1 )( 1 1         (22) This system must be solved for each iteration k via the Least Squares Method, for which the R error function is defined as follows:   2 1 1 1 1)()(                 n i p j k j k j k j k i k ii kAA A V VUAR (23) FEM subsystem replacement techniques for strength problems in variable geometry trusses 02 / 19 / 2007 Macareno 17 To obtain values which minimize the error function the derivatives must be equaled to zero for each parameter. Thus, the following determined compatible equation system is reached of p equations and p variables.   pj A AR k j k ...,1 0 )( 1 1       (24) The solution of this system is the {A}k+1 parameter vector. Set out below are the stop criteria taken into account in the iterative process. 3.2 Stop criteria. To assess the degree of equivalence of the macroelement a value comparing the vectors {U} and {V} must be defined. The following formula is highly appropriate for this case:   n VU n i k ii     1 2 1  (25) With Eq. (25) it can be evaluated the error committed in each iteration. In the general case, as this is a Least Squares method, it will not get a zero error. Therefore, the most interesting stop criterion is that which values stabilization of the error committed. When the error percentage increase in two consecutive iterations is less than a certain q value, the iterative process may be deemed finished. FEM subsystem replacement techniques for strength problems in variable geometry trusses 02 / 19 / 2007 Macareno 18 q k kk      1 (26) In Eq. (26) the absolute value is used because the error is expected to be less in each iteration if the method converges. As an additional criterion it can also be used the one based on the optimization parameters variation norm (27).       ' 2 1 1 1qAAAA p j k j k j kk      (27) However, it must be borne in mind, that it is highly likely there will be a certain number of parameters whose variation has no decisive impact on the objective function value. Therefore, prior to assessing the application of Eq. (27), a detailed study must be performed on parameter sensitivity, as explained in the parameter control section. 3.3 Initial approximation. The selection of initial approximation tasks requires special care, since it has been verified the method used is highly sensitive to said choice, both in process duration and stability. Here the use of genetic algorithm techniques has been chosen. A simple genetic code in real numbers with 5 digits in decimal base indicating the value of optimization variables has been used for this. The variables chosen for the genetic algorithm will be detailed in the results section. 3.4 Parameter control. FEM subsystem replacement techniques for strength problems in variable geometry trusses 02 / 19 / 2007 Macareno 19 Among the optimization methods, this one belong to those of steepest descend or gradient. The difference between the vector obtained as a solution and the previous in each iteration provides a direction vector of maximum gradient. This means that varying the parameters in this direction, a maximum variation of the objective function is achieved. This vector may be distorted in some cases. This is because the variation of some parameters has almost a null impact on the objective function; whereby its value may increase disproportionately without appreciable improvement in the solution. To prevent situations of this kind an increase control has been applied to the parameters. This control consists of applying a maximum value allowed. Once the maximum gradient vector has been obtained, variation of each parameter is calculated. If all values are lower than the percentage permitted, the process continues unchanged. Should said value be exceeded, the vector is scaled taking it to the allowed limit. This action prevents any parameter from growing exaggeratedly yet continues to maintain the maximum variation direction. The maximum variation value allowed is a parameter which may be either variable or fixed. Should it be variable, the most logical would be for it to decrease as the solution approaches optimum, to prevent leaping to another solution. In this study a fixed f value was chosen. Apart from delimiting the parameter maximum percentage variation, sometimes a compromise decision must be taken in relation to the value of some of them. Applicable to the extreme cases commented above of parameters varying disproportionately, without hardly any improvement in the solution. This behavior may be quantified using the information provided by the sensitivity matrix. Each matrix term indicates the variation in elastic strain energy Vi for a unit value increase in the corresponding parameter. When this value is close to zero, it means the energy hardly varies. This FEM subsystem replacement techniques for strength problems in variable geometry trusses 02 / 19 / 2007 Macareno 20 value may be very low for one displacement case and high for another. Therefore, a value which considers all cases must be taken as a reference. This value has been defined as Sj, for each parameter. It must be calculated in each iteration since the sensitivity matrix also varies. n A V S n ij i j            1 2 (28) The value obtained in Eq. (28) includes the objective function sensitivity for all displacement cases in relation to parameter Aj. As previously indicated, this sensitivity enables function invariance evaluation in relation to the parameter in question. FEM subsystem replacement techniques for strength problems in variable geometry trusses 02 / 19 / 2007 Macareno 21 3.5 Flow diagram of the optimization process and term glossary. Fig. 8. Flow diagram of the optimization process FEM subsystem replacement techniques for strength problems in variable geometry trusses 02 / 19 / 2007 Macareno 22 Term glossary. n total number of displacement cases i subindex of current displacement case h number of displacement vector dofs {δ}i displacement vector of case i [δ] vector displacement matrix Ui elastic strain energy of the model for case i {U} vector of elastic strain energies of the model Vi macroelement elastic strain energy for case i {V} vector of macroelement elastic strain energies p total number of parameters j current parameter subindex Aj parameter j {A} parameter vector m increasing percentage of parameters in numerical derivation sensitivity matrix k iteration number R error or objective function ε absolute error q relative error f maximum variation in parameters control Sj parameter j sensitivity         A V FEM subsystem replacement techniques for strength problems in variable geometry trusses 02 / 19 / 2007 Macareno 23 4. Results. This section sets out the optimization process results applied to the macroelement described in section 2. The material used is aluminum with the following characteristics: Young module: E = 70 GPa Shear module: G = 27 GPa Poisson Coefficient: υ = 0.33 Method application conditions: Number of macroelement parameters: p = 8 Number of master dofs: h = 18 Number of displacement cases: n = 60 Increased percentage in the sensitivity matrix calculation: m = 0.01 Maximum Number of iterations: k_max = 15 Value of admissible q: q = 0.5 % Maximum variation in parameter control: f = 10 % Step one is to create the 60 displacement cases, for which vectors {δ} are created randomly. The magnitudes of these displacements were chosen so they were higher than expected for the usual working conditions of the real joint. Thus, it is guaranteed correct macroelement behavior under normal conditions. The random function (29) is centered at zero with values between -10-5 and 10-5 meters. FEM subsystem replacement techniques for strength problems in variable geometry trusses 02 / 19 / 2007 Macareno 24 )12(10 5 Random  (29) Once created the matrix [δ] with the 60 vectors {δ} each of 18 components, vector {U} is calculated. This vector can be obtained in different ways, the most common being calculation of the model elastic strain energy for each displacement case as per the following formula:     i T i iFU   2 1 (30) {F} being the vector comprising components of reaction forces associated to master dofs where displacements are applied. Repeating Eq. (30) for all displacement cases, the vector {U} is obtained. As mentioned previously this vector is only calculated once. The search for initial approximation {A}0 was carried out executing the free Genetic Algorithm code “Pikaia” [19]. The values adopted for the algorithm values are as follows: Number of individuals in the population: 150 Number of generations in the evolution: 400 Encoding: decimal Number of digits to encode e genotype: 5 Crossover probability: 0.85 Mutation mode: adjustable rate based on fitness Initial mutation rate: 0.005 FEM subsystem replacement techniques for strength problems in variable geometry trusses 02 / 19 / 2007 Macareno 25 Minimum mutation rate: 0.0005 Maximum mutation rate: 0.25 Reproduction plan: full generational replacement Elitism: active In this work it was decided to limit the parameter value bearing in mind the physical features of the submodel. On the other hand, the search space was notably reduced taking into account the symmetry extant in the submodel, which implies a certain relationship “a priori” among some parameters. After 400 generations the values obtained for initial approximation (expressed in millimeters) were the following:     T A810.2971.4822.10917.10000.25000.25999.8989.8 0 These are the values initially defining the macroelement; i.e. the values of the 8 radii of the circular section beams forming the same. Now the elastic strain energy Vi can be calculated for each displacement case. By calculating these values for the 60 cases, vector {V} is formed. Now, it can be assessed the degree of equivalence between the macroelement and the submodel, using the Eq. (25), which obtained a value of ε = 0.0292. Fig. 9 represents the quadratic difference values of the U and V energies per displacement case i. This figure confirms the V values are considerably similar to those of U. FEM subsystem replacement techniques for strength problems in variable geometry trusses 02 / 19 / 2007 Macareno 32 References [1] Miura K., Furuya H. and Gokhale D. Variable Geometry truss and its application to deployable truss and space crane arm. Acta Astronautica 1985;12;7:599-607. [2] Gun-shing Chen and Wada B. K. Adaptive Truss Manipulator Space Crane Concept. Journal of Spacecraft and Rocket 1993;30;1:11-5. [3] Stoughton R. S., Tucker J.C., Horner C.G. A Variable Geometry Truss Manipulator For Positioning Large Payloads. In: American Nuclear Society Topical on Robotics and Remote Handling Conference, Monterey, California, USA. 1995:1-8. [4] Hughes P. C. ; Sincarsin W. G. ; Carroll K. A. Trussarm-A variable Geometry Truss Manipulator. Journal of intelligent material systems and structures. 1991;2;2:148-60. [5] Sircasin W. G. and Hughes P. C. Trussarm Candidate Geometries, Dynacon Report 28-611/0401. 1987. [6] Arun V., Reinholtz C. F. and Watson L. T. Enumeration And Analysis Of Variable Geometry Truss Manipulators. Proceedings of ASME Mechanisms Conference 1990:93-8. [7] Shengyang Huang, Natori M.C. and Miura K. Motion Control of Free-Floating Variable Geometry Truss. Part 1: Kinematics. Journal of Guidance, Control and Dynamics 1996;4;19:756-63. FEM subsystem replacement techniques for strength problems in variable geometry trusses 02 / 19 / 2007 Macareno 33 [8] Williams II R.L. and Hexter IV E.R. Maximizing Kinematic Motion For A 3-Dof VGT Module. Journal of Mechanical Design 1998;120;2:333-36. [9] Chen Wu-Jun, Luo Yao-Zhi, Fu Gong-Yi, Gong Jing-Hai and Dong Shi-Lin. A Study on Space Masts Based on Octahedral Truss Family. International Journal of Space Structures 2001;1;16:75-82. [10] Guyan J. Reduction of Stiffness and Mass Matrices. AIAA 1965;3:380. [11] Wilson Edward L. The static condensation algorithm. Journal for numerical methods in engineering 1974;8:198-203. [12] Friswell M. I. Model Reduction Using Dynamic and Iterated IRS Techniques. Journal of Sound and Vibration 1995;186:311-23. [13] Friswell M. I. The Convergence Of The Iterated IRS Method. Journal of Sound and Vibration 1998;211:123-32. [14] Parsons Thomas J., GangaRao Hota V. S. and Peterson William C. Macro-element Analysis. Computer and Structures 1985;20:877-83. [15] Kwan A. S. K. and Pellegrino S. Matrix Formulation of Macro-elements for Deployable Structures. Computer and Structures 1994;50:237-54. FEM subsystem replacement techniques for strength problems in variable geometry trusses 02 / 19 / 2007 Macareno 34 [16] Carpenter Donald L. and Snyder Virgil W. Implementation of Substructuring in Finite Elements. Proceedings of the Third Engineering Mechanics Division Specialty Conference, Austin, Texas 1979:377-79. [17] Fernndez-Bustos, I. , Aguirrebeitia, J., Avilés, R. Angulo, C. Kinematical Synthesis of 1-dof Mechanisms Using Finite Elements and Genetic Algorithms, Finite Elements in Analysis and Design 2005;41;15:1441-63. [18] Avilés, R , Ajuria, G., Amezua, E., Gmez-Garay, V. Finite Element Approach to the Position Problems in Open-loop Variable Geometry Trusses. Finite Elements in Analysis and Design 2000;34;3-4:233-55. [19] Charbonneau, P. An Introduction to Genetic Algorithms for Numerical Optimization. NCAR Technical Note 450+IA (Boulder: National Center for Atmospheric Research);2002. FEM subsystem replacement techniques for strength problems in variable geometry trusses 02 / 19 / 2007 Macareno 35 Appendix. Element stiffness matrices. Fig. 14. Stiffness matrix of elements 1,2,3,4,5 and 6 in local coordinates Fig 15. Stiffness matrix of elements 7 and 8 in local coordinates The row and column corresponding to the released dof θz1 all contain zeros. In this component, these elements do not contribute stiffness to the macroelement. FEM subsystem replacement techniques for strength problems in variable geometry trusses 02 / 19 / 2007 Macareno 36 Figure Legends: Fig. 1. Five-Module VGT Prototype. Fig. 2. Exploded view of the joint. Fig. 3. VGT kinematic and FEM model. Fig. 4. Joint FEM Model. Fig. 5. Macroelement. Fig. 6. Real model and macroelement correspondence. Fig. 7. Released dofs. Fig. 8. Flow diagram of the optimization process. Fig. 9. Difference between U and V energies for the 60 displacement cases. Fig. 10. Error evolution. Fig. 11. Parameter evolution. Fig. 12. Variation percentages of parameters. Fig. 13. Error comparison for several “displacement case groups”. Fig. 14. Stiffness matrix of elements 1,2,3,4,5 and 6 in local coordinates. Fig. 15. Stiffness matrix of elements 7 and 8 in local coordinates. Tables: Table 1. Values of parameters and their relative variation. Table 2. Parameter sensitivity values.