SEMI-ANALYTICAL METHOD FOR THE ANALYSIS AND DESIGN OF CYLINDERS WITH CONTROLLABLE ELASTIC POST-BUCKLING RESPONSE By Ali Imani Azad A DISSERTATION Submitted to Michigan State University in partial fulfillment of the requirements for the degree of Civil Engineering—Doctor of Philosophy 2022 ABSTRACT SEMI-ANALYTICAL METHOD FOR THE ANALYSIS AND DESIGN OF CYLINDERS WITH CONTROLLABLE ELASTIC POST-BUCKLING RESPONSE By Ali Imani Azad Research over the past ten years has generated an increased interest in studying elastic structural instabilities as a useful response for smart applications rather than a failure. Buckling under axial compression is a type of structural instability that can be used for rapid geometric switching and energy harvesting if the resulting deformations are properly controlled, which usually implies the need of external constraints. Short thin-walled cylinders are able to experience multiple buckling events without the need of additional constraints. However, predicting their post-buckling response is challenging due to their high sensitivity to initial imperfections. Cylinders with non-uniform stiffness distribution (NSD), which feature wall regions with varying elastic modulus or thickness, allow controlling the number of elastic buckling events, their location, and the sequence at which they occur. However, designing NSD cylinders for smart applications requires predicting, for each buckling event, the corresponding compressive force and end shortening, the resulting drop in force and strain energy, and the post-buckling stiffness. This dissertation presents a semi-analytical model developed to predict the elastic post- buckling response of NSD cylinders under compression. The model follows three general steps: 1) separate the NSD cylinder into parallel segments, 2) simplify and predict the response of each segment, and 3) integrate the response of individual segments. The first step in predicting the elastic post-buckling response of a cylindrical segment was to simplify its geometry into a cylindrical panel with uniform thickness. Linear springs are connected to the top and bottom of the uniform cylinder to match the stiffness of the simplified segment to the actual one. Based on classical shell theory, the elastic post-buckling response of a cylindrical panel is solved as a boundary value differential equation using the pseudo-arclength method. The analytical model was validated by comparing the predicted post-buckling responses of four cylinders with those from experiments and finite element analyses. The response of cylindrical panels of various dimensions is needed to design NSD cylinders for targeted post-buckling behavior. Thus, the differential equation for a cylindrical panel under axial compression was solved independently of the cylinder’s radius and elastic modulus. These allowed the development of design maps for several parameters, including axial strain and stresses corresponding to the first buckling event. Predictive equations generated by using genetic programming were then used to relate each design parameter to the geometry of the panels. Three NSD cylinders were designed using the developed design maps to validate the proposed approach. One was designed to undergo multiple buckling events under compression at pre- defined end shortenings. A second cylinder was designed to feature a post-buckling force- deformation response that plateaus at a constant force level. The third cylinder was designed to experience the same force drop at each buckling event and for equal axial end shortenings after the first event. Finite element analyses verified that the proposed design procedure yielded NSD cylinders with a post-buckling response close to the desired one, and the ultimate design goal can be achieved by slight modifications to the cylinder geometry. This study advances the knowledge on the elastic buckling and post-buckling response of slender cylindrical shells under axial compression and provides an approach to analyze and design them for a desired far post-buckling response. This work expands the harnessing of elastic instabilities to the area of thin shells buckling under compression, which has received less attention in comparison to other forms of achieving controllable structural instabilities. Copyright by ALI IMANI AZAD 2022 This dissertation is dedicated to My wonderful family for their endless love, advice, and encouragement. My mother Jamileh, my lovely wife Samira, my beautiful daughter Saba, and my late Father Moharram. v ACKNOWLEDGMENTS My Ph.D. program was a long journey with many ups and downs. If it was not for the patience, vision, kindness, and knowledge of my advisor, Professor Rigoberto Burgueño I was not have been able to conclude this journey. When I look back to spring 2017, that I walked into his office and requested to join his research group, there is no doubt about the progress I had under his supervision. He walked me step by step to the goal of my education and did it without losing my independent thinking and confidence. He encouraged me to think critically, dared me to show strength, and demanded to push myself every day. He was the best mentor that I could ever imagine. Professor Burgueño was not only a great researcher, but also the best teacher I have ever seen, and believe me, I saw so many teachers. Professor Burgueño, I cannot express all my appreciation and gratitude for your help and mentoring in words; you supported me in the hardest days of my life. Thank you, Sir. I am also thankful to my committee members, Dr. Thomas Pence, Dr. Weiyi Lu, and Dr. Mahmoodul Haq. I had the opportunity to sit in three classes by Dr. Pence and learned many things from him. He is a passionate teacher, and he is so excited about engineering and science that he turns boring topics into exciting problems. Dr. Weiyi Lu was my mentor when I taught a structural analysis course. I remember when he did not like something about my teaching and wanted to criticize me, he was kind and supportive but direct to the point. I am thankful for the lessons he taught me as a mentor. I am happy that I know Dr. Mahmoodul Haq. He is my teacher, second advisor, mentor, and a caring, kind brother. Whenever I felt disappointed in my life, research, job, or any other things, I just ran to his office and had a short talk with him. His calm attitude and great advice changed the perspective on my life. He has a special place in my heart. vi A special thanks to my lovely, patient, and kind wife, Samira. I could not keep working if she were not supporting me. She was patient for many years and did everything to help me with my Ph.D. Samira took care of all my responsibilities in life and let me focus on my work. I was asked many times how I could keep working for long hours, and the answer was "Samira." She sacrificed everything for me, bared responsibility for our lovely daughter Saba at her young age and let me focus on my work. She is the best! I am grateful to all my friends and colleagues at MSU, especially Shabnam Rajaei and Leila Khalili, for their encouragement when I needed it the most. My gratitude to my colleague Suihan Liu for her collaboration in two research projects. I am also grateful to Nan Hu and Jun Guo for their prior work on the NSD cylinder project. In addition, without financial support from U.S. National Science Foundation and assistantships from Civil Engineering Department at MSU, I could not complete my Ph.D. program. Many thanks to Professor Neeraj Buch for the financial support during the last years of my Ph.D. I also should express my gratitude to professor Buch for the great opportunity of teaching in the Civil Engineering department. Last and not least, I could never accomplish anything without the encouragement and endless support of my mother, Jamile Zamani, and my late father, Moharram Imani Azad. He was my role model, best friend, strongest critic, and the best man I have ever known in my life. May he rest in peace. vii TABLE OF CONTENTS LIST OF TABLES .......................................................................................................................... x LIST OF FIGURES ....................................................................................................................... xi Chapter 1. Introduction ................................................................................................................... 1 1.1. Motivation ...................................................................................................................... 1 1.2. Hypothesis ...................................................................................................................... 4 1.3. Research significance ..................................................................................................... 4 1.4. Background .................................................................................................................... 5 1.4.1. Buckling in cylinders ............................................................................................... 5 1.4.2. Hoop and string stiffeners ....................................................................................... 6 1.4.3. Buckling capacity knockdown factors..................................................................... 7 1.4.4. Imperfection insensitive cylinders ........................................................................... 9 1.4.5. Cylindrical shells with non-uniform stiffness distribution .................................... 12 1.4.6. Research gap .......................................................................................................... 15 1.5. Research plan ............................................................................................................... 16 Chapter 2. A framework to predict the post-buckling response of NSD cylinders ...................... 19 2.1. Dividing NSD cylinder into parallel segments ............................................................ 19 2.1.1. Boundary conditions .............................................................................................. 20 2.1.2. Divider segments ................................................................................................... 22 2.1.3. Equivalent Panels .................................................................................................. 23 2.2. Composition axial responses of separated segments .................................................... 26 2.3. A semi-analytical model for axially compressed cylindrical panels ............................ 29 2.3.1.1. Analytical Solution for cylindrical panel with all edges clamped ................. 30 2.3.1.1.1. The radial deformation functions ............................................................ 31 2.3.1.2. Model verification for post-buckling response .............................................. 44 Chapter 3. NSD buckling segment discussions ............................................................................ 48 3.1. Thickened area and equivalent panel ........................................................................... 48 3.2. Segment radial deformation after the buckling event .................................................. 56 3.3. Segment energy drop from buckling ............................................................................ 58 3.4. Verification of the assumptions.................................................................................... 60 3.4.1. Verification of modeling framework ..................................................................... 63 Chapter 4. Design of NSD cylinders with favored post-buckling response ................................. 68 4.1.1. Design maps for cylindrical panels with specific radius and thickness ................ 69 4.1.2. NSD cylinder with known end shortenings for each buckling event .................... 72 4.2. General Design maps ................................................................................................... 74 4.2.1. Determining acceptable range for panel parameters ............................................. 76 4.2.2. Design Maps .......................................................................................................... 79 4.2.2.1. Axial Strain corresponding to the buckling ................................................... 81 4.2.2.2. Axial stress corresponding to the buckling .................................................... 82 viii 4.2.2.3. Normalized force drop corresponding to the buckling .................................. 83 4.2.2.4. Normalized energy drop corresponding to the buckling................................ 84 4.2.2.5. Normalized secondary stiffness of the panel after buckling .......................... 85 4.2.2.6. Normalized radial deformation after the buckling ......................................... 86 4.2.2.7. Normalized change in the radial deformation after the buckling .................. 87 4.2.2.8. Normalized von mises stress after the buckling ............................................ 88 4.2.2.9. Normalized slope of the von mises stress after the buckling ......................... 89 4.2.3. Discussion.............................................................................................................. 90 4.2.4. Predictive models .................................................................................................. 94 4.3. Design problems ........................................................................................................... 97 4.3.1. NSD cylinder with plateau post-buckling response .............................................. 98 4.3.2. NSD cylinder with equal force drops in buckling events .................................... 108 Chapter 5. Summary and Conclusion ......................................................................................... 112 5.1. Summary .................................................................................................................... 112 5.2. Conclusions ................................................................................................................ 113 5.3. Potential future research ............................................................................................. 114 APPENDIX ................................................................................................................................. 116 BIBLIOGRAPHY ....................................................................................................................... 142 ix LIST OF TABLES Table 2-1: Dimensions of cylindrical panels for verification of the analytical model ................. 46 Table 3-1: Geometry of NSD cylinder segments (all dimensions are in mm).............................. 50 Table 3-2: Geometry of NSD cylinders tested under uniform axial compression ........................ 64 Table 4-1: Dimensions of cylindrical panels to develop empirical models .................................. 78 Table 4-2: Dimensions of the cylindrical panels used for the test in genetic programming ........ 95 Table 4-3: Relation between design parameters and design variables (l/R, b/R, t/R) .................. 96 x LIST OF FIGURES Figure 1-1: (a) Slender beam under compression buckles and deflection keeps growing until damage occurs, (b) the same slender beam under compression but constraint between rigid wall shows several buckling events, (C) axially compressed cylindrical shells with uniform thickness show several buckling events.......................................................................................................... 2 Figure 1-2: Cylindrical shells with nonuniform stiffness distribution (NSD cylinder) limits the buckling in the cylinder to specific areas ........................................................................................ 3 Figure 1-3: Empirical knockdown factors for the design of cylindrical shells ............................... 7 Figure 1-4: (a) section and 3D shape of Aster cylinder, (b) section and 3D face of the corrugated cylinder by Ning and Pellegrino ................................................................................................... 10 Figure 1-5: The thickness of the wall varies in CTS cylinders ..................................................... 12 Figure 1-6: Seeding imperfections by Hu and Burgueño; (a) sinusoidal seeded imperfection, (b) changing the pattern of shell's wall thickness, and (c) radial constraints [60].............................. 13 Figure 1-7: Buckling areas in cylinders with non-uniform stiffness distribution and load drops in the axial response of the cylinders [60] ........................................................................................ 14 Figure 1-8: Schematic research plan (a) predict the post-buckling response of any given NSD cylinder and (b) design an NSD cylinder for a given post-buckling response ............................. 18 Figure 2-1: Wall-thickness pattern along the surface of a cylinder with Non-uniform Stiffness Distribution (NSD); h is the cylinder height, and t is the wall thickness of each area ................. 19 Figure 2-2: Division of the NSD cylindrical shells into parallel independent panels .................. 20 Figure 2-3: Sections in the longitudinal and circumferential direction of the cylinder to investigate the boundary conditions in buckling segments ............................................................................. 21 Figure 2-4: Sequence of radial deformation in the circumference of 50 mm ............................... 22 Figure 2-5: Replacing thick zones (t=t1i) in the panel with an equivalent length of the shell with identical thickness to mid-section (t=tb) ....................................................................................... 24 Figure 2-6: Framework to obtain the post-buckling response of axially compressed NSD cylinders through a semi-analytical solution ................................................................................................ 28 Figure 2-7: (left) Schematics of a cylindrical panel with all edges clamped, (right) Notation for the deformations and stresses in the cylindrical coordinate system [axis, deformation (analytical model), deformation (FE), stress (FE)] ......................................................................................... 30 xi Figure 2-8: Steps of semi-analytical solution for the post-buckling response of a cylindrical shell under compression ........................................................................................................................ 41 Figure 2-9: (a) Schematics of NSD cylindrical segments, (b) schematic distribution of stress prior to buckling event ........................................................................................................................... 42 Figure 2-10: Equilibrium path in snap-back behavior compared to force-displacement from force control and deformation control loading....................................................................................... 44 Figure 2-11: Comparison of force-deformation curves from FE and analytical models for four cylindrical panels with all-clamped boundary conditions ............................................................ 47 Figure 3-1: Cylindrical segments from NSD cylinder in which the stiff area gets bigger from (a) to (c) .............................................................................................................................................. 49 Figure 3-2: Comparing the response of 3 cylindrical segments with the stiffened pattern in Figure 3-1and dimensions in the Table 3-1 (a)the difference in radial deformation of the segment, (b) the difference in force-displacement response, (c) the difference in energy dissipation in the buckling event of each segment, (d) the variation of axial deformation corresponding to the buckling event, and (e) the difference in force drops from buckling in each of the segments ............................... 51 Figure 3-3: Replacing the segment with equivalent panel plus linear springs ............................. 53 Figure 3-4: The deformation of the system of the springs in series when the middle spring buckles. (a) right before buckling, the system is in balance and force in the springs are equal, (b) right after buckling, the total deformation does not change, force in the buckled panel and linear springs in series are not in balance, and (c) the axial deformation in linear springs and buckled panel are adjusted to keep the force in balance ............................................................................................ 57 Figure 3-5: The work lost in the snapback behavior..................................................................... 58 Figure 3-6: The energy dissipated, and the work done in the NSD segment in the moment of buckling......................................................................................................................................... 59 Figure 3-7: The force-displacement curve of segment (a) from finite element analysis compared to the semi-analytical model ............................................................................................................. 60 Figure 3-8: The fall in the internal energy of segment (a) from finite element analysis compared to the semi-analytical model ............................................................................................................. 61 Figure 3-9: The force-displacement curve of segment (b) from finite element analysis compared to the semi-analytical model ......................................................................................................... 61 Figure 3-10: The fall in the internal energy of segment (b) from finite element analysis compared to the semi-analytical model ......................................................................................................... 62 Figure 3-11: The force-displacement curve of segment (c) from finite element analysis compared to semi-analytical model ............................................................................................................... 62 xii Figure 3-12: The fall in the internal energy of segment (c) from finite element analysis compared to the semi-analytical model ......................................................................................................... 63 Figure 3-13: (a) 3D printed cylindrical shell and (b) the axial loading test setup ........................ 65 Figure 3-14: Force-deformation response of axially compressed NSD cylinders as predicted by FE and analytical models compared to experiments .......................................................................... 67 Figure 4-1: Axial response of cylindrical panels from ABS material with R = 50 mm and t = 0.5 mm for varying length (l) and height (b);  and  are the global axial stress and strain of the panel, respectively ................................................................................................................................... 70 Figure 4-2: Contour plot s for cylindrical panels from ABS material with t = 0.5 mm and R = 50 mm for varying length (l) and width (b): (a) buckling strain (0), (b) buckling stress (0), (c) stress drop (), and (d) post-buckling stiffness () ............................................................................... 71 Figure 4-3: (a) Axial response of designed NSD cylinder for local buckling at the axial shortening of 0.24, 0.32, and 0.4 mm, (b) radial deformation of the cylinder at each buckling event ........... 73 Figure 4-4: Axial strain corresponding to the first buckling event for the panels with t/R varies from 0.001to 0.013 ........................................................................................................................ 81 Figure 4-5: Axial stress corresponding to the first buckling event for the panels with t/R varies from 0.001to 0.013 ........................................................................................................................ 82 Figure 4-6: Normalized drop in stress corresponding to the first buckling event for the panels with t/R varies from 0.001to 0.013 ....................................................................................................... 83 Figure 4-7: Normalized drop in the strain energy corresponding to the first buckling event for the panels with t/R varies from 0.001to 0.013 .................................................................................... 84 Figure 4-8: Normalized secondary stiffness after the first buckling event for the panels with t/R varies from 0.001to 0.013 ............................................................................................................. 85 Figure 4-9: Normalized radial deformation in the center of the panel after the first buckling event for the panels with t/R varies from 0.001to 0.013 ........................................................................ 86 Figure 4-10: Normalized slope of radial deformation after the first buckling event for the panels with t/R varies from 0.001to 0.013 ............................................................................................... 87 Figure 4-11: Normalized von-mises stress after the first buckling event for the panels with t/R varies from 0.001to 0.013 ............................................................................................................. 88 Figure 4-12: Normalized slope of von-mises stress after the first buckling event for the panels with t/R varies from 0.001to 0.013 ....................................................................................................... 89 Figure 4-13: Normalized stress vs. strain response of the equivalent panels ............................... 91 xiii Figure 4-14: Schematic plateau post-buckling response of the NSD cylinder ............................. 98 Figure 4-15: Force-displacement response of the designed NSD cylinder for plateau post-buckling from FEA in comparison to one from the semi-analytical model .............................................. 102 Figure 4-16: The shear resulted from D zones deforming less than the center of the segment .. 103 Figure 4-17: Axial deformation distribution before the first buckling in three segments of the designed NSD cylinder ............................................................................................................... 104 Figure 4-18: Force-displacement response of the designed NSD cylinder for plateau post-buckling from modified FEA in comparison to one from a semi-analytical model .................................. 105 Figure 4-19: Axial deformation distribution before the first buckling in segments three of the (a) designed NSD cylinder, and (b)modified NSD cylinder ............................................................ 106 Figure 4-20: Radial deformation in segment three of (a) primarily designed cylinder, and (b) modified cylinder, before the first buckling in segment one ...................................................... 107 Figure 4-21: Schematic post-buckling response of the NSD cylinder where identical force drops occur with the identical difference in end shortening for buckling events ................................. 108 Figure 4-22: Force-displacement response of the designed NSD cylinder for similar load falls in similar distances from FEA in comparison to one from a semi-analytical model ...................... 110 Figure 4-23: Force-displacement response of the designed NSD cylinder for similar load falls in similar distances from FEA in comparison to one from a semi-analytical model ...................... 111 Figure A-1: Axial force-displacement from FE model and experiment ..................................... 118 Figure A-2: Stress distribution in the cylindrical shell for the maximum load (𝛿 = 0.5 𝑚𝑚) .. 119 Figure A-3: Deformation of the cylinder for the maximum load (𝛿 = 0.5 𝑚𝑚) in; (a) radial direction, (b) circumferential direction, and (c) longitudinal direction ...................................... 120 Figure A-4: Sections in the longitudinal and circumferential direction of the cylinder, which has been investigated more in the following figures ......................................................................... 121 Figure A-5: Sequence of radial deformation in the circumference 40 mm ................................ 121 Figure A-6: Sequence of radial deformation in the circumference of 50 mm ............................ 122 Figure A-7: Sequence of radial deformation in the circumference of 60 mm ............................ 122 Figure A-8: Radial deformation along with section A1.............................................................. 123 Figure A-9: Radial deformation along with section A2.............................................................. 124 Figure A-10: Radial deformation along with section A3............................................................ 124 xiv Figure A-11: Radial deformation along with section B1 ............................................................ 125 Figure A-12: Radial deformation along section B2 .................................................................... 125 Figure A-13: Radial deformation along section B3 .................................................................... 126 Figure A-14: Radial deformation along section C1 .................................................................... 126 Figure A-15: Radial deformation along with section C2 ............................................................ 127 Figure A-16: Radial deformation along with section C3 ............................................................ 127 Figure A-17: The longitudinal displacement showed in the circumferential direction for h=40 mm. The radius of the circle here is 1 mm, and it is just to show the magnitude of changes in the longitudinal displacement ........................................................................................................... 128 Figure A-18: The longitudinal displacement showed in the circumferential direction for h=50 mm. The radius of the circle here is 1 mm, and it is just to show the magnitude of changes in the longitudinal displacement ........................................................................................................... 129 Figure A-19: The longitudinal displacement showed in the circumferential direction for h=60 mm. The radius of the circle here is 1 mm, and it is just to show the magnitude of changes in the longitudinal displacement ........................................................................................................... 129 Figure A-20: Axial deformation in section A1 ........................................................................... 130 Figure A-21: Axial deformation in section A2 ........................................................................... 131 Figure A-22: Axial deformation in section A3 ........................................................................... 131 Figure A-23: Axial deformation in section B1 ........................................................................... 132 Figure A-24: Axial deformation in section B2 ........................................................................... 132 Figure A-25: Axial deformation in section B3 ........................................................................... 133 Figure A-26: Axial deformation in section C1 ........................................................................... 133 Figure A-27: Axial deformation in section C2 ........................................................................... 134 Figure A-28: Axial deformation in section C3 ........................................................................... 134 Figure A-29: The longitudinal strain showed in the circumferential direction for h=40 mm. The radius of the circle here is 0.01, and it is just to show the magnitude of changes in the longitudinal strain ............................................................................................................................................ 135 Figure A-30: The longitudinal strain showed in the circumferential direction for h=50 mm. The radius of the circle here is 0.01, and it is just to show the magnitude of changes in the longitudinal strain ............................................................................................................................................ 136 xv Figure A-31: The longitudinal strain showed in the circumferential direction for h=60 mm. The radius of the circle here is 0.01, and it is just to show the magnitude of changes in the longitudinal strain ............................................................................................................................................ 136 Figure A-32: Axial strain in section A1 ...................................................................................... 137 Figure A-33: Axial strain in section A2 ...................................................................................... 137 Figure A-34: Axial strain in section A3 ...................................................................................... 138 Figure A-35: Axial strain in section B1 ...................................................................................... 138 Figure A-36: Axial strain in section B2 ...................................................................................... 139 Figure A-37: Axial strain in section B3 ...................................................................................... 139 Figure A-38: Axial strain in section C1 ...................................................................................... 140 Figure A-39: Axial strain in section C2 ...................................................................................... 140 Figure A-40: Axial strain in section C3 ...................................................................................... 141 xvi Chapter 1. Introduction 1.1. Motivation The viewpoint shift of considering elastic instabilities in thin-walled structural elements to be sources for productive behavior rather than failure is relatively recent [1], [2]. Applications for elastic instabilities in thin-walled elements include stress sensing, smart switches, deployable structures, energy harvesting, and energy dissipation [3], [4], [13]–[16], [5]–[12]. Buckling under axial compression is an unstable response with the potential for switching applications if the structure's deformations are adequately limited. A constrained slender strip between rigid walls is an example of such behavior that was recently used to harvest energy from multiple buckling events [3], [7], [17] (See Figure 1-1). However, it is advantageous to eliminate the need for external constraints and obtain numerous buckling events in the structure. Short thin-walled cylinders can experience several elastic buckling events under axial compression without additional limitations because of their self-confined configuration [18]–[21]. However, thin-walled cylindrical shells are infamous for their unpredictable response to axial compression due to their high sensitivity to minor initial imperfections [22]–[27]. Cylindrical shells have been the subject of numerous studies which tried to predict an acceptable design load to avoid undesirable buckling events. 1 Figure 1-1: (a) Slender beam under compression buckles and deflection keeps growing until damage occurs, (b) the same slender beam under compression but constraint between rigid wall shows several buckling events, (C) axially compressed cylindrical shells with uniform thickness show several buckling events More recently, the concept of cylinders with non-uniform stiffness distribution (NSD) was developed to limit the buckling events in the cylinder to specifically targeted panels [28](See Figure 1-2). This strategy was shown as an effective way of defining the number, location, and order of local buckling events. Yet, no analytical or semi-analytical method were proposed to 2 predict the post-buckling response of NSD cylinders, let alone designing them for a targeted post- buckling response. An example of targeted post-buckling responses is the pre-defined spacing of the buckling events, and the size of the load or energy drop in each buckling event. Figure 1-2: Cylindrical shells with nonuniform stiffness distribution (NSD cylinder) limits the buckling in the cylinder to specific areas The main goal of this research is the development of an analytical framework to predict the post-buckling response of NSD cylinders and to design them for a desirable post-buckling response The analytical framework is based on three steps: I) dividing the cylinder into individual panels 3 assumed to be connected in parallel and series, II) predicting and designing the response of each panel, and III) integrating the response of individual panels. The predictive model and design procedure are expected to be effective for the cylinders which are providing strong boundary conditions for independent panel response. In addition, the proposed semi-analytical model can be used to design novel cylinders, such as cylinders with non-uniform cylindrical stiffeners. Such a cylinder follows the rules of NSD cylinders while increasing the number of buckling events and potentially taking care of the cylinder's post-buckling stiffness. 1.2. Hypothesis The proposed analytical framework in this study predicts the post-buckling response of cylinders by considering three integrated steps: I) dividing the cylinder into individual panels connected in parallel and series, II) predicting and designing the response of the individual panel, and III) integrating the response of individual panels. Cylindrical shells with a targeted post- buckling response can be designed using the analytical framework. 1.3. Research significance The elastic post-buckling response of short thin-walled cylindrical shells is not a commonly considered structure for smart applications because of its unpredictability. This research will provide a powerful tool for I) predicting the elastic post-buckling response of NSD cylinders and II) designing NSD cylinders with a targeted elastic post-buckling response. Being able to design cylinders with predictable elastic post-buckling responses will open new opportunities to develop devices with switchable applications that do not rely on external constraints. 4 1.4. Background 1.4.1. Buckling in cylinders Structural instability is a well-studied field due to the resulting catastrophic effects on traditional structures. However, unstable response in structures has recently regained increased attention due to newly emerged smart applications, such as vibration control, smart switches, additive manufacturing, deployable structures, energy dissipation, and energy harvesting [3], [4], [13]–[16], [5]–[12]. Buckling under axial compression is an unstable response that releases energy with high-rate motion. If this fast deformation is not adequately limited, the large deformations will result in sudden failure of the structure [29]. Axial compression may lead to global buckling for some structures like slender beams and columns. However, controlled buckling, or instability, could be a desired feature rather than a failure limit for specific elements or device design. Managing post- buckling response, which is always associated with a rapid increase of deformations perpendicular to the loading axis, requires external constraints [30], [31]. These constraints confine the unrestrained increase of transverse deformations and allow the formation of new stable states for the buckled structure. For example, constraining the global buckling of slender beams permits novel use of the dynamic unstable transition events to harvest electrical energy for self-powered sensors [3], [7], [17]. It is advantageous to eliminate the need for external constraints when aiming to obtain a post- buckling response with multiple buckling events. Due to their self-confined geometry, short, thin- walled cylindrical shells can experience numerous local elastic buckling events when subjected to axial compression [18]–[21]. Shell structures are widely used in different engineering fields, such as offshore structures and aerospace vehicles. However, their response to axial compression is well 5 known to be highly sensitive to initial imperfections, and their buckling strength never reaches the theoretical limit [22]–[27]. Two approaches are taken to avoid consequences from this uncertainty for traditional uses of cylindrical shells where the buckling is a catastrophic event. The first approach is to lower the design load to the point that there is no chance for the cylinder to buckle in the design forces. A second approach consists of stiffening the cylinder through hoops and stringers, which improves the axial strength of the cylinder and postpones its global buckling to higher demands. These approaches have been taken individually or simultaneously in several studies for different materials, from aluminum to composite materials produced in different ways. However, none of these approaches try to eliminate the sensitivity of the post-buckling response to initial imperfections. 1.4.2. Hoop and string stiffeners One standard method is to stiffen the shell using hoop and string stiffeners to avoid buckling in the cylinder before the expected design load. Most studies smear the stiffeners onto the surface of the cylinder and treat the structure as a cylinder with orthotropic material properties to simplify the calculation [32], [33]. Yet, local buckling in the areas between stiffeners still takes place. The buckling events in the stiffened panels are still sensitive to initial imperfections [34]–[36]. Several studies have proposed semi-analytical models to predict the axial response of a single panels in stiffened cylinders considering flexible boundary conditions [34], [37]. These models showed good accuracy in predicting the pos-buckling response of stiffened panels for a specific initial imperfection. Results from these models were compared to finite element simulations with the same imperfections. However, this is not what happens in stiffened panels, as several researchers have warned about the impact of local buckling events in the stiffened panels on the global buckling of the cylinder [38]–[40]. Although stiffening of a cylinder with hoop and stringers 6 increases its axial load capacity, it does not eliminate the uncertainty in the load capacity regarding the impact of initial imperfections. For several applications, this uncertainty in post-buckling response is accepted, and guidelines address this issue for design purposes. 1.4.3. Buckling capacity knockdown factors To deal with the uncertainty in determining the buckling capacity of cylinders, a knockdown factor (KDF) is defined to decrease their design strength to a reliable level [41]–[44]. The definition of early KDF values was primarily empirical and very conservative [41], [45], [46]. Among these studies, the recommendations by NASA SP-8007 [41] have been used as the design guideline for many years. The knockdown factor in this guideline is based on the lower-bound axial strength from numerous experiments from the 1950s (Figure 1-3). Figure 1-3: Empirical knockdown factors for the design of cylindrical shells Many research studies urged to update NASA’s guideline for modern cylindrical structures manufactured with more accurate technology. NASA launched the SBKF (Shell Buckling Knockdown Factor) project [47]–[50] to adjust cylinder designs with new materials and 7 manufacturing methods. Yet, to the author's knowledge, no design guidelines from this project are available. Finite element analysis is another way to define reliable KDFs for axially compressed cylindrical shells. Since there are infinite possible initial imperfections for cylinders, this approach aims to find the best regime of applying imperfections to develop reliable KDFs using imperfection patterns. The Single Perturbation Load Approach (SPLA), Multi Perturbation Load Approach (MPLA), Geometric Dimple Imperfection (GDI), and Linear Buckling Mode-shaped Imperfection (LBMI) are the most well-known approaches to define imperfection patterns [51], [52]. In an endeavor to design optimum cylinders with desired buckling load, Hao et al. developed a repetitive algorithm [53]. Their algorithm divides the stiffened cylinder into several connected panels, replaces the hoop and spring stiffened panels with equivalent uniform thickness panels, and optimizes the cylinder's weight for specific collapse load. Their results showed promise to decrease the mass of the hardened cylinder while keeping the knockdown factor the same. However, they did not address the boundary condition between panels, and they did not decrease the sensitivity of the cylinder to initial imperfections. The methods mentioned above provide acceptable confidence that the cylinder will not buckle before a specific force. However, thin-walled cylinders for smart applications require knowledge of the exact force, displacement, and location of the multiple local buckling events generated on their surface. During the last century, numerous studies have been conducted to determine methods to predict the post-buckling response of axially compressed short thin-walled cylindrical shells for a given initial imperfection. The approaches can be considered of two kinds. One kind consists of semi- analytical models based on classical shell theory, where two partial differential equations express 8 equilibrium and compatibility conditions in the cylinder [54]. In this approach, radial deformation is simplified as the sum of trigonometric functions that satisfy the boundary conditions of the cylindrical shell. Then, the partial differential equations are transformed into nonlinear algebraic equations using Galerkin's approximation method. The second kind uses the finite element formulation to develop nonlinear equations for the system directly. From this point, both approaches use numerical techniques to solve the corresponding equations. Some researchers have used Newton's method and solved the equations for increments of maximum deflection of the shell or deflection of the shell at its center (only for cylindrical panels). However, this method does not accurately predict the snap-back part of the response. In 1979 Riks [55] developed the arc-length method where the force and axial deformation vary at the same time along a circular arc from the previous step's solution. The original arc-length method had some deficiencies when the structure had a snap-back response [56], [57]. The biggest issue is the possibility of going back instead of going forward. Since the development of the Riks method, several modifications have been developed to solve the problem mentioned above [58], [59]. 1.4.4. Imperfection insensitive cylinders Despite the increased ability to predict the response of an axially compressed cylindrical shell for given initial imperfections, the control of buckling and the post-buckling response of such structures remains unsolved. Researchers have recently tried to control cylindrical shells' buckling and post-buckling response by manipulating the cylinder's shape to trigger specific buckling modes. Jullien and Araar patented Aster shells (Figure 1-4(a)), which consist of the inclusion of multiple arches on the circumference of the cylinder uniformly distributed over the height of the cylinder [60]. Their primary goal was to improve the strength of cylindrical shells under external pressure. 9 Figure 1-4: (a) section and 3D shape of Aster cylinder, (b) section and 3D face of the corrugated cylinder by Ning and Pellegrino Ning and Pellegrino optimized Aster shells by varying the size of the arches to design a cylinder that is insensitive to initial imperfections under axial compression[61], [62]. They showed that increasing the amplitude of the waves boosts the critical buckling load while also decreasing the effect of initial imperfections (Figure 1-4(b)). Yadav and Grasimidis expanded the study on the imperfection insensitivity of the cylinders with wavy sections to the combination of bending and compression [63]. Using finite element analyses, they showed that the buckling force in wavy sections is not changed significantly for several scenarios of initial imperfections. Such a method is an excellent tool to improve the strength of cylinders for structural purposes. For example, Yadav and Grasimidis studied the effectiveness of this type of imperfection insensitivity in wind turbine towers [64]. Their study on two wind turbine towers showed significant insensitivity to initial imperfection for the cylindrical segments of the towers and the conical shells that connect cylindrical pieces. However, none of these studies investigate the post-buckling response of such 10 cylinders, and in several cases material nonlinearity occurred before buckling events [61]. Results from finite element analyses showed that the wavy section of the cylinder is not capable of multi- stable response while it increases the axial strength significantly. In another approach, Cox et al. [65] proposed to "nudge" the cylindrical shell to control the post-buckling response of their cylindrical panels. Using finite element analyses, they showed that nudging the cylindrical shells with the combination of the first few buckling modes has a remarkable positive effect on controlling the post-buckling response of the cylindrical panels without a significant increase in the mass of the shell. Manipulating the stiffness to control the buckling and post-buckling response of cylindrical shells is not limited to cylinders with isotropic materials. As the application of fiber-reinforced composite materials increases, these materials are being used in cylindrical elements. However, the orthotropic nature of composite materials does not banish the sensitivity of the post-buckling response to initial imperfections, just as for thin- walled cylinders from homogenous materials Recently, Lincoln et al. [66] reduced the imperfection sensitivity of a cylinder's axial strength by using the continuous tow shearing (CTS) method to manufacture variable angle composites. The CTS manufacturing method enhances the thickness of the composite cylinder by gradually increasing the fiber angle (Figure 1-5). Thickness variation is the main reason for the observed decline in the sensitivity to initial imperfections. 11 Figure 1-5: The thickness of the wall varies in CTS cylinders 1.4.5. Cylindrical shells with non-uniform stiffness distribution Despite the extensive research on predicting and controlling the strength of axially compressed cylinders (first buckling event), controlling the far post-buckling response of such structures was not appropriately studied. Recently, Hu and Burgueño, motivated by smart applications of structural instabilities, endeavored to manage the post-buckling response of axially compressed cylindrical shells by tailoring the stiffness of the cylinder. They proposed different types of seeded imperfections on cylindrical shells to manipulate stiffness in the desired way [28], [67]–[70]. The approaches included the seeding of sinusoidal geometric imperfections to the cylinder's wall (Figure 1-6(a)), changing the pattern of the shell's wall stiffness Figure 1-6 (b)), and the provision of external and internal radial constraints (Figure 1-6(c)). These approaches were studied using finite element analysis (FEA) and experiments on 3D printed specimens. Results from the experiment and FEA showed that seeding modal imperfections or providing a patterned non-uniform wall stiffness distribution decreased the sensitivity of the post-buckling response to initial imperfections significantly. On the other hand, cylinders with external constraints, internal constraints, or both did not show any advantage over ordinary cylinders regarding the predictability of the buckling events unless one of the other two methods was used 12 to seed imperfection in the cylinder's surface. However, constraints for transverse deformations can increase the number of buckling events due to mode switches in the cylinder [71]. Figure 1-6: Seeding imperfections by Hu and Burgueño; (a) sinusoidal seeded imperfection, (b) changing the pattern of shell's wall thickness, and (c) radial constraints [60] Among the methods for seeding imperfections, changing the shell's wall stiffness showed the most promise over the control of local buckling areas. Hu and Burgueño showed that different patterns of stiffening the shell's wall lead to localized buckling in pre-defined slender regions [71]. However, the stiffening pattern in some of the tailored cylinders allowed deformation continuity between the slender areas, which caused the local buckling events to affect each other. Further, their work showed that providing a specific tailored non-uniform stiffness distribution to the cylinder's wall (NSD cylinder), achieved by different thickening regions of the cylinder, allowed control of the number, location, and sequence of local buckling events in the cylindrical shell [28]. 13 Figure 1-7: Buckling areas in cylinders with non-uniform stiffness distribution and load drops in the axial response of the cylinders [60] Figure 1-7 shows the primary proof that using an NSD cylinder makes it is possible to relate the drops in the axial force to the buckling in thin regions. Also, the response from loading and unloading in Figure 1-7 looks very promising for energy dissipation applications. In a follow-up study, Burgueño and Guo presented a model describing the effect of localized wall thickening on the strain energy released due to the localized buckling events, manifested as load drops in the cylinder's global force-deformation response [72]. The distribution of stiffness in their study divided the cylinder into parallel segments and stiffened areas at the top and bottom of each segment. Also, they showed it does not matter which segments are neighbors, and the post- buckling response will be the same for all the NSD cylinders composed of specific segments. 14 The study by Burgueño and Guo showed promising results regarding the sequence of buckling events and the likely area of buckling in each segment. However, the approach did not permit predicting the exact force, displacement, and force drop for a given thickening pattern. Such a prediction will open the door for many smart applications that require knowing in detail the buckling events in the cylinder. Harvesting electrical energy from the post-buckling response of cylindrical shells is a potential application that recently received more attention [73]–[77]. Several studies tried to formulate and predict the post-buckling response of functionally graded (FG) with piezoelectric layer thin cylindrical shells for energy harvesting. However, their endeavor still lacks predictability since the response of FG cylindrical shells is highly dependent on initial imperfections. NSD cylinders look promising in that area since they provide some degrees of predictability for the FG cylinders [78]. More recently, M. Nazar et al. [79] tried to harvest electrical energy from cylindrical shells with specific stiffness distribution. They attached the piezoelectric coated PVDF film to the surface of the cylinder. Their results showed some jump in voltage during the buckling process, while the post-buckling response from the experiment was far from the predicted response from FE and the experiment on the cylinder without piezoelectric film. They did not study the effect of attached piezoelectric film on a cylinder's body while using finite element analysis to predict the mechanical response of the cylinder. However, it also shows the sensitivity of the cylinder's response to initial imperfections. 1.4.6. Research gap The presented background shows that the buckling capacity of axially compressed cylindrical shells has been studied extensively for structural and mechanical structures. However, the endeavor has concentrated on delaying global buckling or designing cylinders that do not buckle 15 before a particular load. Recently, studies have tailored the cylinder's geometry to eliminate, or reduce, imperfection sensitivity. However, the aim has been to increase the maximum possible axial strength of the cylinder, and thus these efforts did not consider the cylinder's response beyond the first buckling event. Only a few research projects have been conducted to control the post-buckling response of cylinders for smart applications. Among these studies, non-uniformly stiffness distributed (NSD) cylinders are very promising based on previous studies in our research group. We understood that it is possible to treat buckling segments of an NSD cylinder separately since the segments do not affect each other. Yet, there is not an established framework to predict the post-buckling response of an NSD cylinder without modeling the entire cylinder and using the finite element method. More importantly, there is no design procedure for an NSD cylinder with a desirable far post- buckling response. 1.5. Research plan This study aims to develop analysis and design methods that allow gaining control of the post- buckling response of axially compressed NSD cylindrical shells. This control is two-fold: I) predict the response of any given NSD cylinder, and II) design NSD cylinders for a given post-buckling response. The current NSD cylinder concept is based on thickening areas non-uniformly to limit buckling events to specific zones on the cylinder. The following steps will be followed to achieve the objectives of the proposed research: Task I. Develop a semi-analytical framework to predict the post-buckling response of axially compressed NSD cylinders. (Figure 1-8 left) 16 - Developing the framework's structure: outline the steps of the analysis framework to separate the NSD shell into segments, simplify the segments into equivalent panels with uniform thickness, and reintegrate the response of the NSD cylinder. - Use finite-element analysis to predict the post-buckling response of every individual cylindrical panel. Task II. Develop a procedure to design an NSD cylinder for favorite post-buckling responses. (Figure 1-8, right) - Developing design maps that lets the designer select appropriate geometries for equivalent panels in each segment. - Develop a procedure to design an NSD cylinder from design maps. Task III. Design NSD cylinders with a desired post-buckling response. (Figure 1-8, right) - Design NSD cylinders: to show the validity of the design process and study the post- buckling capabilities of NSD cylinders. - Design an NSD cylinder for a targeted post-buckling stiffness. - Design an NSD cylinder for targeted load drops in specific end shortenings in the buckling events. 17 (a) (b) Figure 1-8: Schematic research plan (a) predict the post-buckling response of any given NSD cylinder and (b) design an NSD cylinder for a given post-buckling response 18 Chapter 2. A framework to predict the post-buckling response of NSD cylinders Nan Hu and Burgueño [71] showed that seeding imperfection could force thin-walled cylindrical shells to buckle in a certain way. As one of the strategies, they introduced the concept of cylindrical shells with non-uniformly distributed stiffness (NSD cylinders) to control the post- buckling response. Here, a semi-analytical framework to predict the axial response of NSD cylinders is presented. The framework has three main steps: 1) dividing the cylinder into individual panels assumed to be connected in parallel and series, 2) predicting the axial response of an individual panel with a semi-analytical model, and 3) integrating the response of the individual panels. Figure 2-1: Wall-thickness pattern along the surface of a cylinder with Non-uniform Stiffness Distribution (NSD); h is the cylinder height, and t is the wall thickness of each area 2.1. Dividing NSD cylinder into parallel segments The buckling events in NSD cylinders are limited to slender zones (areas with a thickness equal to 𝑡𝑏 in Figure 2-1), where thicker zones (areas with a thickness equal to 𝑡2 in Figure 2-1) impede interaction between the slender zones. It is assumed that the parallel zones of NSD cylinders are isolated from each other, and thus the cylinder can be divided into parallel panels. Figure 2-2 shows the separated parts of the cylinder and their boundary conditions. These assumptions are based on observations from experiments on NSD cylinders where the deformation on the divider area were 19 visibly negligible and from a detailed study on finite element simulations for a sample NSD cylinder. Figure 2-2: Division of the NSD cylindrical shells into parallel independent panels 2.1.1. Boundary conditions As the cylinder is divided into parallel segments, their boundaries must be restrained in a way that resembles the behavior of the segments in the actual cylinder. The top and bottom boundaries of an individual segment are equal to the boundaries of the actual cylinder. However, the sides of each segment are now experiencing a new condition that must be addressed. There is no concern about the boundaries of segment D in Figure 2-2, as its axial response remains linear during the loading period. On the other hand, the side boundaries of the segments subjected to buckling (A, B, and C in Figure 9) are connected to the thick areas in segment D. Connection to the thicker area provides some degree of rigidity. Here, it is assumed that the thicker area in zone D is at least two times thicker than the body of the cylinder. Doubling the thickness increases the moment of inertia 20 of the section by 8 times and the bending moment from the thin zone during the buckling event has a negligible bending effect on zone D compared to the deformation that would be generated in this zone if it had the same the thickness as the thinner segments. Thus, it is assumed that side boundaries are clamped but can move along the axial direction. Figure 2-3: Sections in the longitudinal and circumferential direction of the cylinder to investigate the boundary conditions in buckling segments A finite element model for an NSD cylinder was developed to examine the assumption of rigid boundaries. The cylinder’s radial deformations during compression were investigated in several sections from the FE model. Figure 2-4 shows the radial deformation at the cylinder’s mid-section. The obtained responses supported the assumption of rigid boundaries on the sides since the radial deformations and rotations at the edge of the stiff regions are both nearly zero. The investigation on boundaries considered the deformations from several sections in the cylinder. These deformations, in addition to details of the finite element model, are shown in Appendix of this document. 21 Figure 2-4: Sequence of radial deformation in the circumference of 50 mm To load the cylinder uniformly during experiments, the top and bottom edges are connected to a thick surface inserted into a stiff grooved plate. This configuration leads to clamped boundary conditions for the top and bottom cylinder edges. Also, the top and bottom of segments A, B, and C are thickened by nature of the NSD cylinder design. Therefore, applying the same displacement on the surface means avoiding bending at the ends, which means that naturally clamped boundary conditions are provided at the top and bottom edges of the cylinder. 2.1.2. Divider segments The cylinder is divided into parallel segments and compressive loading is applied in their longitudinal direction. Segments D, which have a small thin region at their top and bottom ends, are not expected to buckle and, as a result, can be considered as linear axial springs. Since a uniform deformation is applied to the cylinder’s top edge, all the segments experience the same deformation at the top. The total axial load in the cylinder is therefore the sum of forces in the 22 segments, and they can be modeled as parallel springs. The stiffness of the linear spring of segment D is calculated with Equation (1) 1 𝑙𝐷 ℎ − 𝑙𝐷 𝐸𝑡𝐷 𝑡𝑏 𝑏𝐷 = + ⟹ 𝑘𝐷 = (1) 𝑘𝐷 𝐸𝑡𝐷 𝑏𝐷 𝐸𝑡𝑏 𝑏𝐷 𝑙𝐷 𝑡𝑏 + ℎ𝑡𝐷 − 𝑙𝐷 𝑡𝐷 2.1.3. Equivalent Panels Segments with potential buckling regions (A, B, and C in Figure 2-2) have two thick zones at the top and bottom, and a thin zone at the center – the targeted buckling region. Thus, the main problem is to predict/design the post-buckling response of a cylindrical panel with non-uniform thickness under compression. Since there is no fully established model for the panels with non- uniform stiffness distribution, the only way to accurately predict the post-buckling response of such a cylindrical panel is using the finite element method (FEM). Although it is helpful to use the FEM to predict the post-buckling response of the NSD panels, a finite element (FE) model of the exact panel will not help with the design process where the panel is yet to be determined. Also, finite element analysis (FEA) is inefficient for functionally graded cylinders, and the amount of computational time required to deal with contact between material layers is significant. It was thus decided to simplify the NSD segment, which has non- uniform stiffness, with a uniform stiffness panel and then use a semi-analytical model to predict its post-buckling response under compression. Replacing the NSD segment with an equivalent panel with uniform thickness facilitates using FEA for post-buckling response prediction. However, the special geometry of stiffness distribution in the segment makes it totally different from a uniform panel with a regular distribution of radial displacements. The shape of the segment forces the buckling event to occur at the center of the segment even without the explicit definition of an initial imperfection. This is beyond considering larger initial imperfections at the center of the cylindrical panel. Increasing the amplitude of 23 imperfection at the center of the panel smoothens the force-displacement curve of the segment to the point that the force drop due to buckling disappears. Given the considerations noted above, a semi-analytical model that addresses the simulation needs of NSD segments of a cylinder was developed. The proposed method applies to regular panels by considering some adjustments. The special treatment for NSD segments is discussed later in this chapter. First, the NSD segment with an equivalent panel is replaced with one of uniform thickness. The equivalent length of the panel is calculated with the assumption of constant width. The approach to determine the equivalent width of the panel is presented in Chapter 3. Figure 2-5: Replacing thick zones (t=t1i) in the panel with an equivalent length of the shell with identical thickness to mid-section (t=tb) The wall thickness of the top and bottom zones is set equal to the wall thickness of the center zone by modifying its length to simplify the panel with non-uniform thickness. This idea has been used by others to determine the buckling capacity of cylinders with stepwise wall thickness under pressure [80]–[82]. 24 The equivalent length of the thick-walled regions at the top and bottom of the panel should be determined so that the transverse stiffness at the edge of the thin area is unchanged. The equivalent length of the shell with the transformed thickness (top and bottom regions) is calculated in two general steps. First, the ticker regions at the top and bottom of the panels are assumed as one-way cantilevers and their transverse tip displacement, , under an applied tip moment, M, is determined (see Figure 2-5). Second, the length of a cantilever shell with a thickness of the mid-section of the panel (𝑡𝑏 in Figure 2-5) that yields the same deformation  under an equal bending moment M is determined. The bending moment and deformation at the edge of the cantilever elements are related by Equation (2): 3𝐸𝐼1𝑖 3𝐸𝐼𝑏 𝑀= 𝛥= 𝛥 (2) 𝑙12𝑖 𝑙 ′12𝑖 where 𝑖 indicates the segment’s name (A, B, or C), E is the modulus of elasticity, l1i is the initial length of thicker areas at the top or bottom, l’i1 is the equivalent length of the same area, 𝑡𝑏 and 𝑡1𝑖 are the thicknesses of the segment’s mid-section and the thicker areas; M is the boundary's bending moment and ∆ is the edge's deformation in the radial direction. The equivalent length of the thicker zones is then calculated with Equation (3). 𝐼1𝑖 𝐼 ′1𝑖 𝑡𝑏 3 = ′ ⇒ 𝑙 ′1𝑖 = ( )2 𝑙1𝑖 . (3) 𝑙1𝑖 𝑙 1𝑖 𝑡1𝑖 Equation (3) shows a simple solution for finding the equivalent length for the NSD cylindrical segments. NSD cylindrical segments do not fulfill the assumption of uniform thicker areas at the top and bottom of a cantilever beam behavior. However, these assumptions were considered acceptable to develop the proposed framework. Chapter 3 presents a detailed discussion on the calculation of the panel’s equivalent length. 25 2.2. Composition axial responses of separated segments The geometrical modifications described in Section2.1 permit modeling cylindrical NSD segments as uniform cylindrical panels under compression, which simplifies determining an analytical solution to their post-buckling response. However, the equivalent shell panel has a different axial stiffness than the original one. To correct the axial stiffness deviation, springs are added in series to the equivalent panel to match the initial panel stiffness. This spring is to respond linearly under axial compression since the equivalent panel is responsible for the nonlinear part of the response. The stiffness of this spring is calculated by matching the axial stiffness of the original NSD cylinder. Considering the dimensions of the cylinder from Figure 2-2 and the length of the equivalent panel as 𝑙1′ 𝑖 , the stiffness of the linear spring is calculated with Equation (4): 1 1 1 = − 𝑘𝑖1 𝑘𝑠𝑒𝑔𝑚𝑒𝑛𝑡𝑖 𝑘𝑝𝑎𝑛𝑒𝑙𝑖 1 ℎ − 2𝑙1𝑖 2𝑙1𝑖 − 2𝑙2 2𝑙2 = + + 𝑘𝑠𝑒𝑔𝑚𝑒𝑛𝑡𝑖 𝐸𝑡𝑏 𝑏 𝐸𝑡1𝑖 𝑏 𝐸𝑡1𝑖 𝑏1 + 𝐸𝑡𝑏 (𝑏 − 𝑏1 ) (4) 𝐸𝑡𝑏 𝑏 𝑘𝑝𝑎𝑛𝑒𝑙𝑖 = ′ 𝑙1𝑖 𝑏𝐸𝑡1𝑖 𝑡𝑏 (𝑏1 𝑡1𝑖 + 𝑏𝑡𝑏 − 𝑏1 𝑡𝑏 ) 𝑘𝑖1 = 𝑏𝑡𝑏 ((𝑙1′ 𝑖 − 𝑙2 )𝑡1𝑖 − 𝑙1𝑖 𝑡𝑏 ) − 𝑏1 (𝑡1𝑖 − 𝑡𝑏 )(𝑙1′ 𝑖 𝑡1𝑖 + 𝑙1𝑖 𝑡𝑏 ) The total axial displacement in an individual panel is equal to the sum of the displacements from the linear springs and the displacement from the nonlinear response of the panel. Similarly, the axial force corresponding to a given displacement in the cylinder is the sum of axial forces in all segments (A, B, C, and three D panels) at that displacement. In other words, the segments contribute to the response of the cylinder as resisting elements in parallel, while different regions 26 (i.e., different thicknesses) within an individual segment behave like elements in series. Figure 2-6 shows the framework followed to predict the elastic post-buckling response of the NSD cylindrical shell in Figure 2-1. The stiffness of all the equivalent linear springs in the cylinder are known, except the axial nonlinear response of the equivalent cylindrical panels 𝐾𝐴𝑛𝑙 , 𝐾𝐵𝑛𝑙 , 𝐾𝐶𝑛𝑙 . In Section 0, a classical method is used to develop a relationship for the nonlinear response of an individual cylindrical panel with all-clamped boundary conditions. 27 Figure 2-6: Framework to obtain the post-buckling response of axially compressed NSD cylinders through a semi-analytical solution 28 Two aforementioned simplifications are the basis for the proposed framework: 1) diving the NSD cylinder into parallel isolated panels, and 2) replacing the non-uniformly stiffened segments with an equivalent uniform panel. The first assumption was validated for several NSD cylindrical shells by Guo et al. [72], while the second simplification had not been evaluated before this work. It seems inevitable to simplify a non-uniform stiffened cylindrical segment with one with uniform stiffness. However, the stiffened parts of the segments have several variables that could affect the length of the equivalent uniform panel. For instance, the ratio of 𝑏1 ⁄𝑏 , and 𝑙1 /𝑙2 in Figure 2-1 are two geometric factors that potentially affect the length of equivalent panel. A study of the geometry of the stiffened area and its impact on the post-buckling response of the segment is presented in Chapter 3. 2.3. A semi-analytical model for axially compressed cylindrical panels In Section 2.1 the concept of separating panels and transforming them into equivalent ones with uniform thickness was discussed. This section presents the analytical solution for the elastic post-buckling response of a single panel with clamped boundaries at the top and bottom edges and a guided roller on the sides. The resulting force-displacement curves are compared to FEA results for individual panels with equal boundary conditions. The final results from the semi-analytical model and FEA are compared to experimental results for several cases. The similarity between all three results validates the idea of separating panels of the NSD cylinder and recomposing the response from that obtained for individual panels. At first, the semi-analytical model for the cylindrical panel looks impractical when there is access to FEA programs like ABAQUS or ANSYS [87]. Nevertheless, in reality, because of the specific shape of the thickened areas in NSD cylinders, specific buckling deformations are being 29 dictated to the individual panels, which are very hard to predict by an FEA software. Also, the semi-analytical model is very useful for situations without access to commercial FEA programs. But as mentioned in Section 0, most importantly, the semi-analytical solution is a better option when we are dealing with a functionally graded cylindrical panel where FEA programs need significant time to solve contact between layers. 2.3.1.1. Analytical Solution for cylindrical panel with all edges clamped The semi-analytical model assumes that the lateral deformation of the cylindrical panel is a combination of pre-defined shapes that are compatible with the boundary conditions. Figure 2-7 (left) shows the panel, the boundary conditions, and the direction of the applied compressive displacement. Notations for the coordinate system, deformations, and stresses in each direction for the analytical and FE models are shown in Figure 2-7 (right). Figure 2-7: (left) Schematics of a cylindrical panel with all edges clamped, (right) Notation for the deformations and stresses in the cylindrical coordinate system [axis, deformation (analytical model), deformation (FE), stress (FE)] 30 According to classical shell theory, the deformation of cylindrical shells should satisfy the following equilibrium and compatibility equations based on von Karman’s theory [21]: 𝑓,𝑥𝑥 𝐷𝛻 4 𝑤 − − [𝑓,𝑦𝑦 (𝑤,𝑥𝑥 + 𝑤,𝑥𝑥∗ ∗ ) − 2𝑓,𝑥𝑦 (𝑤,𝑥𝑦 + 𝑤,𝑥𝑦 ∗ ) + 𝑓,𝑥𝑥 (𝑤,𝑦𝑦 + 𝑤,𝑦𝑦 )] = 0 𝑅 4 𝜕4 𝜕4 𝜕4 (𝛻 = 4 + 2 2 2 + 4 ) (5) 𝜕 𝑥 𝜕 𝑥𝜕 𝑦 𝜕 𝑦 𝜕 𝜕 𝐸𝑡 3 , ( ,𝑥 = 𝜕𝑥 , ,𝑦 = 𝜕𝑦) , 𝐷 = 12(1−𝜈2 ) 𝑤,𝑥𝑥 𝛻 4 𝑓 − 𝐸𝑡 (𝑤,𝑥𝑦 2 − 𝑤,𝑥𝑥 𝑤,𝑦𝑦 − ∗ + 2𝑤,𝑥𝑦 𝑤,𝑥𝑦 ∗ − 𝑤,𝑥𝑥 𝑤,𝑦𝑦 ∗ − 𝑤,𝑦𝑦 𝑤,𝑥𝑥 )=0 (6) 𝑅 In these equations, 𝑅 is the radius of the panel, 𝑤 is the radial deformation, 𝑤 ∗ is the initial radial imperfection, and 𝑓 is the Airy stress function. The Airy stress function definition considers the stresses in the x-direction, y-direction, and x-y plane to be given by: 𝑁𝑥 = 𝑓,𝑦𝑦 𝑁𝑦 = 𝑓,𝑥𝑥 𝑁𝑥𝑦 = −𝑓,𝑥𝑦 (7) The in-plane deformations 𝑢 and 𝑣 are related to out-of-plane deformation 𝑤 and ∗ , and the airy stress function is the following way: 𝐸𝑡 1 𝑤 1 𝑓,𝑦𝑦 = 2 [𝑢,𝑥 + 𝑤,𝑥2 + 𝑤,𝑥 𝑤,𝑥∗ + 𝜈 (𝑣,𝑦 − + 𝑤,𝑦2 + 𝑤,𝑦 𝑤,𝑦∗ )] 1−𝜈 2 𝑅 2 𝐸𝑡 1 𝑤 1 𝑓,𝑥𝑥 = 2 [𝜈 (𝑢,𝑥 + 𝑤,𝑥2 + 𝑤,𝑥 𝑤,𝑥∗ ) + 𝑣,𝑦 − + 𝑤,𝑦2 + 𝑤,𝑦 𝑤,𝑦∗ ] (8) 1−𝜈 2 𝑅 2 𝐸𝑡 −𝑓,𝑥𝑦 = [𝑢 + 𝑣,𝑥 + 𝑤,𝑥 𝑤,𝑦 + 𝑤,𝑥 ∗ 𝑤,𝑦 + 𝑤,𝑥 𝑤,𝑦∗ ] 2(1 + 𝜈) ,𝑦 2.3.1.1.1. The radial deformation functions The solution for Equations.(5) and (6) is based on taking the radial deformation w and the radial initial imperfection w* as functions of coordinates x and y with unknown coefficients in a way that they satisfy the boundary conditions of the panel. It is then possible to find the unknown 31 coefficients by minimizing the integration of Equations.(5) and (6) with respect to them. Considering a clamped boundary condition in all edges with the panel sides allowed to move axially, w and w*should satisfy Equation (9). 𝑤 = 𝑤∗ = 0 (for 𝑦 = 0 & 𝑏) , 𝑓𝑥𝑦 = 0 (for 𝑦 = 0 & 𝑏) (9) Any equation with the following shape can be applied to the Equation (9). 𝑖1 𝜋𝑥 𝑖2 𝜋𝑥 𝑗1 𝜋y 𝑗2 𝜋y 𝑤 = ∑ ∑ ∑ ∑ 𝑎𝑖1𝑖2 𝑗1 𝑗2 (sin ) (sin ) (sin ) (sin ) 𝑙 𝑙 b b 𝑖1 𝑖2 𝑗1 𝑗2 (10) 𝑖1 𝜋𝑥 𝑖2 𝜋𝑥 𝑗1 𝜋y 𝑗2 𝜋y 𝑤 ∗ = ∑ ∑ ∑ ∑ 𝜇𝑖1𝑖2𝑗1 𝑗2 𝑡 (sin ) (sin ) (sin ) (sin ) 𝑙 𝑙 b b 𝑖1 𝑖2 𝑗1 𝑗2 The above equations can be transformed into a cosine form as: 1 (𝑖2 − 𝑖1 ) 𝜋𝑥 (𝑗2 − 𝑗1 ) 𝜋𝑦 𝑤 = ∑ ∑ ∑ ∑ 𝑎𝑖1𝑖2 𝑗1 𝑗2 (𝑐𝑜𝑠 ( ) . 𝑐𝑜𝑠 ( ) 4 𝑙 𝑏 𝑖1 𝑖2 𝑗1 𝑗2 (𝑖2 − 𝑖1 ) 𝜋𝑥 (𝑗2 + 𝑗1 ) 𝜋𝑦 − 𝑐𝑜𝑠 ( ) . 𝑐𝑜𝑠 ( ) 𝑙 𝑏 (𝑖2 + 𝑖1 ) 𝜋𝑥 (𝑗2 − 𝑗1 ) 𝜋𝑦 − 𝑐𝑜𝑠 ( ) . 𝑐𝑜𝑠 ( ) 𝑙 𝑏 (𝑖2 + 𝑖1 ) 𝜋𝑥 (𝑗2 + 𝑗1 ) 𝜋𝑦 + 𝑐𝑜𝑠 ( ) . 𝑐𝑜𝑠 ( )) 𝑙 𝑏 32 1 (𝑖2 − 𝑖1 )𝜋𝑥 (𝑗2 − 𝑗1 )𝜋𝑦 𝑤 ∗ = ∑ ∑ ∑ ∑ 𝜇𝑖1𝑖2𝑗1 𝑗2 𝑡(𝑐𝑜𝑠 ( ) . 𝑐𝑜𝑠 ( ) 4 𝑙 𝑏 𝑖1 𝑖2 𝑗1 𝑗2 (𝑖2 − 𝑖1 )𝜋𝑥 (𝑗2 + 𝑗1 )𝜋𝑦 − 𝑐𝑜𝑠 ( ) . 𝑐𝑜𝑠 ( ) 𝑙 𝑏 (𝑖2 + 𝑖1 )𝜋𝑥 (𝑗2 − 𝑗1 )𝜋𝑦 (𝑖2 + 𝑖1 )𝜋𝑥 − 𝑐𝑜𝑠 ( ) . 𝑐𝑜𝑠 ( ) + 𝑐𝑜𝑠 ( ).𝐶 𝑙 𝑏 𝑙 (𝑗2 + 𝑗1 ) 𝜋𝑦 = 𝑐𝑜𝑠 ( )) 𝑏 One more transformation of the equations gives us the following equations for radial deformation and initial imperfections. 1 (𝑖2 − 𝑖1 )𝜋𝑥 (𝑗2 − 𝑗1 )𝜋𝑦 𝑤 = ∑ ∑ ∑ ∑ 𝑎𝑖1𝑖2 𝑗1 𝑗2 (𝑐𝑜𝑠 ( + ) 8 𝑙 𝑏 𝑖1 𝑖2 𝑗1 𝑗2 (𝑖2 − 𝑖1 )𝜋𝑥 (𝑗2 − 𝑗1 )𝜋𝑦 (𝑖2 − 𝑖1 )𝜋𝑥 (𝑗2 + 𝑗1 )𝜋𝑦 + 𝑐𝑜𝑠 ( − ) − 𝑐𝑜𝑠 ( + ) 𝑙 𝑏 𝑙 𝑏 (𝑖2 − 𝑖1 )𝜋𝑥 (𝑗2 + 𝑗1 ) 𝜋𝑦 (𝑖2 + 𝑖1 )𝜋𝑥 (𝑗2 − 𝑗1 )𝜋𝑦 − 𝑐𝑜𝑠 ( − ) − 𝑐𝑜𝑠 ( + )−𝐶 𝑙 𝑏 𝑙 𝑏 (𝑖2 + 𝑖1 ) 𝜋𝑥 (𝑗2 − 𝑗1 )𝜋𝑦 (𝑖2 + 𝑖1 ) 𝜋𝑥 (𝑗2 + 𝑗1 ) 𝜋𝑦 = 𝑐𝑜𝑠 ( − ) + 𝑐𝑜𝑠 ( + ) 𝑙 𝑏 𝑙 𝑏 (𝑖2 + 𝑖1 ) 𝜋𝑥 (𝑗2 + 𝑗1 ) 𝜋𝑦 + 𝑐𝑜𝑠 ( − )) 𝑙 𝑏 33 1 (𝑖2 − 𝑖1 ) 𝜋𝑥 (𝑗2 − 𝑗1 ) 𝜋𝑦 𝑤 ∗ = ∑ ∑ ∑ ∑ 𝜇𝑖1𝑖2𝑗1 𝑗2 𝑡(𝑐𝑜𝑠 ( + ) 8 𝑙 𝑏 𝑖1 𝑖2 𝑗1 𝑗2 (𝑖2 − 𝑖1 ) 𝜋𝑥 (𝑗2 − 𝑗1 ) 𝜋𝑦 (𝑖2 − 𝑖1 ) 𝜋𝑥 (𝑗2 + 𝑗1 ) 𝜋𝑦 + 𝑐𝑜𝑠 ( − ) − 𝑐𝑜𝑠 ( + ) 𝑙 𝑏 𝑙 𝑏 (𝑖2 − 𝑖1 ) 𝜋𝑥 (𝑗2 + 𝑗1 ) 𝜋𝑦 (𝑖2 + 𝑖1 ) 𝜋𝑥 (𝑗2 − 𝑗1 ) 𝜋𝑦 − 𝑐𝑜𝑠 ( − ) − 𝑐𝑜𝑠 ( + ) 𝑙 𝑏 𝑙 𝑏 (𝑖2 + 𝑖1 ) 𝜋𝑥 (𝑗2 − 𝑗1 )𝜋𝑦 (𝑖2 + 𝑖1 ) 𝜋𝑥 (𝑗2 + 𝑗1 ) 𝜋𝑦 − 𝑐𝑜𝑠 ( − ) + 𝑐𝑜𝑠 ( + ) 𝑙 𝑏 𝑙 𝑏 (𝑖2 + 𝑖1 ) 𝜋𝑥 (𝑗2 + 𝑗1 ) 𝜋𝑦 + 𝑐𝑜𝑠 ( − )) 𝑙 𝑏 In the literature, there are only a few research studies that considered clamped boundaries and tried to solve the fundamental shell equations using pre-defined radial deformations. Even those studies oversimplified the problem by considering the deformation equation as follows [34], [83]. 𝑚 𝑛 𝑖1 𝜋𝑥 2 𝑖2 𝜋y 2 (11) 𝑤 = ∑ ∑ 𝑎𝑖1 𝑖2 (sin ) (sin ) 𝑙 b 𝑖1 =1 𝑖2 =1 This equation implies that 𝑖1 and 𝑖2 and 𝑗1, and 𝑗2 are identical, respectively. If we reshape Equation (11) based on this assumption, w can be written as follow: 𝑚 𝑛 1 2𝑖2 𝜋𝑦 2𝑖1 𝜋𝑥 2𝑖1 𝜋𝑥 2𝑖2 𝜋𝑦 𝑤 = ∑ ∑ 𝑎𝑖1 𝑖2 (2 − 2𝑐𝑜𝑠 ( ) − 2𝑐𝑜𝑠 ( ) + 𝑐𝑜𝑠 ( + ) 8 𝑏 𝑙 𝑙 𝑏 𝑖1 =1 𝑖2 =1 2𝑖1 𝜋𝑥 2𝑖2 𝜋𝑦 + 𝑐𝑜𝑠 ( − )) 𝑙 𝑏 By considering the deformation as Equation (11), we lose the chance of having both negative and positive deformations in a panel in one deformation mode, which is not necessarily a correct deformed shape. The other problem with this assumption is that we have only one multiplying factor for any of the cosine expressions. For instance, if we consider 𝑐𝑜𝑠(2𝜋𝑥⁄𝑙 + 10𝜋𝑦⁄𝑏), it only appears in the expression of 𝑎15 . These two limitations make the response inaccurate, and we need 34 a large number of expressions to find a decent post-buckling response. However, we can assume the radial deformation in a way that covers the old assumption. The following expression for the radial displacement and initial imperfections can satisfy the boundary conditions. Λ 𝜋𝑥 𝑚𝑖 𝜋𝑥 𝜋y 𝑛𝑖 𝜋y 𝑤 = ∑ 𝑎𝑖 (sin ) (sin ) (sin ) (sin ) 𝑙 𝑙 b b 𝑖=1 (12) Λ ∗ 𝜋𝑥 𝑚𝑖 𝜋𝑥 𝜋y 𝑛𝑖 𝜋y 𝑤 = ∑ 𝜇𝑖 𝑡 (sin ) (sin ) (sin ) (sin ) 𝑙 𝑙 b b 𝑖=1 Here, x and y are the in-plane coordinates of the panel, 𝑚𝑖 , 𝑛𝑖 are the number of the waves in the x and y directions (respectively), 𝑎𝑖 is the unknown coefficient for the radial deformation, 𝜇𝑖 is the known coefficient for the initial imperfection, and t is the panel thickness. A larger number of modes 𝑚𝑘 , and 𝑛𝑘 considered in the solution is expected to provide a more accurate force- displacement response. However, increasing the number of modes in the deformation functions increases the computational cost of the solution. In Equation (12), we can write: 𝜋𝑥 𝑚𝑖 𝜋𝑥 1 (𝑚𝑖 − 1) 𝜋𝑥 (𝑚𝑖 + 1) 𝜋𝑥 (sin ) (sin ) = (𝑐𝑜𝑠 ( ) − 𝑐𝑜𝑠 ( )) 𝑙 𝑙 2 𝑙 𝑙 𝜋𝑦 𝑛𝑖 𝜋𝑦 1 (𝑛𝑖 − 1) 𝜋𝑦 𝑛𝑖 + 1) 𝜋𝑦 (sin ) (sin ) = (𝑐𝑜𝑠 ( ) − 𝑐𝑜𝑠 ( )) 𝑏 𝑏 2 𝑏 𝑏 Then we can write 𝑤 as follows: 35 Λ 1 (𝑚𝑖 − 1) 𝜋𝑥 (𝑛𝑖 − 1) 𝜋𝑦 (𝑚𝑖 − 1) 𝜋𝑥 (𝑛𝑖 + 1) 𝜋𝑦 𝑤 = ∑ 𝑎𝑖 (𝑐𝑜𝑠 ( ) . 𝑐𝑜𝑠 ( ) − 𝑐𝑜𝑠 ( ) . 𝑐𝑜𝑠 ( ) 4 𝑙 𝑏 𝑙 𝑏 𝑖=1 (𝑚𝑖 + 1) 𝜋𝑥 (𝑛𝑖 − 1) 𝜋𝑦 − 𝑐𝑜𝑠 ( ) . 𝑐𝑜𝑠 ( ) 𝑙 𝑏 (𝑚𝑖 + 1) 𝜋𝑥 (𝑛𝑖 + 1) 𝜋𝑦 + 𝑐𝑜𝑠 ( ) . 𝑐𝑜𝑠 ( )) 𝑙 𝑏 Again, considering 1 𝑐𝑜𝑠(𝛼). 𝑐𝑜𝑠(𝛽) = (𝑐𝑜𝑠(𝛼 + 𝛽) + 𝑐𝑜𝑠(𝛼 − 𝛽)) 2 the function for w ca be reconfigured as follows: Λ 1 (𝑚𝑖 − 1) 𝜋𝑥 (𝑛𝑖 − 1) 𝜋𝑦 (𝑚𝑖 − 1) 𝜋𝑥 (𝑛𝑖 − 1) 𝜋𝑦 𝑤 = ∑ 𝑎𝑖 (𝑐𝑜𝑠 ( + ) + 𝑐𝑜𝑠 ( − ) 8 𝑙 𝑏 𝑙 𝑏 𝑖=1 (𝑚𝑖 − 1) 𝜋𝑥 (𝑛𝑖 + 1) 𝜋𝑦 (𝑚𝑖 − 1) 𝜋𝑥 (𝑛𝑖 + 1) 𝜋𝑦 − 𝑐𝑜𝑠 ( + ) − 𝑐𝑜𝑠 ( − ) 𝑙 𝑏 𝑙 𝑏 (𝑚𝑖 + 1) 𝜋𝑥 (𝑛𝑖 − 1) 𝜋𝑦 (𝑚𝑖 + 1) 𝜋𝑥 (𝑛𝑖 − 1)𝜋𝑦 − 𝑐𝑜𝑠 ( + ) − 𝑐𝑜𝑠 ( − ) 𝑙 𝑏 𝑙 𝑏 (𝑚𝑖 + 1) 𝜋𝑥 (𝑛𝑖 + 1) 𝜋𝑦 (𝑚𝑖 + 1) 𝜋𝑥 (𝑛𝑖 + 1) 𝜋𝑦 + 𝑐𝑜𝑠 ( + ) + 𝑐𝑜𝑠 ( − )) 𝑙 𝑏 𝑙 𝑏 The same procedure is applied to 𝑤 ∗ 36 Λ 1 (𝑚𝑖 − 1) 𝜋𝑥 (𝑛𝑖 − 1) 𝜋𝑦 (𝑚𝑖 − 1) 𝜋𝑥 (𝑛𝑖 − 1) 𝜋𝑦 𝑤 ∗ = ∑ 𝜇𝑖 𝑡(𝑐𝑜𝑠 ( + ) + 𝑐𝑜𝑠 ( − ) 8 𝑙 𝑏 𝑙 𝑏 𝑖=1 (𝑚𝑖 − 1) 𝜋𝑥 (𝑛𝑖 + 1) 𝜋𝑦 (𝑚𝑖 − 1) 𝜋𝑥 (𝑛𝑖 + 1) 𝜋𝑦 − 𝑐𝑜𝑠 ( + ) − 𝑐𝑜𝑠 ( − ) 𝑙 𝑏 𝑙 𝑏 (𝑚𝑖 + 1) 𝜋𝑥 (𝑛𝑖 − 1) 𝜋𝑦 (𝑚𝑖 + 1) 𝜋𝑥 (𝑛𝑖 − 1)𝜋𝑦 − 𝑐𝑜𝑠 ( + ) − 𝑐𝑜𝑠 ( − ) 𝑙 𝑏 𝑙 𝑏 (𝑚𝑖 + 1) 𝜋𝑥 (𝑛𝑖 + 1) 𝜋𝑦 (𝑚𝑖 + 1) 𝜋𝑥 (𝑛𝑖 + 1) 𝜋𝑦 + 𝑐𝑜𝑠 ( + ) + 𝑐𝑜𝑠 ( − )) 𝑙 𝑏 𝑙 𝑏 Having 𝑤 and 𝑤 ∗ in the above forms, their second derivatives follow the same form with different coefficients. Λ 𝜋2 (𝑚𝑖 − 1) 𝜋𝑥 (𝑛𝑖 − 1) 𝜋𝑦 (𝑚𝑖 − 1) 𝜋𝑥 (𝑛𝑖 − 1) 𝜋𝑦 𝑤,𝑥𝑥 = ∑ 2 𝑎𝑖 (−(𝑚𝑖 − 1)2 (𝑐𝑜𝑠 ( + ) + 𝑐𝑜𝑠 ( − ) 8𝑙 𝑙 𝑏 𝑙 𝑏 𝑖=1 (𝑚𝑖 − 1) 𝜋𝑥 (𝑛𝑖 + 1) 𝜋𝑦 (𝑚𝑖 − 1) 𝜋𝑥 (𝑛𝑖 + 1) 𝜋𝑦 − 𝑐𝑜𝑠 ( + ) − 𝑐𝑜𝑠 ( − )) 𝑙 𝑏 𝑙 𝑏 (𝑚𝑖 + 1) 𝜋𝑥 (𝑛𝑖 − 1)𝜋𝑦 (𝑚𝑖 + 1) 𝜋𝑥 (𝑛𝑖 − 1) 𝜋𝑦 + (𝑚𝑖 + 1)2 (𝑐𝑜𝑠 ( + ) + 𝑐𝑜𝑠 ( − ) 𝑙 𝑏 𝑙 𝑏 (𝑚𝑖 + 1) 𝜋𝑥 (𝑛𝑖 + 1) 𝜋𝑦 (𝑚𝑖 + 1) 𝜋𝑥 (𝑛𝑖 + 1) 𝜋𝑦 − 𝑐𝑜𝑠 ( + ) − 𝑐𝑜𝑠 ( − ))) 𝑙 𝑏 𝑙 𝑏 Λ 𝜋2 (𝑚𝑖 − 1) 𝜋𝑥 (𝑛𝑖 − 1) 𝜋𝑦 (𝑚𝑖 − 1) 𝜋𝑥 (𝑛𝑖 − 1)𝜋𝑦 𝑤,𝑦𝑦 = ∑ 2 𝑎𝑖 (−(𝑛𝑖 − 1)2 (𝑐𝑜𝑠 ( + ) + 𝑐𝑜𝑠 ( − ) 8𝑏 𝑙 𝑏 𝑙 𝑏 𝑖=1 (𝑚𝑖 + 1) 𝜋𝑥 (𝑛𝑖 − 1)𝜋𝑦 (𝑚𝑖 + 1)𝜋𝑥 (𝑛𝑖 − 1) 𝜋𝑦 − 𝑐𝑜𝑠 ( + ) − 𝑐𝑜𝑠 ( − )) 𝑙 𝑏 𝑙 𝑏 (𝑚𝑖 − 1)𝜋𝑥 (𝑛𝑖 + 1)𝜋𝑦 (𝑚𝑖 − 1)𝜋𝑥 (𝑛𝑖 + 1)𝜋𝑦 + (𝑛𝑖 + 1)2 (𝑐𝑜𝑠 ( + ) + 𝑐𝑜𝑠 ( − ) 𝑙 𝑏 𝑙 𝑏 (𝑚𝑖 + 1)𝜋𝑥 (𝑛𝑖 + 1)𝜋𝑦 (𝑚𝑖 + 1)𝜋𝑥 (𝑛𝑖 + 1)𝜋𝑦 − 𝑐𝑜𝑠 ( + ) − 𝑐𝑜𝑠 ( − ))) 𝑙 𝑏 𝑙 𝑏 37 Λ 𝜋2 (𝑚𝑖 − 1)𝜋𝑥 (𝑛𝑖 − 1)𝜋𝑦 𝑤,𝑥𝑦 =∑ 𝑎𝑖 (−(𝑚𝑖 − 1)(𝑛𝑖 − 1)(𝑐𝑜𝑠 ( + ) 8𝑙𝑏 𝑙 𝑏 𝑖=1 (𝑚𝑖 − 1)𝜋𝑥 (𝑛𝑖 − 1)𝜋𝑦 − 𝑐𝑜𝑠 ( − )) + (𝑚𝑖 − 1)(𝑛𝑖 𝑙 𝑏 (𝑚𝑖 − 1)𝜋𝑥 (𝑛𝑖 + 1)𝜋𝑦 (𝑚𝑖 − 1)𝜋𝑥 (𝑛𝑖 + 1)𝜋𝑦 + 1)(𝑐𝑜𝑠 ( + ) − 𝑐𝑜𝑠 ( − )) 𝑙 𝑏 𝑙 𝑏 (𝑚𝑖 + 1)𝜋𝑥 (𝑛𝑖 − 1)𝜋𝑦 + (𝑚𝑖 + 1)(𝑛𝑖 − 1)(𝑐𝑜𝑠 ( + ) 𝑙 𝑏 (𝑚𝑖 + 1)𝜋𝑥 (𝑛𝑖 − 1)𝜋𝑦 − 𝑐𝑜𝑠 ( − )) − (𝑚𝑖 + 1)(𝑛𝑖 𝑙 𝑏 (𝑚𝑖 + 1) 𝜋𝑥 (𝑛𝑖 + 1) 𝜋𝑦 (𝑚𝑖 + 1) 𝜋𝑥 (𝑛𝑖 + 1) 𝜋𝑦 + 1)(𝑐𝑜𝑠 ( + ) − 𝑐𝑜𝑠 ( − ))) 𝑙 𝑏 𝑙 𝑏 As we put these expressions in Equation (5), it is obvious that 𝛻 4 𝑓 becomes the product of several (𝑐𝑜𝑠((𝑚𝑖1 ± 1)πx/l + (𝑛𝑖1 ± 1)𝜋𝑦/𝑏) × 𝑐𝑜𝑠((𝑚𝑖2 ± 1)πx/l + (𝑛𝑖2 ± 1) πy⁄𝑏)) which can be written as follow: 2mΛ +2 2𝑛Λ +2 2 𝑖𝜋𝑥 𝑗(−1)𝑘 𝜋𝑦 𝛻4𝑓 = ∑ ∑ ∑ 𝑓𝑖𝑗𝑘 𝑐𝑜𝑠 ( + ) (13) 𝑙 𝑏 𝑖=0 𝑗=0 𝑘=1 where 𝑚Λ and 𝑛Λ are the maximum number of contributed modes in Equation (12) for the x and y directions, and 𝑓𝑖𝑗𝑘 is the multiplier for the cosine function in 4f. A computer code developed in the software Mathematica to solve Equation (13). The result is the function F, which has similar cosine expressions as 4f but different multipliers to the cosine terms (Fijk). 2mΛ +2 2𝑛Λ +2 2 𝑖𝜋𝑥 𝑗(−1)𝑘 𝜋𝑦 𝐹(𝑥, 𝑦) = ∑ ∑ ∑ F𝑖𝑗𝑘 𝑐𝑜𝑠 ( + ) (14) 𝑙 𝑏 𝑖=0 𝑗=0 𝑘=1 38 For any given set of {{𝑚1 , 𝑛1 }, … , {𝑚Λ , 𝑛Λ }} all the multipliers 𝐹𝑖𝑗𝑘 are known. However, the Airy stress function can also contain functions (y) and (x) with 4(y) = 4(x) = 0. Then, it should be expressed as follows: 𝑓 = 𝜆(𝑦) + θ(x) + 𝐹(𝑥, 𝑦) (15) Since the side edges of the panel are movable in the direction of the applied axial force (see Figure 2-7(a)), the integration of Nx over the width of the panel should be equal to the applied force (P). Thus, 𝑏 1 𝑃 −𝑃 = ∫ 𝑓,𝑦𝑦 𝑑𝑦 ⟹ 𝜆(𝑦) = − ( ) ( ) 𝑦 2 (16) 0 2 𝑏 Also, the moveable edges on the sides of the panel only can move parallel to the x-direction, and as a result the change in distance between these edges should be zero during the loading process. This condition is fulfilled in an average sense as [84]: 𝑙 𝑏 𝑉 = ∫ ∫ 𝑣,𝑦 𝑑𝑦𝑑𝑥 = 0 0 0 (17) 𝑁𝑦 − 𝜈𝑁𝑥 𝑤 1 2 𝑣,𝑦 = + − 𝑤,𝑦 − 𝑤,𝑦 𝑤,𝑦∗ 𝐸𝑡 𝑅 2 Thus (x) can be expressed as: 1 𝑃 2 1 𝑙 𝑏 𝑤 1 𝜃(𝑥) = 𝜈 ( ) ( ) 𝑥 − ∫ ∫ (𝑓,𝑥𝑥 − 𝜈𝑓,𝑦𝑦 ) + 𝐸. 𝑡 ( − 𝑤,𝑦2 − 𝑤,𝑦 𝑤,𝑦∗ ) 𝑑𝑦𝑑𝑥 (18) 2 𝑏 𝑙𝑏 0 0 𝑅 2 Having the Airy stress function related to the radial deformation allows solving Equation (5) numerically for any axial deformation. The axial deformation itself is calculated from Equation (16). 𝑙 𝑑 = ∫ 𝑢,𝑥 𝑑𝑥 (19) 0 39 𝑁𝑥 − 𝜈𝑁𝑦 1 2 𝑢,𝑥 = − 𝑤,𝑥 − 𝑤,𝑥 𝑤,𝑥∗ 𝐸𝑡 2 By replacing f and w in Equation (5) and using Galerkin’s approximation method for each unknown in Equation (12), an equation is produced. It means that any 𝑚𝑖 , between 1 and 𝑚Λ and any 𝑛𝑖 between 1 and 𝑛Λ : 𝑙 𝑏 𝜋𝑥 𝑚𝑖 𝜋𝑥 𝜋y 𝑛𝑖 𝜋y 𝑓,𝑥𝑥 ∫ ∫ (sin ) (sin ) (sin ) (sin ) (𝐷𝛻 4 𝑤 − 𝑥=0 𝑦=0 𝑙 𝑙 b b 𝑅 ∗ ∗ (20) − [𝑓,𝑦𝑦 (𝑤,𝑥𝑥 + 𝑤,𝑥𝑥 ) − 2𝑓,𝑥𝑦 (𝑤,𝑥𝑦 + 𝑤,𝑥𝑦 ) ∗ + 𝑓,𝑥𝑥 (𝑤,𝑦𝑦 + 𝑤,𝑦𝑦 )])𝑑𝑥𝑑𝑦 = 0 The number of equations is equal to Λ, which is one equation less than the number of unknowns in Equation (12). The unknowns are 𝑎1 to 𝑎Λ and 𝑃. Axial deformation 𝑑 can be calculated from having these unknown variables. Thus, we have Λ + 1 unknowns and Λ equations. The last equation comes from assuming that all the unknown variables are function of the curve length (s). 𝑎1 = 𝑎1 (𝑠), … , 𝑎Λ = 𝑎Λ (𝑠), 𝑃 = 𝑃(𝑠) As a result, the axial end shortening d is a function of s. Then, based on general continuation theory, these variables should be continuous with respect to s, which means: Λ 2 𝛿𝑎𝑖 2 𝑙(1 − 𝜈 2 ) 𝛿𝑝 𝛿𝑑 2 ∑( ) + ( ) +( ) =1 (21) 𝛿𝑠 𝐸𝑡𝑏 𝛿𝑠 𝛿𝑠 𝑖=1 In Equation (21), 𝛿𝑝⁄𝛿𝑠 is multiplied by the factor of 𝑙(1 − 𝜈 2 )⁄𝐸𝑡𝑏 to scale it to a reasonable range close to other variables. Without this scaling factor the solution process will be very slow or will not converge. These equations are solved using the pseudo arclength method [85]. The method is simply Newton’s method with consideration of Equation (21). Since Equation (21) requires a solution for specific 𝛿𝑠 at each increment, this solution is called pseudo arclength method. 40 Figure 2-8: Steps of semi-analytical solution for the post-buckling response of a cylindrical shell under compression Now we know how to solve equations for any given combination of mode shapes, but the question is which modes we should consider in the calculation. We previously used Equation (10) to introduce radial deformation and changed 𝑖1 , 𝑖2 , 𝑗1 , and 𝑗2 from 1 to 4 [86]. This left us with 100 equations to solve, which still did not match the post-buckling response of panels and cylinders captured from FE analysis very well. In addition, the resulting response from such an assumption did not solve the imperfection sensitivity of the cylindrical panels. However, by taking a deeper look into the stress distribution in the NSD segments, we can see that only certain modes can be activated during axial loading. 41 (a) (b) Figure 2-9: (a) Schematics of NSD cylindrical segments, (b) schematic distribution of stress prior to buckling event Figure 2-9 shows a representative stress distribution in the NSD panels prior to the buckling event. It shows significant stress concentrations in the middle part of the panel, which is the result of a shorter thin part in the mid-width of the NSD segment. The big difference between stresses in the mid-section and edges ensures that buckling modes with the maximum radial deformation occur at the center of the panel. Therefore, it is reasonable to consider only the modes with the potential of maximum deformation at the panel center. All the modes with 𝑚𝑖 and 𝑛𝑖 as odd numbers in Equation (9) fulfill this requirement. The semi-analytical model was implemented in three programs in Mathematica Wolfram software. The first program finds the Airy stress function by assuming the form of radial deformation. The second program develops nonlinear equations from Equation (6) by using Galerkin’s method. The third program solves these equations and finds the axial force and deformation, radial deformation in the center of the panel, and external energy of the panel. 42 In the rest of this study, the semi-analytical model uses Equation (12) for the radial deformation of the panel with only 49 deformations modes (𝑚𝑖 & 𝑛𝑖 ∈ {1,3,5,7,9,11,13}). Also, since the compression at the center of the panel is very high, the first buckling mode looks very probable, and thus a large initial imperfection coefficient (0.001h) is applied to the first mode, and other modes are kept to 10% of that coefficient. This way, the first mode dominates the response over the others, which based on Figure 2-9 looks reasonable. To continue, we find out about the validity of such an assumption. To better show how the current semi-analytical model improves the prediction of post-buckling response of NSD segments, the model’s results are compared with the responses obtained from FE analyses, experiments, and a semi-analytical model previously developed by the author [86]. Using the semi-analytical model with the pseudo-arc length method (Equation (21)) provides the path for the panel’s force-displacement response. If buckling occurs in the panel, the force- displacement curve shows snap-back behavior. It means that after buckling both the force and displacement decrease to take the force-deformation response to a new equilibrium position. FEA and experiments cannot capture this behavior since it is only a mathematical explanation. Rather, these methods show a force plateau or a drop in force at a constant displacement for load control and displacement control, respectively. In this study, the experiments and finite element analyses were conducted in displacement control. Using displacement-controlled loading helps to capture the load-deformation of the system even when the slope of the force-displacement response of the cylindrical panel is zero or negative after buckling in the panel. 43 Figure 2-10: Equilibrium path in snap-back behavior compared to force-displacement from force control and deformation control loading 2.3.1.2. Model verification for post-buckling response Four finite element models with different geometries were prepared using the software ABAQUS [87] to verify the semi-analytical model. The geometry of the panels is shown in Table 2-1. Four node SR4 shell elements represent the cylindrical panels, and a linear elastic material was assigned with an elastic modulus of 1300 N/mm2 and Poisson’s ratio of 0.2. The top and bottom edges of the panels are clamped, while the side boundaries are fixed for radial and circumferential deformations and outward rotations but deform in the axial direction. Loading was in the form of a uniform axial displacement, equal to 1.5% of the panel’s length, applied to one end of the panel. To eliminate the effect of vibrations from the buckling events, a dynamic implicit solver with the quasi-static application was used with a maximum time step of 0.01 s. Results of finite element analyses were compared to the semi-analytical model using the general equation for 𝑤 and 𝑤 ∗ (see Equation (10)) [86]. To have an acceptable response from a panel, numerous cosine modes are needed. However, as described in Section 0, the model for 𝑤 44 and 𝑤 ∗ were improved as given in Equation (11). For the NSD cylinders, where stresses concentrate in the center of the panel, 𝑚𝑖 and 𝑛𝑖 are odd numbers (in this case 1,3,5,7,9,11,13). However, for the panels without stiffened areas concentrated in the center, radial deformations from Equation (11) will not necessarily match with correct answer. Thus, comparisons between FE analyses and the semi-analytical model are only done for the radial deformation from Equation (10). To validate the efficiency of Equation (11) for predicting the behavior of an NSD segment, results from the semi-analytical model to FEA solutions and experimental results are provided in Chapter 3. In the semi-analytical model, the initial imperfections follow Equation (10), while in the FE model two options are available to apply the initial imperfection: 1) conducting a linear eigenvalue analysis and importing modes of interest to the post-buckling model, or 2) considering the imperfection at the node’s coordinates. In this study, the initial imperfections were applied at the coordinates of the nodes based on the first buckling mode with 𝜇1111 = 0.05 and all other 𝜇𝑖1𝑖2𝑗1𝑗2 in Equation (10) equal to zero. The panel side edges were clamped but with the ability of axial motion along y = 0, y = b, and x = l (see Figure 2-7). A uniform axial compressive displacement equal to 0.015h was applied to the top edge of the panel. For each panel in Table 2-1, four semi-analytical models were developed for ( =  = 1) to ( =  = 4) in Equation (12). The material properties and initial imperfections in the semi- analytical models are identical to those in the FE models. 45 Table 2-1: Dimensions of cylindrical panels for verification of the analytical model Length Width Radius Thickne Cases (l) (b) (R) ss (t) mm mm mm mm (a) 80 80 50 0.5 (b) 80 40 50 0.5 (c) 40 60 50 0.5 (d) 60 40 50 0.5 Figure 2-11 shows that the analytical model’s accuracy improves as the number of modes contributing to the radial deformation equation increases, especially for larger cylindrical panels. That is, more modes should be used in the analytical model when evaluating bigger panels in order to obtain an accurate prediction of the post-buckling response. For all cases, the analytical model with  =  = 4 accurately predicts the displacement corresponding to the first buckling event. However, the panel’s post-buckling stiffness predicted by the analytical model is higher than the response calculated by FEA. While increasing the number of involved modes in the radial deformation ( and ) leads to improved accuracy, it also increases analysis time since the number of nonlinear equations escalates significantly. For example, by increasing  and  from 4 to 5, the number of nonlinear equations to solve rises from 100 to 225, which significantly increases computational time. 46 Figure 2-11: Comparison of force-deformation curves from FE and analytical models for four cylindrical panels with all-clamped boundary conditions 47 Chapter 3. NSD buckling segment discussions Chapter 2 presented the concept of replacing the parallel segments of the NSD cylinder with equivalent panels and then matching the stiffness of the equivalent system with the original segment by adding linear springs at the top and bottom of the equivalent panel. The equivalent panel was previously defined as a cylindrical panel with a uniform thickness that had the same width as the original segment but a different height. However, this simplification neglects the effect of other dimensions of the stiffened area in the NSD segment. In this Chapter, the effect of those dimensions is discussed, and the equivalent panel dimensions are modified to consider the width and the depth of the stiffened zones in a simplified way. Then, relations between the force and deformation in the segment and the equivalent panel are developed for later use in designing NSD cylinders with a desired post-buckling response. Finally, the validity of the assumptions in this section is examined against FE analyses for several given NSD segments. 3.1. Thickened area and equivalent panel The stiffened areas in NSD cylinders in previous works of Hu and Burgueño and Guo, Liu and Burgueño are assumed as two rectangular areas on top of each other[2], [70], [72], [86], [88]. This approach was taken to 1) concentrate the axial stresses in the center of the thin zone, and 2) avoid modeling and production complications. For the first reason, it is sound to have the shortest thin zone in the mid area of the segment. When the thin area is shorter, it attracts more stress, and therefore buckling more likely starts there. However, the second reason (simplicity) is responsible for the rectangular shape of the stiffened areas since a triangular shape or curve with a peak at the center also fulfills the goal of concentrating axial stresses in the mid buckling zone. 48 Here, for the same very valid reasons, the shape of the thickened areas is kept the same as in previous works. In the proposed framework, an equivalent panel has replaced the thickened area. However, the effect on other dimensions of the thickened area, such as the width of the thickened area in comparison to the width of the segment or the height of the thickened area in the bottom and top of the cylinder, were not considered. For instance, the framework takes the same equivalent length for all the thick/stiff segments in Figure 3-1 since there is no impact in the definition of the equivalent length from the width of the thickened area [86]. Having the same equivalent length for equivalent panels means that these segments share the same buckling force and drop in force after buckling. The end shortening corresponding to the buckling event and drop in energy of the panel is not independent of the linear springs at the top and bottom of the panel, and thus they can be different for the segments with the same equivalent panel length. (a) (b) (c) Figure 3-1: Cylindrical segments from NSD cylinder in which the stiff area gets bigger from (a) to (c) A finite element model for cylindrical segments like that shown in Figure 3-1 was developed using the software ABAQUS to show how the cylindrical segments with the same equivalent length behave differently under axial compression. The length, width, radius, thickness of thin and 49 thick parts, and the maximum length of the thick area are the same for both segments. The only difference is the width of the thick area. Dimensions of the cylindrical segments are presented in Table 3-1. The material is acrylonitrile butadiene styrene (ABS) with an elastic modulus (𝐸) equal to 1300 N/mm2 and Poison’s (𝜈) ratio of 0.3. The equivalent length of the panel (el) using Equation (2) in the framework for all the panels is calculated as follows: 3 𝑡𝑏 3 0.5 2 𝑙 ′1𝑖 = ( )2 𝑙1𝑖 . = ( ) × 60 = 11.54 𝑚𝑚 → 𝑒𝑙 = 11.54 + 40 = 51.54 𝑚𝑚 𝑡1𝑖 1.5 Table 3-1: Geometry of NSD cylinder segments (all dimensions are in mm) R h h1 b b1(a) b1(b) b1(c) l1(a) l1(b) l1(c) tb t1 50 100 40 65 24 32 40 25 20 15 0.5 1.5 The comparison of results in Figure 3-1 shows that despite the proposed model considering all these segments as equal (with the exception of deformation corresponding to buckling event and energy drop), the dimensions of the thick area has a significant role in the desired design factors from buckling of the segment such as energy drop, force drop, and radial deformation of the segment right after buckling. 50 (a) (b) (c) (d) (e) Figure 3-2: Comparing the response of 3 cylindrical segments with the stiffened pattern in Figure 3-1and dimensions in the Table 3-1 (a)the difference in radial deformation of the segment, (b) the difference in force-displacement response, (c) the difference in energy dissipation in the buckling event of each segment, (d) the variation of axial deformation corresponding to the buckling event, and (e) the difference in force drops from buckling in each of the segments 51 In the following, the response of the equivalent cylindrical panel using the semi-analytical model is developed. Comparing the results of segments (a), (b), and (c) with the response from the semi-analytical model for the equivalent panel shows that Case (c) is very close to the assumption of the equivalent panel. The buckling load from the semi-analytical model for Case (c) is 222 N, while the FEA for the same case gives a force of 214.4 N right before the buckling event. The semi-analytical model shows the biggest difference for Case (a). Also, the radial deformation of Cases (a), (b), and (c) in Figure 3-2 shows that Case (c) has a wider deformed area than Case (b); and Case (b) has a bigger deformed area than Case (a). The differences in the cases with equal equivalent length raise the idea of considering an equivalent width for the buckling panel as well as equivalent length. This way, the effect of size of the thickened region on the segment will show up in the equivalent panel. This panel response is responsible for all the nonlinear behavior in the segment. However, we should consider the additional linear springs surrounding this buckling panel to match the linear stiffness of the segment, like what is shown in Figure 3-3. To simplify the system, it is assumed that: 1) the top and bottom of the equivalent panel deform uniformly, and 2) the remaining width of the segment is modeled as a linear spring. Here, the stiffness of the linear springs is calculated for any given dimensions of the cylindrical segment and any possible equivalent panel. These calculations are necessary for the design of NSD cylinders for a given post-buckling response. 52 Figure 3-3: Replacing the segment with equivalent panel plus linear springs The stiffness of the linear springs on top and bottom of the equivalent panel are calculated to provide the same axial stress in the center of the panel. Thus, 2 𝑙′ ℎ − ℎ1 ℎ1 + ′ = ′ + ′ 𝑘1 𝐸𝑏 𝑡𝑏 𝐸𝑏 𝑡1 𝐸𝑏 𝑡𝑏 (22) 2Eb′𝑡1 𝑡𝑏 𝑘1 = (ℎ1 − l′)𝑡1 + (h − ℎ1 )𝑡𝑏 The width of the panel is calculated in a way that the total axial force in the section of the equivalent panel becomes identical to the axial force in the section of the segment. It is worth noting that this assumption is a simplification of the situation in the segment. Of course, this assumption neglects several complexities in the behavior of the segment. To define the width of equivalence, the following assumptions were made: 53 • The axial stress under the central 𝑏1 width of the segment is constant and equal to 𝜎𝑐 (stress in the center). • The axial stress in the remaining 𝑏 − 𝑏1 width of the segment is constant and equal to 𝜎𝑠 (stress on the side). • The slope of change in the stiffness is 1⁄𝑙1. • The width of the equivalent panel is calculated in a way that the axial force in the panel and segment match each other. • If the width of the segment is larger than 𝑏1 + 2𝑙1 , then the panel's width is calculated for the segment with 𝑏1 + 2𝑙1 width, and the remaining part of the segment is treated as segments with a linear response (segments D) For the case of 𝑏 > 𝑏1 + 2𝑙1, the stiffness of the remaining part is calculated as follows: 1 ℎ1 + 2𝑙1 ℎ − (ℎ1 + 2𝑙1 ) 2ℎ1 𝑡1 + 4𝑙1 𝑡1 + ℎ𝑡𝑏 − ℎ1 𝑡𝑏 − 2𝑙1 𝑡𝑏 = ′ + = 𝑘2 𝐸(𝑏 − 𝑏 )𝑡𝑏 𝐸(𝑏 − 𝑏 ′ )𝑡1 𝐸(𝑏 − 𝑏 ′ )𝑡1 𝑡𝑏 2 2 𝐸(𝑏 − 𝑏 ′ )𝑡1 𝑡𝑏 𝑘2 = (23) 2ℎ1 𝑡1 + 4𝑙1 𝑡1 + 2ℎ𝑡𝑏 − 2ℎ1 𝑡𝑏 − 4𝑙1 𝑡𝑏 The stress at the center of the panel is calculated as follows: 𝑑𝑓 𝛿 𝜎𝑐 = = 𝑑𝑏𝑡𝑏 𝑑𝑏 𝑡𝑏 + 2 𝑑𝑏𝑡𝑏 (24) 𝑘1𝑎 𝑘1𝑏 The stress at 𝑏1 ⁄2 + 𝑙1 from the center of the panel is calculated as: 𝑓 𝛿 𝜎𝑠 = = (𝑏 − 𝑏1 )𝑡𝑏 𝑑𝑏 𝑡𝑏 2𝑑𝑏 𝑡𝑏 (25) + 𝑘1𝑐 𝑘1𝑑 The stiffness of springs in Equations (24) and (25) are expressed as follows: 54 𝐸 𝑑𝑏 𝑡𝑏 2𝐸 𝑑𝑏 𝑡1 𝐸 𝑑𝑏 𝑡𝑏 𝑘1𝑎 = , 𝑘1𝑏 = , 𝑘1𝑐 = , ℎ1 ℎ − ℎ1 (ℎ1 + 2𝑙1 ) (26) 2𝐸 𝑑𝑏 𝑡1 𝑘1𝑑 = ℎ − ℎ1 − 2𝑙1 Then: 𝑑𝑓 𝛿 𝛿 𝐸𝑡1 𝛿 𝜎𝑐 = = = = 𝑑𝑏𝑡𝑏 𝑑𝑏 𝑡𝑏 + 2 𝑑𝑏𝑡𝑏 ℎ1 (ℎ − ℎ1 )𝑡𝑏 ℎ1 𝑡1 + (ℎ − ℎ1 )𝑡𝑏 (27) 𝑘1𝑎 𝑘1𝑏 𝐸 + 𝐸 𝑡1 𝑑𝑓 𝛿 𝛿 𝜎𝑠 = = = 𝑑𝑏 𝑡𝑏 𝑑𝑏 𝑡𝑏 + 2𝑑𝑏 𝑡𝑏 ℎ1 + 2𝑙1 (ℎ − ℎ1 − 2𝑙1 )𝑡𝑏 𝑘1𝑐 𝑘1𝑑 + 𝐸 𝐸𝑡1 (28) 𝐸𝑡1 𝛿 = (ℎ1 + 2𝑙1 )𝑡1 + (ℎ − ℎ1 − 2𝑙1 )𝑡𝑏 ℎ1 𝑡1 + (ℎ − ℎ1 )𝑡𝑏 𝜎𝑠 = 𝜎𝑐 (29) (ℎ1 + 2𝑙1 )𝑡1 + (ℎ − ℎ1 − 2𝑙1 )𝑡𝑏 For the case when 𝑏 ≤ 2𝑙1 + 𝑏1 , 𝑏′ can be calculated as follow: 𝑏 ′ 𝑡𝑏 𝜎𝑐 = 𝑏1 𝑡𝑏 𝜎𝑐 + (𝑏 − 𝑏1 )𝑡𝑏 𝜎𝑠 (𝜎𝑠 + 𝜎𝑐 ) 𝑏 ′ 𝑡𝑏 𝜎𝑐 = 𝑏1 𝑡𝑏 𝜎𝑐 + (𝑏 − 𝑏1 )𝑡𝑏 2 𝑏 + 𝑏1 (𝑏 − 𝑏1 ) ℎ1 𝑡1 + (ℎ − ℎ1 )𝑡𝑏 𝑏′ = + 2 2 (ℎ1 + 2𝑙1 )𝑡1 + (ℎ − ℎ1 − 2𝑙1 )𝑡𝑏 ℎ1 𝑡1 + (ℎ − ℎ1 )𝑡𝑏 𝑏 ′ = 𝑏1 + (𝑏 − 𝑏1 ) (30) (ℎ1 + 2𝑙1 )𝑡1 + (ℎ − ℎ1 − 2𝑙1 )𝑡𝑏 55 3.2. Segment radial deformation after the buckling event In Section 3.1 the buckling segment has been simplified to three springs in series, two linear at top and bottom with a stiffness of 𝑘1 , and one nonlinear spring in the middle, which has the stiffness of 𝑘𝑠 prior to buckling and 𝑘𝑠 ’ after buckling. It is assumed that a deformation 𝛿 is applied at the end of the system and that it increases to the point that the middle spring “buckles.” At that point, the system loses some force since the middle spring cannot resist more force. If there were no other springs in the system, the force drop would happen at the exact same end shortening as for the equivalent panel. However, having the spring 𝑘1 makes the problem more complex. Now, the system of springs is not in balance as the force in the 𝑘1 spring dropped and then the deformation must change to produce the corresponding balancing force. The spring 𝑘1 must then extend, and considering the constant total deformation, the spring in the middle should compress more. Only then will both springs be able to balance the force 𝑃𝑏 . Figure 3-4 illustrates this process in detail. Total end shortening and the portion of the deformation in parallel springs are shown in Equation (31): 𝛿𝑡 = 𝛿𝑠 ′ + 𝛿′1 𝛿𝑠′ = 𝛿𝑡 − 𝑃0 ⁄𝑘1 (31) If 𝛿𝑠 ’ < 𝛿0 then 𝛿𝑠 = 𝛿𝑠 ’ When the middle spring (equivalent panel) buckles, then: 𝛿𝑡 − 𝛿𝑠 (𝛿𝑠 − 𝛿0 )𝑘𝑠 ’ + (𝑃0 − Δ𝑃0 ) = 𝑘1 2 𝑘1 𝛿𝑡 + 𝑘𝑠′ 𝛿0 − (𝑃0 − Δ𝑃0 ) 𝛿𝑠 = 2 (33) 𝑘1 ′ 2 + 𝑘𝑠 𝛿𝑡 − 𝛿𝑠 𝛿1 = 2 56 Also, the force decline in the buckling event instead of Δ𝑃0 should be calculate as following. Δ𝑃1 = 𝑃0 − 𝑃𝑏 = Δ𝑃0 + (𝛿𝑠 − 𝛿0 )𝑘𝑠 ’ (32) Figure 3-4: The deformation of the system of the springs in series when the middle spring buckles. (a) right before buckling, the system is in balance and force in the springs are equal, (b) right after buckling, the total deformation does not change, force in the buckled panel and linear springs in series are not in balance, and (c) the axial deformation in linear springs and buckled panel are adjusted to keep the force in balance 57 3.3. Segment energy drop from buckling The force-deformation response of the equivalent panel shows a drop in the internal energy when it snaps back. The extra applied work transforms into kinetic energy and causes rapid deformation in the panel [89],[90]. The amount of energy drop in the equivalent panel is calculated by multiplying the difference in the deformation by the average of the force in each step of the incremental numerical solution, which is equal to the area under the snapback part of the force- displacement curve (Figure 3-5). Figure 3-5: The work lost in the snapback behavior The solid red area in Figure 3-5 schematically shows the amount of work lost during the snapback response. However, this is not the only energy loss in the segment. Since the force in the segment drops suddenly due to buckling, the linear springs at the top and bottom of the panel lose some energy while the equivalent panel is in stable mode and use some of this energy for more deformation. The total energy fall in the segment is calculated as follows: Δ 𝐸 = Δ 𝐸𝑠 + 2Δ 𝐸1 − Δ 𝐸𝑠 ’ (34) 58 where Δ 𝐸𝑠 is the dissipated energy in the panel, and it is calculated using the presented semi- analytical model. Δ 𝐸1 is the dissipated energy in the top and bottom linear springs. and Δ 𝐸𝑠 ’ represents the work done by the equivalent panel after buckling to keep the balance of force and displacement. (See Figure 3-6) Figure 3-6: The energy dissipated, and the work done in the NSD segment in the moment of buckling The values of Δ 𝐸1 and Δ 𝐸𝑠 ’are calculated as follows: 1 𝛿𝑡 − 𝛿𝑠 2 Δ 𝐸1 = ( ) 𝑘1 ( ) 2 2 (35) 1 Δ 𝐸𝑠 ′ = ( ) 𝑘𝑠′ (𝛿𝑠 − 𝛿0 )2 2 where 𝛿0 , and 𝑘’𝑠 are found from the presented semi-analytical model, and 𝛿𝑠 is calculated from Equation (33). 59 3.4. Verification of the assumptions To examine the validity of the proposed modeling approach, the equivalent widths for panels (a), (b), and (c) in Figure 3-1 were determined. To do so, five analytical models for the same equivalent length but different panel widths ( ’= 40 mm, ’=45 mm, ’=50 mm, ’=55 mm, and ’=60 mm) were developed. Then, for cases (a), (b), and (c), the response of the segment was calculated to see which equivalent width gave the closest response to results to the FEA, which is presented in Figure 3-2. For Case (a): 𝑙 ′ = 51.54 𝑚𝑚, 𝑏 = 65 𝑚𝑚, 𝑏1 = 24 𝑚𝑚, 𝑙1 = 25 𝑚𝑚, 𝐸 = 1300 𝑀𝑃𝑎, then using Equation (30), we find 𝑏’ equal to 50.35 mm. Figure 3-7: The force-displacement curve of segment (a) from finite element analysis compared to the semi-analytical model 60 Figure 3-8: The fall in the internal energy of segment (a) from finite element analysis compared to the semi-analytical model For Case (b): 𝑙 ′ = 51.54 𝑚𝑚, 𝑏 = 65 𝑚𝑚, 𝑏1 = 32 𝑚𝑚, 𝑙1 = 20 𝑚𝑚, 𝐸 = 1300 𝑀𝑃𝑎, then using Equation (30) we find 𝑏’ equal to 54.84 mm. Figure 3-9: The force-displacement curve of segment (b) from finite element analysis compared to the semi-analytical model 61 Figure 3-10: The fall in the internal energy of segment (b) from finite element analysis compared to the semi-analytical model For Case (c): 𝑙 ′ = 51.54 𝑚𝑚, 𝑏 = 65 𝑚𝑚, 𝑏1 = 40 𝑚𝑚, 𝑙1 = 15 𝑚𝑚, 𝐸 = 1300 𝑀𝑃𝑎, then using Equation (30) 𝑏’ is equal to 58.75 mm. Figure 3-11: The force-displacement curve of segment (c) from finite element analysis compared to semi-analytical model 62 Figure 3-12: The fall in the internal energy of segment (c) from finite element analysis compared to the semi-analytical model By interpolation, it can be said that the length of the equivalent panel for Case (a) should be (54.54 mm), and for Case (b) it can be (52.54 mm). First, these differences are not very significant. Then, in the process of design, we can use the previously proposed relation for the equivalent panel’s length. Second, the equivalent length was determined based on the buckling force. The reason is that the analytical model is known to be accurate when it comes to predicting the axial strain and force corresponding to the buckling event, while the level of accuracy is not the same for the drop in force. 3.4.1. Verification of modeling framework Four non-uniform stiffness (NSD) cylinders were 3D printed and axially compressed to verify the predictive model for post-buckling response presented in Chapter 2. The geometry of the cylinders, with reference to Figure 2-1, is presented in Table 3-2. 63 Table 3-2: Geometry of NSD cylinders tested under uniform axial compression Test tb t1A t1B t1c t2 l1 ’1A ’1B ’1C l3 b1 b2 Unit (mm) CYL-1 0.5 1 1.3 4 1 30 10.56 7.15 1.33 60 35 65 CYL-2 0.5 1 1.5 4 1 30 10.56 5.75 1.33 60 35 65 CYL-3 0.5 1 2 4 1 30 10.56 3.75 1.33 60 35 65 CYL-4 0.5 1 2.5 4 1 30 10.56 2.65 1.33 60 35 65 The NSD cylinders were 3D printed from acrylonitrile butadiene styrene (ABS+) material using a Fortus 250MC printer (Stratasys Ltd., MN). All cylinders had an effective length of 100 mm and an inner radius of 50 mm. The cylinders were provided with 10 mm wide thickened edges at the top and bottom for connection to test platens to approximate clamped boundary conditions (the platens had grooves to fit the cylinder’s edges and constrain rotation and radial displacements.) Also, to avoid stress concentrations at the interfaces between thicker and thinner areas in the main cylinder body, the thickened areas had a stepwise change in thickness. The tests were conducted on a universal testing frame (Instron 5982) under displacement control. Loading consisted of a uniform quasi-static ramp of axial shortening up to 0.5 mm. 64 Figure 3-13: (a) 3D printed cylindrical shell and (b) the axial loading test setup For each of the cylinders in Table 3-2, a finite element model was prepared using the software ABAQUS. The modeling details were as described in Section 2.3.1.2, with two differences. First, the initial imperfection was simulated by using the cylinder’s first buckling mode with an amplitude of 0.05t. Second, the elastic modulus for the analytical model was determined by matching the initial stiffness of the experiment and finite element models, which resulted in an elastic modulus equal to1300 N/mm2. The force-deformation post-buckling response for each cylinder in Table 3-2 was predicted by using the proposed semi-analytical model presented in Chapter 2 for two cases: 1) by using Equation (10) for radial deformation and 2) by using Equation (12) for radial deformation . The framework presented in Figure 2-6 was followed whereby the NSD cylinder is separated into 65 panels, the equivalent length of each panel is calculated from Equation (3), and the post-buckling response of the A, B, and C equivalent panels is calculated as presented in Section 0. The response of a parallel system of linear and nonlinear springs is determined by adding the forces from all springs for equal displacement increments to the system. The axial post-buckling responses of the cylinders from experiments, finite element models, and the proposed semi-analytical model are plotted in Figure 3-14. A comparison of the responses in Figure 3-14 shows that the difference between corresponding displacements to the buckling events between the analytical model, FE model, and experiment is smaller than 0.02 mm, which is less than 4% of the total deformation. Interestingly, in general, the analytical model better matches the experimental data for the displacement at buckling events. For all cases, the response from the analytical model is stiffer than the experimental response and finite element solution. This is consistent with the results presented in Figure 2-11, which also shows the analytical post-buckling response of cylindrical panels being stiffer than that obtained from the FE model. It can be observed that the experimental traces in Figure 3-14 initially show a lower stiffness with a stiffening response after approximately 0.1 mm of total deformation. This behavior is primarily attributed to settlement effects of the test set up; that is, the imperfect contact between the interacting surfaces during loading. Another source of discrepancy for this response can be manufacturing imperfections, particularly at the transitions between thin and thick regions in the cylinder wall. 66 (a) (b) (c) (d) Figure 3-14: Force-deformation response of axially compressed NSD cylinders as predicted by FE and analytical models compared to experiments It can also be observed that the finite element analyses predict the magnitude of the force drops at the buckling events much better than the analytical model. The differences between the analytical and experimental responses can have several sources. The first and most important reason is the extensive level of assumptions made in the analytical model, and the second reason is that the initial imperfections in the 3D printed cylinder are different from what was assumed in the analytical model. 67 Chapter 4. Design of NSD cylinders with favored post-buckling response The main goal of this dissertation is to design NSD cylinders for a predefined elastic post- buckling response. This goal can be accomplished by following some backward steps in the proposed framework in Chapter 2. However, it is necessary to formulate the desired elastic post- buckling responses and then approach the design procedure as an optimization problem. The formulation for linear springs involved in the cylinder’s response is presented in Chapters 2 and 3 of this dissertation. Also, the equivalent length and width of the buckling panels are formulated based on the assumption presented in Chapter 3. Considering that the post-buckling responses of the panels are formulated as a function of their length, width, and thickness, it then becomes possible to mathematically connect the dimensions of the cylinder to its post-buckling response. As a result, a series of equations will represent the NSD cylinder behavior. Of course, the segment dimensions should add up to the geometry of the cylinder, which leads to several constraints for the optimization problem of designing the NSD cylinder. The only missing part of the design procedure is developing equations for the post-buckling response using a statistical approach. A pool of responses for numerous panels with varying geometry is needed to provide an empirical predictive model for any post-buckling response that enables the selection of panels that fulfill the design goal. Of course, the selected panels from design maps represent the equivalent panels described in the framework in Figure 2-6. It must be noted that the dimensions of the thickened areas on the top and bottom of the buckling panels affect some design factors, such as its stiffness and the axial shortening corresponding to the buckling event. After developing design maps, designing the cylinder with any desired elastic post- buckling response becomes an optimization problem that can be solved in different ways. The 68 optimization methods for each design parameter are not part of this study, and for each of these goals available routines in MATLAB were used. The design method just outlined was used in a recent study to obtain NSD cylinders with a specific radius and wall thickness [86]. In that study, the post-buckling responses were presented in contour maps that related each post-buckling response to the panel’s geometry. Here, such contour plots are called design maps. 4.1.1. Design maps for cylindrical panels with specific radius and thickness Design maps for buckling panels for different response parameters of NSD cylinders with a particular radius, elastic modulus, and body thickness were developed. They are: - The deformation corresponding to the buckling event in the panel. - The force corresponding to the buckling event. - The secondary stiffness of the panel. - Size of the force drop due to the buckling event. The availability of design maps for the characteristic features enables a designer to find panels with a desirable post-buckling response, such as a panel that buckles at a specific deformation or releases a particular amount of force from the buckling event. Due to the potential use of the NSD cylinders in energy dissipation systems, another valuable characteristic of a cylindrical panel is the dissipated hysteretic energy from the multi-buckling response. The semi-analytical model for cylindrical panels presented in Section 2.3 was used in this design. Design maps for cylindrical panels in a minimal range of dimensions and specific materials were developed. A parametric evaluation was conducted first, considering the post-buckling response of 16 cylindrical panels with changes in width and length from 20 mm to 80 mm. The radius and thickness of all panels are 50 mm and 0.5 mm, respectively. The material is ABS+ with 69 an elastic modulus (E) of 1300 N/mm2and Poisson’s ratio () of 0.3. The simulation was conducted by considering the contribution of 64 coupled modes (𝑚Λ = 𝑛Λ = 13). The axial stress-strain responses for these cases are provided in Figure 4-1. As it can be seen in the figure, the post-buckling response of the panels does not follow a linear pattern. However, after buckling, it is assumed that it is linear to simplify the response. Figure 4-1: Axial response of cylindrical panels from ABS material with R = 50 mm and t = 0.5 mm for varying length (l) and height (b);  and  are the global axial stress and strain of the panel, respectively Design contours were developed for four response parameters: a) buckling strain, b) buckling stress, c) stress drop, and d) secondary stiffness based on information from figure 4-1. The contour 70 curves were created by fitting a surface plot to the 16 data points from the parametric evaluation (for each design parameter). These plots can help select the buckling panels' dimensions based on the desired response parameter. (a) (b) (c) (d) Figure 4-2: Contour plot s for cylindrical panels from ABS material with t = 0.5 mm and R = 50 mm for varying length (l) and width (b): (a) buckling strain (0), (b) buckling stress (0), (c) stress drop (), and (d) post-buckling stiffness () 71 4.1.2. NSD cylinder with known end shortenings for each buckling event The design of an NSD cylinder with a targeted elastic post-buckling response under uniform axial compression using the design maps in Section 0 was showcased in a preliminary study [86]. The cylinder is to have an inner radius of 50 mm, a height of 100 mm, and the thickness of the buckling areas equal to 0.5 mm. The material is linear elastic with a modulus of 1300 N/mm2 and Poisson’s ratio of 0.2. The cylinder must experience three local buckling events when the global axial shortening reaches 0.24 mm, 0.32 mm, and 0.4 mm. Also, b, l1, l2, and b1 (see Figure 1-1) in all three segments are 60 mm, 30 mm, 10 mm, and 30 mm, respectively. The displacement corresponding to the buckling in segment A (Figure 2-2) is calculated as follows: 𝜎0 𝑏𝑡𝑏 𝛿0𝐴 = 𝜖0 𝑙1′𝐴 + 2 = 0.24 𝑚𝑚 (36) 𝑘1𝐴 where  is the end shortening corresponding to buckling in segment A. Also, 𝑙’1𝐴 , and 𝑘1𝐴 are expressed in Equations (3) and (4). In Equation (36), 0 is a function of l ( ’1A), and b; ’1A is a function of l1, 𝑡𝑏 , and t1A; 0 is a function of l ( ’1A) and b; and k1A is a function of b, b1, l1, l2, 𝑡𝑏 , and t1A. Thus,  is a function of six independent variables, from which five are predefined here and only t1A remains unknown. One can design a cylindrical panel that buckles under a 0.24 mm axial shortening by following a trial-and-error process. However, this is time-consuming. Thus, a numerical minimization problem using a genetic algorithm method was formulated using the software Mathematica [91]. The objective was to minimize (- 0.24)2. The same process was used to find t1B and t1C numerically minimize (- 0.32)2 and (C- 0.40)2. The obtained solutions were: t1A= 3.89 mm, t1B = 1.30 mm, and t1C =0.9 mm. 72 In this case, segments D do not affect the design objective, and they just should be stiff enough to separate the buckling events in segments A, B, and C. A thickness of 1.0 mm for t2 and a length of 60 mm for l3 separated the buckling segments well in previous experiments (Section 3.4.1), and thus the same values were used here. A finite element model for the NSD cylinder with the designed dimensions was developed in ABAQUS to evaluate the design. The axial response of the designed NSD cylinder from the finite element analysis is shown in Figure 4-3(a), and contour maps of the radial deformations are shown in Figure 4-3(b). The buckling events happen in the designed order, and the axial shortening corresponding to the buckling events matches well the design target and the analysis results. This example shows that it is possible to design NSD cylinders with a targeted post-buckling response using the presented analytical model. (a) (b) Figure 4-3: (a) Axial response of designed NSD cylinder for local buckling at the axial shortening of 0.24, 0.32, and 0.4 mm, (b) radial deformation of the cylinder at each buckling event 73 4.2. General Design maps Design maps provided in the previous section can help design an NSD cylinder with a desired post-buckling response of a cylinder with a specific radius and elastic modulus. However, the problem is that the dimensions of a cylindrical panel and cylinders (thickness, width, and length) have no limitations. It is beneficial to find the response of a cylindrical panel to compression in a closed-form solution for Equations (4) and (5). This is impossible for w (radial deformation) with more than one sinusoidal mode. However, it is possible to solve those equations independently from the radius (R) of the cylindrical panel if we know the ratios of length, width, and thickness to the radius. Thus, if l/R=, b/R= , and t/R=  are defined, R and E can be removed from the shell governing equations. Then, the ratio of radial deformation w to R can be expressed as follows: Λ 𝜋𝑥′ 𝑚𝑖 𝜋𝑥′ 𝜋y′ 𝑛𝑖 𝜋y′ 𝑤 = 𝑅 ∑ 𝜍𝑖 (sin ) (sin ) (sin ) (sin ) 𝜉1 𝜉1 𝜉2 𝜉2 𝑖=1 (37) Λ 𝜋𝑥′ 𝑚𝑖 𝜋𝑥′ 𝜋y′ 𝑛𝑖 𝜋y′ 𝑤 ∗ = 𝑅 ∑ 𝜇𝑖 𝜉3 (sin ) (sin ) (sin ) (sin ) 𝜉1 𝜉1 𝜉2 𝜉2 𝑖=1 where 𝜍𝑖 is 𝑎𝑖 ⁄𝑅 . We can further define 𝑤 = 𝑅𝑤 ̂ and 𝑤 ∗ = 𝑅𝑤̂ ∗ . Expressing the radial deformations this way leads to: 𝑤̂,𝑥′𝑥′ 𝑤̂,𝑥′𝑦′ 𝑤 ̂ ,𝜉2𝜉2 𝑤,𝑥𝑥 = , 𝑤,𝑥𝑦 = , 𝑤,𝑦𝑦 = , 𝑅 𝑅 𝑅 4 4 𝛻̂ 𝑤̂ 4 𝜕4 𝜕4 𝜕4 𝛻 𝑤= , ̂ = (𝛻 +2 + ) 𝑅 𝜕4 𝑥′ 𝜕2 𝑥′𝜕2 𝑦′ 𝜕4 𝑦′ ∗ ∗ ∗ ∗ 𝑅𝑤 ̂,𝑥′𝑥′ ∗ 𝑅𝑤 ̂,𝑥′𝑦′ ∗ 𝑅𝑤 ̂,𝑦′𝑦′ 𝑤,𝑥𝑥 = , 𝑤,𝑥𝑦 = , 𝑤,𝑦𝑦 = 𝑅 𝑅 𝑅 Rewriting Equation (5) results in: 74 𝛻̂ 4 𝑓̂𝐸 𝑤 2 ̂,𝑥′𝑦′ 𝑤 ̂,𝑥′𝑥′ 𝑤̂,𝑦′𝑦′ 𝑤 ̂,𝑥′𝑥′ 2𝑤 ̂,𝑥′𝑦′ 𝑤 ∗ ̂,𝑦′𝑦′ 𝑤̂,𝑥′𝑥′ 𝑤 ∗ ̂,𝑦′𝑦′ 𝑤 ̂,𝑦′𝑦′ 𝑤 ∗ ̂,𝑥′𝑥′ − 𝐸𝜉3 𝑅 ( 2 − − 2 + − − )=0 𝑅 𝑅 𝑅2 𝑅 𝑅2 𝑅2 𝑅2 𝛻4𝑓 𝑅 where 𝛻 ̂ 4 𝑓̂ = . Then, Equation (5) now can be written for variables 𝜉1 , 𝜉2 , 𝜉3 as follows: 𝐸 𝛻̂ 4 𝑓̂ − 𝜉3 (𝑤 2 ̂,𝑥′𝑦′ −𝑤̂,𝑥′𝑥′ 𝑤 ̂,𝑦′𝑦′ − 𝑤 ̂,𝑥′𝑥′ + 2𝑤 ̂,𝑥′𝑦′ 𝑤 ∗ ̂,𝑦′𝑦′ −𝑤̂,𝑥′𝑥′ 𝑤 ∗ ̂,𝑦′𝑦′ −𝑤 ̂,𝑦′𝑦′ 𝑤 ∗ ̂,𝑥′𝑥′ )=0 (38) By solving Equation (38), we can say: 3 𝐸𝑡3 𝐸𝜉3 𝑅3 𝑓,𝑥𝑥 = 𝑅𝐸𝑓̂,𝑥′𝑥′ , 𝑓,𝑥𝑦 = 𝑅𝐸𝑓̂,𝑥′𝑦′ , 𝑓,𝑦𝑦 = 𝑅𝐸𝑓̂,𝑦′𝑦′ , 𝐷 = = 12(1 − 𝜈2 ) 12(1 − 𝜈2 ) Then by replacing these expressions in Equation (6), we have: 3 𝐸 𝜉3 𝑅 3 𝛻̂ 4 𝑤 ̂ 12(1 − 𝜈 2 ) − 𝐸𝑅𝑓̂,𝑥′𝑥′ 𝑅3 𝐸𝑅𝑓̂,𝑦′𝑦′ (𝑤 ̂ ,𝑥′𝑥′ + 𝑤̂ ,∗𝑥′𝑥′ ) 2𝐸𝑅𝑓̂,𝑥′𝑦′ (𝑤 ̂ ,𝑥′𝑦′ + 𝑤 ∗ ̂,𝑥′𝑦′ ) 𝐸𝑅𝑓̂,𝑥′𝑥′ (𝑤 ̂ ,𝑦′𝑦′ + 𝑤̂ ,∗𝑦′𝑦′ ) −[ − + ]=0 𝑅 𝑅 𝑅 Factoring out E and R, the equation turns to be: 𝜉3 3 ̂ − 𝑓̂,𝑥′𝑥′ 𝛻̂ 4 𝑤 12(1 − 𝜈 2 ) (39) − [𝑓̂,𝑦′𝑦′ (𝑤̂,𝑥′𝑥′ + 𝑤 ∗ ̂,𝑥′𝑥′ ) − 2𝑓̂,𝑥′𝑦′ (𝑤 ̂,𝑥′𝑦′ + 𝑤 ∗ ̂,𝑥′𝑦′ ) + 𝑓̂,𝑥′𝑥′ (𝑤 ̂,𝑦′𝑦′ + 𝑤 ∗ ̂,𝑦′𝑦′ )] =0 Instead of solving Equations (5) and (6), we can solve equations (38) and (39), which are independent of R and E. Solving these equations, we can find 𝜍𝑖 in Equation (37) and 𝛿 ⁄𝐸 = 𝑃⁄𝐸𝑡𝑏 . Using these, we see the axial deformation ratio to the radius (𝛿 ⁄𝑅 ). This form responds to any given cylindrical panel in the format of 𝛿 ⁄𝐸 vs. 𝛿 ⁄𝑙 . The next step for designing an NSD cylinder for a desired post-buckling response is to define a numerical range for variables 𝜉1 , 𝜉2 , 𝜉3 that cover most of the possible design cases. 75 4.2.1. Determining acceptable range for panel parameters Cylinders with a wall thickness of less than 10% of the inner diameter are considered thin walled. The term “thin-walled” refers to the activation of the hoop and longitudinal stress under internal pressure, which is not the case for this research. This research focuses on thin cylinders under compression. In 1968 NASA conducted a series of experiments to define knockdown factor for buckling in cylinders [41]. They concentrated most of their experiments with 𝜉3 between 0.01 and 0.001. However, they had some experiments for 𝜉3 around 0.0003. This small ratio makes more sense for use in aerospace structures where the radius is significantly large, and thus a value of 𝜉3 =0.0003 can be a realistic thickness. Further, in the NASA report there are several studies that concentrated on 𝜉3 between 0.02 and 0.001 in their experiments on isotropic cylindrical shells under compression [92]–[101]. The goal here is to cover an as large area as possible for NSD cylinders that are not expected to have a large radius. A value of 0.001 for 𝜉3 means that a cylinder with a radius of 100 mm will have a thickness of 0.1 mm, which is the current limit for commercial polymer 3D printers. Thus, the minimum value of 𝜉3 was set to 0.001 based on manufacturing limitations. The largest 𝜉3 in the NASA study was 0.01. To cover the possibility of having a cylinder with a slightly larger value of 𝜉3 , the maximum number in this study is fixed to 0.013. Defining b/R (𝜉2 ) does not need further discussion since the circumference of a cylinder is at most 2𝜋𝑅, which is approximately 6.28𝑅. Then, considering that NSD cylinders have at least two NSD segments and two dividing segments, with the width of the dividing segment being half the width of the buckling segment, it leads to a width of approximately 2R for the buckling segment. 76 For the smallest ratio here, 0.4R is selected. The ratio of l/R (𝜉1 ) is defined based on its ratio with 𝜉2 from 0.1 to 10, and thus it is between 0.4R and 2R. Based on these limits, a set of 125 semi-analytical models was developed for five thicknesses of the panels in Table 4-1. All cases considered a Poison’s ratio of 0.3. Poison’s ratio has a minor effect on the response of the cylindrical panels under compression, and 0.3 is a reasonable value for any material unless we are dealing with incompressible material, which is beyond the scope of this study. 77 Table 4-1: Dimensions of cylindrical panels to develop empirical models l/R 0.4 0.8 1.2 1.6 2.0 b/R t/R 0.001 A-1-1 A-1-2 A-1-3 A-1-4 A-1-5 0.004 A-2-1 A-2-2 A-2-3 A-2-4 A-2-5 0.4 0.007 A-3-1 A-3-2 A-3-3 A-3-4 A-3-5 0.010 A-4-1 A-4-2 A-4-3 A-4-4 A-4-5 0.013 A-5-1 A-5-2 A-5-3 A-5-4 A-5-5 0.001 B-1-1 B-1-2 B-1-3 B-1-4 B-1-5 0.004 B-2-1 B-2-2 B-2-3 B-2-4 B-2-5 0.8 0.007 B-3-1 B-3-2 B-3-3 B-3-4 B-3-5 0.010 B-4-1 B-4-2 B-4-3 B-4-4 B-4-5 0.013 B-5-1 B-5-2 B-5-3 B-5-4 B-5-5 0.001 C-1-1 C-1-2 C-1-3 C-1-4 C-1-5 0.004 C-2-1 C-2-2 C-2-3 C-2-4 C-2-5 1.2 0.007 C-3-1 C-3-2 C-3-3 C-3-4 C-3-5 0.010 C-4-1 C-4-2 C-4-3 C-4-4 C-4-5 0.013 C-5-1 C-5-2 C-5-3 C-5-4 C-5-5 0.001 D-1-1 D-1-2 D-1-3 D-1-4 D-1-5 0.004 D-2-1 D-2-2 D-2-3 D-2-4 D-2-5 1.6 0.007 D-3-1 D-3-2 D-3-3 D-3-4 D-3-5 0.010 D-4-1 D-4-2 D-4-3 D-4-4 D-4-5 0.013 D-5-1 D-5-2 D-5-3 D-5-4 D-5-5 0.001 E-1-1 E-1-2 E-1-3 E-1-4 E-1-5 0.004 E-2-1 E-2-2 E-2-3 E-2-4 E-2-5 2.0 0.007 E-3-1 E-3-2 E-3-3 E-3-4 E-3-5 0.010 E-4-1 E-4-2 E-4-3 E-4-4 E-4-5 0.013 E-5-1 E-5-2 E-5-3 E-5-4 E-5-5 78 4.2.2. Design Maps For all the cases in Table 4-1 loading increases at a constant rate until a change in stiffness is observed. The stiffness difference does not mean that the panel experienced a buckling event, and buckling is recognized only if there is a significant force drop in a single deformation increment. The results from the semi-analytical model show several snap-back responses similar to the response that has been shown in Figure 2-10 (equilibrium path). However, knowing that the loading process is controlled by deformation, we can change the responses to Figure 2-10 (deformation control). Four types of responses were obtained for a given panel: • Force-displacement • Applied energy • Radial deformation ate the panel center • Maximum von-mises stress. Obtaining the force and displacement of the panel provides the ability to design a cylinder for three types of behavior: 1) the end shortening regarding the buckling event, 2) the drop in force during the buckling event, and 3) the post-buckling stiffness. From recording applied energy to the panel, one can find the drop in the energy from each buckling event. The radial deformation in the center of the panel is of interest to possibly make contact between a shell and a rigid wall or make contact between two different cylindrical panels. The maximum von-Mises stress provides three essential design factors that are seldom talked about in the literature about buckling in cylindrical shells; these are: I) the von-Mises stress before the buckling event, 2) the jump in the von-Mises stress after the buckling event, and 3) the rate of stress after buckling. Considering these essential factors lead to design maps for nine desired responses: 79 • Axial strain corresponding to the buckling event • Normalized axial force corresponding to buckling • Normalized force-drop due to buckling • Normalized stiffness of the panel after buckling • Maximum normalized von-Mises stress in the panel after the buckling event • Strain energy reduction in the buckling event • Normalized radial deformation at the panel center in the buckling event. 80 4.2.2.1. Axial Strain corresponding to the buckling Figure 4-4: Axial strain corresponding to the first buckling event for the panels with t/R varies from 0.001to 0.013 81 4.2.2.2. Axial stress corresponding to the buckling Figure 4-5: Axial stress corresponding to the first buckling event for the panels with t/R varies from 0.001to 0.013 82 4.2.2.3. Normalized force drop corresponding to the buckling Figure 4-6: Normalized drop in stress corresponding to the first buckling event for the panels with t/R varies from 0.001to 0.013 83 4.2.2.4. Normalized energy drop corresponding to the buckling Figure 4-7: Normalized drop in the strain energy corresponding to the first buckling event for the panels with t/R varies from 0.001to 0.013 84 4.2.2.5. Normalized secondary stiffness of the panel after buckling Figure 4-8: Normalized secondary stiffness after the first buckling event for the panels with t/R varies from 0.001to 0.013 85 4.2.2.6. Normalized radial deformation after the buckling Figure 4-9: Normalized radial deformation in the center of the panel after the first buckling event for the panels with t/R varies from 0.001to 0.013 86 4.2.2.7. Normalized change in the radial deformation after the buckling Figure 4-10: Normalized slope of radial deformation after the first buckling event for the panels with t/R varies from 0.001to 0.013 87 4.2.2.8. Normalized von mises stress after the buckling Figure 4-11: Normalized von-mises stress after the first buckling event for the panels with t/R varies from 0.001to 0.013 88 4.2.2.9. Normalized slope of the von mises stress after the buckling Figure 4-12: Normalized slope of von-mises stress after the first buckling event for the panels with t/R varies from 0.001to 0.013 89 4.2.3. Discussion The colorful design maps in the previous section show the distribution of design parameters in the defined range. These maps can be directly used to prepare equations for design, which are provided in the next section, and they provide constructive insight into the design process, which is discussed here. The strain maps for 𝑡𝑏 ⁄𝑅 = 0.001 and 0.004 do not show a pattern of change regarding the length and width of the panels. However, the difference between the maximum and minimum buckling strain values is very close to each other and the average buckling strain. This difference does not exceed 10% of the buckling strain in the worst case (𝑡𝑏 ⁄𝑅 = 0.001). For thicker panels, the buckling strain shows an inverse relation with the length and width of the panel. However, the higher buckling strain or normalized buckling stress for shorter and narrower panels is not clear since there is no significant load or energy drop. Thus, the points of buckling for these panels were judged by the author, and a slight difference can be a consequence from this judgment. Equal buckling strain and normalized stress in the panels with the same initial imperfections were experienced before by the author for a minimal set of data for specific radius and elastic modulus [86]. Figure 4-13 indicates the similarity between buckling stresses and strains for the 125 studied cases. This figure shows how increased thickness-to-radius ratio increases the average buckling strain and normalized stress. It also indicates that the gap between maximum and minimum buckling strain increases with panel thickness. 90 Figure 4-13: Normalized stress vs. strain response of the equivalent panels Such a similarity in the results is not observed for the force and energy drop regarding the buckling events. However, the design maps for these two parameters suggest that to see a sudden 91 load decrease in the NSD cylinder, one should avoid narrow, thick, and short equivalent panels. The design maps also suggest that the thicker the panel is, the range of the width and height with energy or force drop reduces. The design maps for force and energy drop answer the question about the number of buckling events in one NSD cylinder. The narrowest equivalent panel with force-drop has 𝑏/𝑅 = 0.5 for 𝑡/𝑅 = 0.001, which is the best case for force-drop. Also, to keep the buckling panels separated we assume the narrowest D zone with 𝑏𝑑 ⁄𝑅 is equal to 0.3. Then one buckling panel takes the width of the cylinder equal to 0.8 times the radius. Knowing that the circumference of a cylinder is 2𝜋𝑅, or roughly 6.28R, then, at best, only eight buckling events are possible with force drops in the segment. However, even in that case, since the load drop magnitude is small, a force drop in the total cylinder response is not observable. Thus, it can be confirmed that there are no NSD cylinders with nine buckling events without the buckling occurring in the divider regions. The secondary stiffness of the panels was harder to define since there is no clear distance after the first buckling to measure the force and displacement, and thus identifying the secondary stiffness relies heavily on visual judgment. In this study, the post-buckling (i.e., secondary) stiffness was obtained by recording the force at a displacement that was twice of the value at which the buckling event occurred. In some cases, other buckling events can occur between these two points. Also, the response after buckling is not linear, particularly for cases without additional buckling events and large end shortenings. These last two noted issues can cause errors in the design process. The normalized radial deformation follows a similar pattern to the force and energy drops due to buckling. This deformation is captured in the center of the panel right after buckling occurs. The goal of capturing the central deformation of the panels is to potentially design NSD cylinders that 92 can interact, via contact, with each lateral constraints such as concentric cylinders, additional panels, or rigid walls. The radial deformation of the already buckled panel keeps changing in the far post-buckling response. This far post-buckling radial deformation can be necessary for the design of cylinders with maximum deformation is needed. This situation can be experienced when the goal is to design a cylinder that touches another surface at its end after shortening. This can be a design objective for application like a smart switch or making contact between the surfaces of two close cylindrical panels. A valuable design parameter is the maximum von Mises stress in the panel after the buckling event. Here, this stress is normalized by the uniform axial stress in the panel right before the buckling. The reason for providing such a parameter is to avoid designing an NSD cylinder with the expectation of several elastic buckling events while in reality, the panel that buckled first is not acting elastic anymore. To design an elastic NSD cylinder, we want to avoid such behavior. The normalized von Mises stress is not enough to predict the maximum stress far after buckling. To monitor the cylinder’s strength, the rate of increase in the maximum von Mises stress is given by another design map (Figure 4-12). This map has low level of accuracy, similar to the secondary stiffness or even lower than that. Calculating the von Mises stress based on classic shell theory is based on several simplifications already. These simplifications along with the assumptions of the radial deformations make reading the von Mises stress far after first buckling more inaccurate. However, it may be generally said that the rate of increase in von Mises stress decreases with the thickness of the cylindrical panel. Using this slope from Figure 4-12 helps finding a rough approximation of the axial compressive strain that causes failure in the cylindrical shell. 93 4.2.4. Predictive models To design an NSD cylinder, each desired post-buckling response should have an equation that relates it to the dimensions of the equivalent panel. The dimensions of equivalent panels are formulated in Chapter 3 and the design maps cover the post-buckling of the possible range of panel size. It then becomes possible to relate the post-buckling response of the cylinder to its geometry. For any desired post-buckling response, based on the equations provided in Chapter 3, we need a relation between variables (𝑡𝑏 ⁄𝑅 , 𝑙′⁄𝑅 , 𝑏′⁄𝑅 ) and several design factors in the maps. Thus, each of these design factors need to be expressed as a function of the variables. For instance, the secondary stiffness and force drop from the design maps are needed to predict the force drop in one segment during the buckling event. Another example is energy loss. To predict the energy loss in one segment, an energy drop map and secondary stiffness map are needed. The author could not define any physically meaningful equation that expresses the relation between variables and the design maps. However, the data set behind the maps cover a very vast range of responses. It seems logical that if an equation is provided for each map that fits the data well, then it can predict the response of any panel in the range of the variables with reasonable accuracy. Empirical models have been used for many years by researchers and engineers by fitting an equation to a data set when no closed-form solution is available. Several techniques for such fitting include genetic algorithms, genetic programming, and neural networks. Several software packages are available for each of these methods. Here, the genetic programming package, GPTIPS2, under the platform of software MATLAB was used [102]. The genetic programming package provides equations that we can later apply to the optimization problem when designing the NSD cylinder. 94 In this study, for any of the design parameters there are 125 data points. However, the GPTIPS2 package needs another data set to test the predictive model. Then, the design parameters for another 25 models were developed to serve as a test set. The dimensions of the panels were as given in Table 4-2. Table 4-2: Dimensions of the cylindrical panels used for the test in genetic programming Case l/R b/R t/R Name T-1 0.62 0.168 0.011 T-2 0.133 0.71 0.007 T-3 0.129 0.85 0.013 T-4 0.131 0.152 0.005 T-5 0.160 0.124 0.004 T-6 0.136 0.181 0.010 T-7 0.123 0.145 0.009 T-8 0.144 0.128 0.002 T-9 0.169 0.51 0.006 T-10 0.83 0.190 0.010 T-11 0.130 0.95 0.009 T-12 0.185 0.179 0.007 T-13 0.185 0.144 0.008 T-14 0.161 0.172 0.003 T-15 0.55 0.123 0.007 T-16 0.89 0.118 0.006 T-17 0.163 0.175 0.001 T-18 0.154 0.198 0.004 T-19 0.123 0.179 0.003 T-20 0.85 0.89 0.002 T-21 0.198 0.75 0.010 T-22 0.185 0.71 0.012 T-23 0.154 0.113 0.002 T-24 0.177 0.156 0.009 T-25 0.173 0.139 0.001 The dimensions in Table 4-2 were randomly generated in the range of the design maps. Thus, the GPTIPS algorithm was used with 125 training data sets and 25 test data sets. The resulting 95 equation for each design parameter is provided in Table 4-3. The table shows the 𝑅 2 the factor for each of the design parameters and the predictive equation. The 𝑅 2 is one factor that shows how good a fit is in the predictive equation. Table 4-3: Relation between design parameters and design variables (l/R, b/R, t/R) response Predictive equation 𝑹𝟐 65.5𝜉33 0.0000494𝜉1 + 0.0000494𝜉2 + 0.602𝑥3 + − 0.0395𝜉1 𝜉3 𝜉22 𝜖0 0.99 𝜉33 (32.7𝜉3 − 18.8𝜉1 + 32.7(𝜉33 ) + 69.8) + 0.000027𝜉13 + 𝜉23 − 0.000268 𝜎0 66.7𝜉33 0.553𝜉3 − 0.0000442𝜉1 + 3 − 0.0325𝜉15 𝜉32 − 𝜉2 (0.0162𝜉14 𝜉32 0.99 𝐸 𝜉2 − 0.000075) − 0.000111 7.23(𝜉1 − 𝜉3 )(𝜉2 + 2𝜉3 ) 0.0000253(𝜉2 + 2𝜉3 )2 − 9.01𝜉1 + 9.01𝜉2 + 27𝜉3 𝜉3 0.000542(1.48𝜉2 + 2.23𝜉3 ) − − 15.7𝜉3 Δ𝜎0 (𝜉16 )(𝜉2 − 𝜉3 ) 0.89 𝐸 − 0.005𝜉1 𝜉2 (𝜉1 − 𝜉3 )3 0.00000678(1.16𝜉2 + 2.33𝜉3 + 7.36) + 𝜉3 (𝜉1 + 𝜉3 )(𝜉1 − 𝜉3 ) 0.297𝜉3 (𝜉22 + 7.13) + − 0.112 (𝜉13 )(2𝜉1 + 𝜉2 ) 3.29𝜉1 + 3.39 𝜉3 − 0.319 𝜉1 (2 𝜉2 + 2 𝜉3 ) 𝜉2 0.103𝜉1 1.01 𝜉3 − 0.0000131 (𝜉1 − 4.66) ( + 𝜉23 ) − − 𝜉3 𝜉2 𝜉2 𝐸𝑛 0.93 𝐸𝑙𝑏𝑡 30(𝜉1 + 𝜉3 + 𝜉1 𝜉3 ) 3.39 𝜉1 𝜉3 0.82 𝜉1 𝜉3 − − + 𝜉2 + 𝜉22 + 10 𝜉2 𝜉23 + 0.00325 (𝜉2 − 6)(𝜉1 − 𝜉2 )(𝜉1 + 𝜉2 − 6) + 0.0397 96 Table 4-3 (cont’d) 0.00231𝜉23 − 0.00231𝜉22 − 813𝜉1 𝜉2 𝜉32 + 383 𝜉1 𝜉22 𝜉32 + 7.7𝜉1 𝜉3 (2𝜉2 − 2𝜉3 + 𝜉12 𝜉3 ) − 0.704 𝜉1 𝜉3 (9𝜉22 − 2𝜉2 + 7.1) 𝑤0 0.9 𝑅 5955𝜉33 (𝜉1 + 2 𝜉3 ) 0.000634 (2𝜉1 − 6.38)(𝜉2 − 𝜉3 )2 + + 𝜉2 𝜉1 𝜉2 + 0.00816 1.38(𝜉1 − 𝜉23 ) 12.1𝜉3 − 0.105𝜉2 − 0.00863𝜉1 − 72.1𝜉3 − 72.1𝜉2 + 300𝜉2 𝜉3 + 72.1𝜉23 0.0173𝜉1 0.105𝜉3 𝑘𝑠2𝑙 − + + 0.000221(2𝜉1 + 𝜉23 − 4.06)3 𝜉2 𝑥1 0.78 𝐸𝑏𝑡 𝜉3 47.1𝜉2 𝜉32 12.1𝜉3 (𝜉1 − 𝜉2 )2 − 270𝜉2 (𝜉32 ) (𝜉2 − ) − + 𝜉1 𝜉1 3.76 𝜉1 (𝜉2 + ) 𝜉2 + 0.243 23.1 0.000646 + 0.00261(𝜉2 + 𝜉13 − 5.33)3 + 1.13𝜉1 + 5.34 𝜉3 6.09𝜉3 𝑣𝑚 − 0.611𝜉1 (𝜉2 − 6.32 𝜉3 ) − 4.5𝜉1 − 2.35 0.87 𝐸 6.74 𝜉2 (𝜉1 + 5.53𝜉3 ) + + 696𝜉2 𝜉32 (𝜉2 − 5.14) 𝜉1 + 2𝜉2 + 3366𝜉23 𝜉33 (𝜉2 − 5.11) (𝜉2 − 5.14) − 3.57 The 𝑅 2 factors for the buckling strain and normalized buckling stress are less than 0.01, shy of 1.0. This describes the same point as what has been said in the discussion. 4.3. Design problems This section presents the design of four cylinders for specific post-buckling responses. Two examples are provided for controlling the secondary stiffness of the cylinder. This can be used in structural or mechanical systems to control the behavior of a system from sudden loading. A very 97 general example can be a building in which the designer aims to actively manage its stiffness after certain deformation and redistribute forces in the building’s stories. The other example is about designing a cylinder for given energy dissipation. This type of cylinder can be used to dissipate energy from axial loading. One smart application could be to design devices for specific energy dissipation. Controlling the amount of force drop in the buckling event can be a purpose for a designer of an NSD cylinder. Finally, a cylinder with a cylindrical stiffener is another possible design that can increase the number of buckling events. 4.3.1. NSD cylinder with plateau post-buckling response The goal is to design a cylinder with three buckling events where all the buckling events happen at the same force by using the design maps provided in the previous section. The schematic force- displacement of such behavior is shown in Figure 4-14. Figure 4-14: Schematic plateau post-buckling response of the NSD cylinder 98 The plateau post-buckling response means the buckling events should occur at the same force but at different displacements. The equation to relate the force-displacement at the first buckling follows the following formulation. 𝛿01 𝛿01 𝑃1 = 𝑃01 + + + 3𝛿01 𝑘𝐷 1 1 1 1 (40) + + 𝑘𝑠12 𝑘12 𝑘𝑠13 𝑘13 The stiffnesses of the linear springs in this equation comes from Chapter 3. These stiffnesses are a function of the geometry of the NSD cylinder. 𝑃01 𝑃01 = 𝜎01 𝑏1′ 𝑡𝑏 and 𝛿01 = 𝜖01 𝑙1′ + (41) 𝑘11 Also, using genetic programming (GP) 𝜎01 is defined as a function of 𝑙’, 𝑏’, and 𝑡𝑏 . After the first buckling, the system's stiffness adjusts because the first buckled panel shows secondary stiffness instead of its primary stiffness. 𝛿02 − 𝛿01 𝛿02 𝑃2 = 𝑃01 − Δ𝑃1 + + 𝑃02 + + 3𝛿02 𝑘𝐷 1 1 1 1 + + 𝑘𝑠21 𝑘11 𝑘𝑠13 𝑘13 (42) 𝑃02 = 𝜎02 𝑏2′ 𝑡𝑏 𝑃02 𝛿02 = 𝜖02 𝑙2′ + 𝑘12 the second buckling, the system's stiffness adjusts because the first buckled panel and second buckled panel show secondary stiffness instead of their primary stiffness. 𝛿03 − 𝛿01 𝛿03 − 𝛿02 𝛿03 𝑃3 = 𝑃01 − Δ𝑃1 + + 𝑃02 − Δ𝑃2 + + + 3𝛿03 𝑘𝐷 1 1 1 1 1 1 + + + 𝑘𝑠21 𝑘11 𝑘𝑠22 𝑘12 𝑘𝑠13 𝑘13 (43) 𝑃03 = 𝜎03 𝑏3′ 𝑡𝑏 𝑃03 𝛿03 = 𝜖03 𝑙3′ + 𝑘13 99 These equations define the force at the buckling events. Now, if we decide to design a cylinder that shows a bilinear post-buckling response, the force and displacement corresponding to buckling events should be located in a line. This condition mathematically means: 𝑃3 − 𝑃2 𝑃2 − 𝑃1 = =𝜆 (44) 𝛿03 − 𝛿02 𝛿02 − 𝛿01 Then, based on the need of the design, we can define this secondary slope of force- displacement of the NSD cylinder. Here, in the first case, we assume that 𝜆 = 0, which means we are looking for a case with a force plateau after the first buckling event. Of course, this is not the only restriction for the design. We want the buckling occurs in three different end shortenings, which means: 𝛿01 < 𝛿02 < 𝛿03 Since we use an algorithm to find the optimum design, it is probable that the program just selects very close-end shortenings, which is not desirable for our design purpose. Then, we should decide how close these buckling events can be close to each other. Here, it is assumed that 𝛿02 − 𝛿01 ⁄ℎ > 0.0005 and 𝛿03 − 𝛿02 ⁄ℎ > 0.0005. It is also assumed that the material of the NSD cylinder is PLA, with E =1300 N/mm2 and  = 0.3. The geometry of the cylindrical segments should be compatible with the dimensions of the cylinder. Then: 𝑏01 + 𝑏02 + 𝑏03 + 3𝑏𝐷 = 2𝜋𝑅 0 < 𝑏11 < 𝑏01 , 0 < 𝑏12 < 𝑏02 , 0 < 𝑏13 < 𝑏03 ℎ01 + 2𝑙11 < ℎ , ℎ01 > 0 𝑎𝑛𝑑 𝑙11 > 0 ℎ02 + 2𝑙12 < ℎ , ℎ02 > 0 𝑎𝑛𝑑 𝑙12 > 0 ℎ03 + 2𝑙13 < ℎ , ℎ03 > 0 𝑎𝑛𝑑 𝑙13 > 0 0 < 𝑙𝐷 < ℎ 100 There are also restrictions regarding the production of the designed cylinder. For instance, it is desired to have a radius equal to 50 mm, a height of 120 mm. and a wall thickness of 0.6 mm. These are due to the available test setup we have for the 3D printed cylinders and the available 3D printer, which not have adequate quality for thicknesses less than 0.6 mm. Another important constraint for the dimensions of the cylinder is the size of the equivalent panels, which should remain in the range of 0.4 and 2.0 times the radius of the cylinder. Based on these conditions, a Mathematica code was developed to find the geometrical unknowns to minimize (𝑃3 − 𝑃2 )2 + (𝑃2 − 𝑃1 )2 + (𝑃3 − 𝑃1 )2 . Numerous variables result in having an unlimited number of designs. Then, the designer should provide some initial input. For example, using the same process by just changing the thickness of the body and length of the buckling zone, three cylinders were designed to show a plateau after the first buckling event. Initial dimensions for the cylinder were as follows: ℎ = 120𝑚𝑚, 𝑅 = 50 𝑚𝑚, 𝑡𝑏 = 0.6 𝑚𝑚 ℎ01 = ℎ01 = ℎ01 = 70 𝑚𝑚 𝑙11 = 15 𝑚𝑚, 𝑙12 = 15 𝑚𝑚 , 𝑙13 = 15 𝑚𝑚, 𝑙𝑑 = 90 𝑚𝑚, 𝑡𝑑 = 1.2 𝑚𝑚, 𝑡𝑏𝑑 = 0.6 𝑚𝑚 Using the noted optimization algorithm, the dimensions were determined as follows: 𝑏01 = 90 𝑚𝑚, 𝑏02 = 70𝑚𝑚, 𝑏03 = 55 𝑚𝑚 𝑏11 = 67.53 𝑚𝑚, 𝑏12 = 54.97 𝑚𝑚, 𝑏13 = 37.26 𝑚𝑚 𝑡11 = 3.95 𝑚𝑚, 𝑡12 = 1.85 𝑚𝑚 , 𝑡13 = 1.23 𝑚𝑚 A finite element model was developed based on these dimensions. The designer expects that segment one buckles before segment two and segment two buckles before segment three. The design results expect the buckling force for all these segments to be 1573 N. 101 Figure 4-15: Force-displacement response of the designed NSD cylinder for plateau post- buckling from FEA in comparison to one from the semi-analytical model The finite element model shows two buckling events in the same load instead of three. The third buckling occurs right after the second and obviously under a lower axial load. The problem is that the third segment buckles earlier than expected. The author believes that the reason for this premature buckling is the boundary conditions of the third segment. In the proposed semi- analytical model, the side boundaries of the buckling segment are assumed to be free in the axial direction. It means elements D axially deform in harmony with the center of the buckling segment. In reality, this assumption is not accurate since the stiffness of segment D and that of the buckling zone are significantly different. 102 Figure 4-16: The shear resulted from D zones deforming less than the center of the segment The difference in deformation causes a shear force applied to the stiffened area limits. In the hardened areas with high thickness and a relatively large distance between the center and the edge, the effect of this shear is neglected. However, when the thickened area is relatively thin, like segment three of the current design, this shear generates initial imperfections, and the segment buckles prematurely. Figure 4-17 shows the axial deformations in all three segments of the designed NSD cylinder. Segments one and two are not affected significantly in the center, while segment three shows the concentration in the middle. This concentration of deformations is the reason for the premature buckling. 103 Figure 4-17: Axial deformation distribution before the first buckling in three segments of the designed NSD cylinder To fix the problem, the D areas surrounding buckling area three were divided in half. The external half remains the same. For the internal half, the thickness of the top and bottom of zone D was adjusted to 1.2 mm (the thickest of zone D). This adjustment causes the deformation of zone D, and the buckling panel to be close to each other. A finite element model based on this adjustment was developed, and the results are shown in Figure 4-18. The finite element analysis results for the designed cylinder show an approximately flat post- buckling response with the following forces for each buckling event: 𝛿01 = 0.475 𝑚𝑚, 𝐹01 = 1472 𝑁 𝛿02 = 0.547 𝑚𝑚, 𝐹02 = 1481 𝑁 𝛿03 = 0.619 𝑚𝑚, 𝐹03 = 1454 𝑁 If the first and last buckling force peaks are connected with a line, the slope of the line will be: 1454 − 1472 𝑁 𝐾𝑠𝑒𝑐𝑜𝑛𝑑𝑎𝑦 = = −125 0.619 − 0.475 𝑚𝑚 104 The initial stiffness of the cylinder is calculated as: 1472 𝐾𝑖𝑛𝑖𝑡𝑖𝑎𝑙 = = 3098.95 0.475 It makes the secondary stiffness -4% of the initial stiffness, which is close to the design goal. 𝐾𝑠𝑒𝑐𝑜𝑛𝑑𝑎𝑟𝑦 125 ⟹ =− = −0.04 𝐾𝑖𝑛𝑖𝑡𝑖𝑎𝑙 3098.95 Figure 4-18: Force-displacement response of the designed NSD cylinder for plateau post- buckling from modified FEA in comparison to one from a semi-analytical model 105 Figure 4-19: Axial deformation distribution before the first buckling in segments three of the (a) designed NSD cylinder, and (b)modified NSD cylinder The axial deformation of segment 3 in the model with 0.6 mm at the top and bottom is compared to axial deformation in the same segment in the adjusted segment. Figure 4-19 confirms the concentration of axial deformation in the center of the initially designed cylinder and a more even distribution of axial deformation in the modified NSD cylinder. The radial deformation of segment three in the modified and originally designed NSD cylinders are compared. As shown in Figure 4-20, the radial deformation of the primarily designed segment three is almost twice the radial deformation of the modified cylindrical segment. Of course, significant initial imperfections result in premature buckling in the segment. 106 Figure 4-20: Radial deformation in segment three of (a) primarily designed cylinder, and (b) modified cylinder, before the first buckling in segment one The other interesting difference between the force-displacement responses from the FEA and the designed response using the semi-analytical model is the size of the force drop in the second and third buckling events. Since the comparison of the semi-analytical model and FEA for the individual segment showed good agreement, it is believed that the size of this force drop in the FE model is related to later buckling events in segments one and two. In this example, the width of segments one and two are very large, and oversized panels are subjected to numerous buckling events. Considering numerous buckling events and stress, strain, and stress drops for each of those events does not help the design process. The complication level rises to the extent that, with the proposed method, it would not be possible to achieve a design. Thus, the assumption of one buckling event per segment is maintained and therefore the aim is to select the dimensions in a way consistent with this assumption as much as possible. 107 4.3.2. NSD cylinder with equal force drops in buckling events Control over the deformation corresponding to buckling events and the size of the load drop could be one desirable response. For instance, an NSD cylinder with a constant spacing between buckling events and constant load drop is designed here. The schematic force-displacement response for such a cylinder is demonstrated in Figure 4-21. Figure 4-21: Schematic post-buckling response of the NSD cylinder where identical force drops occur with the identical difference in end shortening for buckling events 108 Considering equations in Chapter 3, the size of the force falls off follows: Δ𝑃1 = Δ𝑃0 + (𝛿𝑠 − 𝛿0 )𝑘𝑠 ’ (45) In the equation above Δ𝑃0 , 𝑘𝑠′ , and 𝛿0 are calculated from design maps, and 𝛿𝑠 is calculated from the following equation: 𝑘1 𝛿𝑡 + 𝑘𝑠′ 𝛿0 − (𝑃0 − Δ𝑃0 ) 𝛿𝑠 = 2 (46) 𝑘1 ′ 2 + 𝑘𝑠 The requirement of equal force drops in the cylinder is expressed as: Δ𝑃11 = Δ𝑃12 = Δ𝑃13 = Δ𝑃14 (47) 𝛿𝑡4 − 𝛿𝑡3 = 𝛿𝑡3 − 𝛿𝑡2 = 𝛿𝑡2 − 𝛿𝑡1 Such that for any given segment: 2𝑃0 𝛿𝑡 = 𝛿0 + (48) 𝑘1 Here, we design an NSD cylinder with four buckling segments with a radius of 65 mm, a height of 100 mm, and a body thickness of 0.6 mm. Initial dimensions for the cylinder are as follows: ℎ = 100𝑚𝑚, 𝑅 = 65𝑚𝑚, 𝑡𝑏 = 0.6 𝑚𝑚 ℎ01 = ℎ02 = ℎ03 = ℎ04 = 50 𝑚𝑚 𝑙11 = 𝑙12 = 𝑙13 = 𝑙14 = 15 𝑚𝑚 𝑙𝑑 = 90 𝑚𝑚, 𝑡𝑑 = 1.2 𝑚𝑚, 𝑡𝑏𝑑 = 1.2 𝑚𝑚 𝑏01 = 𝑏02 = 𝑏03 = 𝑏04 = 58 𝑚𝑚 Using the optimization algorithm, the final design dimensions are: 𝑏11 = 46.76 𝑚𝑚, 𝑏12 = 43.22 𝑚𝑚, 𝑏13 = 38.32𝑚𝑚, 𝑏14 = 29.27 𝑚𝑚 𝑡11 = 3.70 𝑚𝑚, 𝑡12 = 2.49 𝑚𝑚 , 𝑡13 = 1.87 𝑚𝑚, 𝑡14 = 1.51 𝑚𝑚 109 Figure 4-22: Force-displacement response of the designed NSD cylinder for similar load falls in similar distances from FEA in comparison to one from a semi-analytical model In Figure 4-22, the blue color shows the design using the semi-analytical model. The distance between buckling events and the size of the drops were designed as: 𝑑𝛿 = 0.02 𝑚𝑚, 𝑑𝑝 = 40 𝑁 The brown line shows the response of the designed cylinder from finite element analysis, which yielded the following results: 𝑑𝛿1 = 0.017 𝑚𝑚 , 𝑑𝛿2 = 0.019 𝑚𝑚 , 𝑑𝛿3 = 0.034 𝑚𝑚 𝑑𝑝1 = 58 𝑁 , 𝑑𝑝2 = 55 𝑁 , 𝑑𝑝3 = 54 𝑁 , 𝑑𝑝4 = 61 𝑁 In this case, the designed cylinder meets the demand of having identical force drops in the buckling events. However, the last buckling occurs at a deformation 0.014 mm bigger than desired. This error can be addressed by slight modification in the fourth segment’s stiffened areas. Here, we only change the thickness of the thickened areas to the following values. 𝑡11 = 3.70 𝑚𝑚, 𝑡12 = 2.49 𝑚𝑚 , 𝑡13 = 1.85 𝑚𝑚, 𝑡14 = 1.67 𝑚𝑚 110 Axial response from finite element analysis of an NSD cylinder with these thicknesses was compared to designed post-buckling response. The distance between buckling events gets much closer to the design while the drop in force still is higher than the designed force drop. From the comparisons in Chapter 3 of this document, we expect our model to underestimate the size of the force drop. Figure 4-23: Force-displacement response of the designed NSD cylinder for similar load falls in similar distances from FEA in comparison to one from a semi-analytical model Figure 4-23 compares the axial post-buckling response of the modified finite element model (brown line) and the designed response (blue line). The distance between buckling events and force drops in the response from the finite element model are measured as follows. 𝑑𝛿1 = 0.017 𝑚𝑚 , 𝑑𝛿2 = 0.02 𝑚𝑚 , 𝑑𝛿3 = 0.018 𝑚𝑚 𝑑𝑝1 = 53 𝑁 , 𝑑𝑝2 = 52 𝑁 , 𝑑𝑝3 = 57 𝑁, 𝑑𝑝4 = 60 𝑁 This comparison shows the effectiveness of the proposed design method and design tables. However, the proposed model significantly underestimates the force drop in the buckling events. 111 Chapter 5. Summary and Conclusion 5.1. Summary This study developed a semi-analytical model to predict the elastic post-buckling response of cylinders with non-uniform stiffness distribution (NSD) under compression. NSD cylinders show promise for predictability of their axial compression response and may find applications to provide structures or devices with functionalities such as load-limit control, energy dissipation or energy harvesting. A challenge to realize the potential of NSD cylinders is to adequately understand their elastic post-buckling response so that it may be controlled and designed for. The developed semi- analytical model takes three steps to predict the post-buckling response: • Separate the cylinder into parallel segments, • Simplify and predict the response of each part, and • Integrate the response of individual segments. To predict the post-buckling response of a cylindrical segment, its geometry is simplified to a cylindrical panel with uniform thickness connected in series with linear springs at the top and bottom of the panel. The side boundary of the panels is assumed fixed except for the ability to move in the axial direction. Based on the assumed boundary conditions and classical shell theory, the elastic post-buckling response of a cylindrical panel is solved as a boundary value differential equation using the pseudo-arc length method. The classic differential equation of the axial compression of the cylindrical panels was solved independent of the cylinder radius and elastic modulus of the material. The results were then normalized by the radius and elastic modulus for a wide range of the cylinder dimensions. These results allowed the development of design maps for several post-buckling responses such as axial strain and stresses corresponding to the first buckling event, force and energy drops from the 112 buckling event, the secondary (or post-buckling) stiffness of the panel, the radial deformation at the panel center, and the maximum von Mises stress in the panel. Three cylinders were designed, and finite element analyses validated the targeted post-buckling responses to show the effectiveness of the design maps and the design procedure. One NSD cylinder was designed to undergo several buckling events under compression at pre-defined end shortenings. A second NSD cylinder was designed to feature a post-buckling force-deformation response that plateaus at a constant force level. The third cylinder was designed to experience the same force drop at each buckling event and in identical axial end shortenings after the first event. 5.2. Conclusions The main goal of this dissertation was to predict the elastic post-buckling response NSD cylinders (cylinder with non-uniform stiffness distribution) and to design NSD cylinders for a desired elastic post-buckling response. Both goals were achieved, and several validating examples were provided for each. In the process of meeting the major goal of designing NSD cylinders, a semi-analytical model was developed for single cylindrical panel with unique stress distribution conditions. This model in combination with the idea of linear springs in series provided a very clear understanding about what happens to an NSD cylindrical segment under axial compression and why the load drops and energy drops observed from FEA models are different than those from single uniform panel. In Addition, design maps which provided for the panels independent from elastic modulus and radius of the cylinder provide useful information about several post-buckle response of an individual panel. The information includes the geometries that one can expect the buckling occurs in the panel, the maximum amount of drop in the force during a buckling in the panel, the maximum drop in the energy for different geometries, and the expected radial deformation in the 113 center of the panel. This information are not only useful for designing NSD cylinders but also they can be used in the process of developing new ideas to use cylindrical panel’s nonlinear response. The presented research provides an insight into the snap back behavior of slender structures under axial compression and how to approach designing such a structure for a desired far post- buckling response. The proposed framework, concept of decomposing NSD cylindrical segments into linear and nonlinear springs in series, proposed semi-analytical model for NSD equivalent panels, and idea of developing maps for nonlinear design parameters improve the ability to design smart structures relying on structural instabilities. This work expands the harnessing of structural instabilities to the area of thin-shell buckling phenomena which received less attention in comparison to other instabilities with less imperfection sensitivity. 5.3. Potential future research The predictive and design models presented in this dissertation paved the road to using NSD cylinders for smart applications. Future research on this topic could follow three paths. The first direction is to improve the predictive model and the proposed model's accuracy. The improvement can be achieved by replacing the fixed boundary conditions with flexible boundaries and developing design maps for different levels of flexibility along the edges. Other valuable research could target the dividing panels and investigate their effect on the post-buckling response, as there is a lack of knowledge on efficient dimensions for the dividing segments in the NSD cylinders. The numerical solution's efficiency is another area worth of more attention. Research in this area can help make the solution faster and consider facilitate considering more radial deformation modes, which improves the quality of the response. The size of the force drop, and energy drop in the proposed semi-analytical model is less than what finite element analysis 114 suggests. Addressing the source of this difference can also be a reasonable extension of this research. The second research path can be expanding the NSD cylinders in different ways. For example, an NSD cylinder with cylindrical stiffeners was preliminarily evaluated in a parallel study for possible stiffness regain by contact mechanism. It was found that it is impossible to make contact between the face-to-face cylindrical panels. However, it is worth investigating other possibilities of making NSD cylinders without manipulating the thickness of the segments. Using different materials with the same thickness may be as suitable approach worth further research. On the other hand, the current research only considered isotropic materials, and studying the possibility of designing cylinders from orthotropic materials can be a valuable future research topic. The third direction for future research is implementing the proposed NSD cylinder design for smart applications. One example can be energy harvesting by connecting piezoelectric PVDF film to the buckling areas of the NSD cylinders. The effect of the PVDF layer should be considered, which is possible by using equivalent elastic modulus and thickness for functionally graded materials. 115 APPENDIX 116 I prepared a finite element model for the cylinder shown in Figure A-1: Axial force- displacement from FE model and experiment, the cylinder is assumed to be made from acrylonitrile butadiene styrene (ABSPlus) with 1.6 GPa Young’s modulus and 0.2 Poisson’s ratio. The shells were modeled with four-node quadrilateral finite membrane strain elements with reduced integration (S4R). The displacement of the edges in the radial direction has been restrained for both ends. The bottom edge of the cylinder is fixed for axial displacement, while axial displacement has been applied to the upper edge. Both top and bottom edges are fixed for rotating around the circumference of the section. All the shell elements are aligned by their internal face. An Implicit dynamic analysis with Quasi-static application has been conducted to compare the force-displacement with the results from the experiment. The force-displacement from FE analysis shows good agreement with the one from the experiment. (Figure A-1) Similar to the experiment, the sequence of the buckling events is in zones 1, 2, and 3. This sequence shows that buckling first occurs between the stiffest boundary conditions. The behavior of the polymers does not follow the simple von misses stress criteria. However, in the absence of several influential factors, the von misses stress still could be the best option to estimate the status of material regarding failure. 117 Figure A-1: Axial force-displacement from FE model and experiment Here, the maximum stress appears around the first buckled zone. After the third buckling event, the maximum von misses stress turns to 28 MPa, close to ABSPlus maximum tensile strength (32 MPa). This very primary investigation shows that it is reasonable to neglect the damage in the material for the first three buckling events. 118 In contrast, further, loading could damage the material and defy the linear elastic behavior assumptions in the finite element model. Figure A-2 b-d shows the in-plane at the maximum load (𝛿 = 0.5 𝑚𝑚). (a) (b) (c) (d) Figure A-2: Stress distribution in the cylindrical shell for the maximum load (𝛿 = 0.5 𝑚𝑚) S11 shows the axial stress in the cylinder, S22 shows the stress in circumference of the cylinder, and S12 shows the shear stress. 119 (a) (b) (c) Figure A-3: Deformation of the cylinder for the maximum load (𝛿 = 0.5 𝑚𝑚) in; (a) radial direction, (b) circumferential direction, and (c) longitudinal direction Figure A-3 indicates the deformations of the cylinder where U1 is radial deformation, U2 is circumferential deformation, and U3 is axial deformation. Figure A-3 (a) shows the biggest radial deformations in the center of zones 1, 2, and 3. However, it is hard to understand the function of the deformations without emphasizing specific sections in the radial plates (ℎ = 𝑐𝑜𝑛𝑠𝑡𝑎𝑛𝑡) and normal to radial plates (𝜃 = 𝑐𝑜𝑛𝑠𝑡𝑎𝑛𝑡). Thus, we picked three circumferences (h=40, 50, and 60 mm) and 9 longitudinal sections (A1-A3, B1-B3, C1-C3) and plotted the radial and longitudinal deformations for six steps of loading (𝛿 = 0.25, 0.30, 0.35, 0.40, 0.45, 0.50 𝑚𝑚). 120 Figure A-4: Sections in the longitudinal and circumferential direction of the cylinder, which has been investigated more in the following figures Figure A-5: Sequence of radial deformation in the circumference 40 mm 121 Figure A-6: Sequence of radial deformation in the circumference of 50 mm Figure A-7: Sequence of radial deformation in the circumference of 60 mm 122 Figure A-5 to Figure A-7 show that the deformations in the thicker parts of the cylinder are negligible in comparison to deformations of the thinner areas. In addition, the thicker areas in their boundary with thin regions do not show significant rotations. Thus, clamped boundary conditions with the ability to move in the axial direction of the cylinder can be a reasonable assumption. Figure A-8: Radial deformation along with section A1 123 Figure A-9: Radial deformation along with section A2 Figure A-10: Radial deformation along with section A3 124 Figure A-11: Radial deformation along with section B1 Figure A-12: Radial deformation along section B2 125 Figure A-13: Radial deformation along section B3 Figure A-14: Radial deformation along section C1 126 Figure A-15: Radial deformation along with section C2 Figure A-16: Radial deformation along with section C3 Figure A-8 to Figure A-16 show the continuity in the rotation of the thin and thick areas in their mutual boundaries. The thicker areas on top and bottom of the thin areas (blue lines in Figure 127 A-8 to Figure A-16) show significant rotational flexibility, unlike the thicker areas on circumferential boundaries of thin zones (blue arcs in Figure A-5 to Figure A-7). Also, comparing the deformation (rotation and displacement) of blue lines (thicker areas on top and bottom on the buckling panel) in Figure A-9, Figure A-12, and Figure A-15 shows that the thicker the boundary, the deformation in the radial direction decreases. Figure A-17 to Figure A-28 present the longitudinal deformation of the cylinder under the uniform axial displacement. The first three figures (Figure A-17 to Figure A-19) assumed a circle with a 1 mm radius and showed the longitudinal deformation in the radial direction of these circles. There is a small difference between deformations in the axial (longitudinal) direction in different angular locations. This small dependence of the axial deformation on the circumferential location could be seen in Figure A-3(c) as well, where the contours of axial deformations are approximately parallel to each other. Figure A-17: The longitudinal displacement showed in the circumferential direction for h=40 mm. The radius of the circle here is 1 mm, and it is just to show the magnitude of changes in the longitudinal displacement 128 Figure A-18: The longitudinal displacement showed in the circumferential direction for h=50 mm. The radius of the circle here is 1 mm, and it is just to show the magnitude of changes in the longitudinal displacement Figure A-19: The longitudinal displacement showed in the circumferential direction for h=60 mm. The radius of the circle here is 1 mm, and it is just to show the magnitude of changes in the longitudinal displacement 129 Figure A-20 to Figure A-28 show the same deformation but in the longitudinal sections. From these figures, I expect to understand the relation between axial displacement and longitudinal coordinates. Figure A-20: Axial deformation in section A1 130 Figure A-21: Axial deformation in section A2 Figure A-22: Axial deformation in section A3 131 Figure A-23: Axial deformation in section B1 Figure A-24: Axial deformation in section B2 132 Figure A-25: Axial deformation in section B3 Figure A-26: Axial deformation in section C1 133 Figure A-27: Axial deformation in section C2 Figure A-28: Axial deformation in section C3 134 Figure A-20 to Figure A-28 indicates that the thinner areas show more deformations than thick areas, as expected. However, the longitudinal deformation of the thicker area is not negligible. Also, the deformation in higher levels of loads shows some sinusoidal pattern in thin areas. Figure A-29: The longitudinal strain showed in the circumferential direction for h=40 mm. The radius of the circle here is 0.01, and it is just to show the magnitude of changes in the longitudinal strain 135 Figure A-30: The longitudinal strain showed in the circumferential direction for h=50 mm. The radius of the circle here is 0.01, and it is just to show the magnitude of changes in the longitudinal strain Figure A-31: The longitudinal strain showed in the circumferential direction for h=60 mm. The radius of the circle here is 0.01, and it is just to show the magnitude of changes in the longitudinal strain 136 Figure A-32: Axial strain in section A1 Figure A-33: Axial strain in section A2 137 Figure A-34: Axial strain in section A3 Figure A-35: Axial strain in section B1 138 Figure A-36: Axial strain in section B2 Figure A-37: Axial strain in section B3 139 Figure A-38: Axial strain in section C1 Figure A-39: Axial strain in section C2 140 Figure A-40: Axial strain in section C3 Figure A-29 to Figure A-40 suggest the longitudinal strain shows approximately sinusoidal relation with both x and y coordinates. 141 BIBLIOGRAPHY 142 BIBLIOGRAPHY [1] P. M. Reis, “A Perspective on the Revival of Structural (In) Stability With Novel Opportunities for Function: From Buckliphobia to Buckliphilia,” Journal of Applied Mechanics, Transactions ASME. 2015, doi: 10.1115/1.4031456. [2] N. Hu and R. Burgueño, “Buckling-induced smart applications: Recent advances and trends,” Smart Materials and Structures. 2015, doi: 10.1088/0964-1726/24/6/063001. [3] N. Lajnef, R. Burgueño, W. Borchani, and Y. Sun, “A concept for energy harvesting from quasi-static structural deformations through axially loaded bilaterally constrained columns with multiple bifurcation points,” Smart Mater. Struct., 2014, doi: 10.1088/0964- 1726/23/5/055005. [4] J. Shim, P. Wang, and K. Bertoldi, “Harnessing instability-induced pattern transformation to design tunable phononic crystals,” Int. J. Solids Struct., 2015, doi: 10.1016/j.ijsolstr.2014.12.018. [5] M. H. Hasan, F. M. Alsaleem, and H. M. Ouakad, “Novel threshold pressure sensors based on nonlinear dynamics of MEMS resonators,” J. Micromechanics Microengineering, 2018, doi: 10.1088/1361-6439/aab515. [6] M. Alturki and R. Burgueño, “Displacement-based design procedure for ductile reinforced concrete bridge pier walls,” ACI Struct. J., 2017, doi: 10.14359/51689497. [7] S. Liu, A. I. Azad, and R. Burgueño, “Energy harvesting from quasi-static deformations via bilaterally constrained strips,” J. Intell. Mater. Syst. Struct., p. 1045389X1878645, Jul. 2018, doi: 10.1177/1045389X18786456. [8] S. Liu, A. I. Azad, and R. Burgueño, “Architected materials for tailorable shear behavior with energy dissipation,” Extrem. Mech. Lett., vol. 28, pp. 1–7, 2019, doi: 10.1016/j.eml.2019.01.010. [9] Y. Hu, P. Kim, and J. Aizenberg, “Harnessing structural instability and material instability in the hydrogel-actuated integrated responsive structures (HAIRS),” Extrem. Mech. Lett., 2017, doi: 10.1016/j.eml.2017.02.003. [10] H. Yuk and X. Zhao, “A New 3D Printing Strategy by Harnessing Deformation, Instability, and Fracture of Viscoelastic Inks,” Adv. Mater., 2018, doi: 10.1002/adma.201704028. [11] J. Shim et al., “Harnessing instabilities for design of soft reconfigurable auxetic/chiral materials,” Soft Matter, 2013, doi: 10.1039/c3sm51148k. [12] K. Liu, J. Wu, G. Paulino, H. Q.-S. reports, and undefined 2017, “Programmable deployment of tensegrity structures by stimulus-responsive polymers,” nature.com, Accessed: Nov. 16, 2021. [Online]. Available: https://www.nature.com/articles/s41598- 143 017-03412-6. [13] M. Konaković-luković et al., “Rapid Deployment of Curved Surfaces via Programmable Auxetics,” ACM Trans. Graph, vol. 37, p. 13, 2018, doi: 10.1145/3197517.3201373. [14] J. Liu, W. Liu, A. Pantula, Z. Wang, D. H. Gracias, and T. D. Nguyen, “Periodic buckling of soft 3D printed bioinspired tubes,” Extrem. Mech. Lett., vol. 30, p. 100514, Jul. 2019, doi: 10.1016/J.EML.2019.100514. [15] Z. Zhao, K. Wang, L. Zhang, L. C. Wang, W. L. Song, and D. Fang, “Stiff reconfigurable polygons for smart connecters and deployable structures,” Int. J. Mech. Sci., vol. 161–162, p. 105052, Oct. 2019, doi: 10.1016/J.IJMECSCI.2019.105052. [16] X. Tan, B. Wang, S. Chen, S. Zhu, and Y. Sun, “A novel cylindrical negative stiffness structure for shock isolation,” Compos. Struct., vol. 214, pp. 397–405, Apr. 2019, doi: 10.1016/J.COMPSTRUCT.2019.02.030. [17] R. Burgueño and N. Lajnef, “Energy harvesting devices for low frequency applications.” Google Patents, 2016. [18] T. Kobayashi, Y. Mihara, and F. Fujii, “Path-tracing analysis for post-buckling process of elastic cylindrical shells under axial compression,” Thin-Walled Struct., vol. 61, pp. 180– 187, 2012, doi: 10.1016/j.tws.2012.05.018. [19] T. Von Karman and H.-S. TSEIN, “The buckling of thin cylindrical shells under axial compression,” J. Spacecr. Rockets, vol. 40, no. 6, pp. 898–907, 2003. [20] N. Yamaki and S. Kodama, “Postbuckling behavior of circular cylindrical shells under compression,” Int. J. Non. Linear. Mech., vol. 11, no. 2, pp. 99–111, 1976, doi: 10.1016/0020-7462(76)90008-1. [21] N. Yamaki, Elastic stability of circular cylindrical shells, vol. 27. Elsevier, 1984. [22] C. R. Calladine, “Understanding imperfection-sensitivity in the buckling of thin-walled shells,” Thin-Walled Struct., 1995, doi: 10.1016/0263-8231(95)00013-4. [23] G. J. Simitses, “Buckling and Postbuckling of Imperfect Cylindrical Shells: A Review,” Appl. Mech. Rev., vol. 39, p. 1517, 1986, doi: 10.1115/1.3149506. [24] C. Y. Song, J. G. Teng, and J. M. Rotter, “Imperfection sensitivity of thin elastic cylindrical shells subject to partial axial compression,” Int. J. Solids Struct., 2004, doi: 10.1016/j.ijsolstr.2004.05.040. [25] J. W. Hutchinson, “Buckling and initial postbuckling behavior of oval cylindrical shells under axial compression,” J. Appl. Mech. Trans. ASME, 1964, doi: 10.1115/1.3601175. [26] R. C. Tennyson, “Buckling modes of circular cylindrical shells under axial compression,” AIAA J., 1969, doi: 10.2514/3.5419. 144 [27] L. Wullschleger and H. R. Meyer-Piening, “Buckling of geometrically imperfect cylindrical shells - Definition of a buckling load,” Int. J. Non. Linear. Mech., vol. 37, no. 4–5, pp. 645– 657, 2002, doi: 10.1016/S0020-7462(01)00089-0. [28] N. Hu and R. Burgueño, “Cylindrical Shells with Tunable Postbuckling Features Through Non-Uniform Patterned Thickening Patches,” Int. J. Struct. Stab. Dyn., vol. 18, no. 02, p. 1850026, 2018, doi: 10.1142/S0219455418500268. [29] D. Bushnell, “Buckling of Shells-Pitfall for Designers,” https://doi.org/10.2514/3.60058, vol. 19, no. 9, pp. 1183–1226, May 2012, doi: 10.2514/3.60058. [30] J. T. Miller, T. Su, J. Pabon, N. Wicks, K. Bertoldi, and P. M. Reis, “Buckling of a thin elastic rod inside a horizontal cylindrical constraint,” Extrem. Mech. Lett., 2015, doi: 10.1016/j.eml.2015.03.002. [31] J. Xiao, Y. Chen, X. Lu, B. Xu, X. Chen, and J. Xu, “Three dimensional buckling beam under cylindrical constraint,” Int. J. Mech. Sci., 2019, doi: 10.1016/j.ijmecsci.2018.10.041. [32] O. K. Bedair, “The Elastic Behaviour of Multi-Stiffened Plates under Uniform Compression,” Thin-Walled Struct., vol. 27, no. 4, pp. 311–335, 1997, doi: 10.1016/S0263- 8231(96)00032-8. [33] O. K. Bedair, “A contribution to the stability of stiffened plates under uniform compression,” Comput. Struct., vol. 66, no. 5, pp. 535–570, 1998, doi: 10.1016/s0045- 7949(97)00102-8. [34] P. Buermann, R. Rolfes, J. Tessmer, and M. Schagerl, “A semi-analytical model for local post-buckling analysis of stringer- and frame-stiffened cylindrical panels,” Thin-Walled Struct., vol. 44, no. 1, pp. 102–114, 2006, doi: 10.1016/j.tws.2005.08.010. [35] J. Singer, M. Baruch, and O. Harari, “On the stability of eccentrically stiffened cylindrical shells under axial compression,” Int. J. Solids Struct., 1967, doi: 10.1016/0020- 7683(67)90001-7. [36] J. Singer, J. Arbocz, and C. D. Babcock, “Buckling of imperfect stiffened cylindrical shells under axial compression,” AIAA J., 1971, doi: 10.2514/3.6125. [37] C. Lynch, A. Murphy, M. Price, and A. Gibson, “The computational post buckling analysis of fuselage stiffened panels loaded in compression,” Thin-Walled Struct., 2004, doi: 10.1016/j.tws.2004.04.002. [38] M. Ahmer Wadee and M. Farsi, “Cellular buckling in stiffened plates,” Proc. R. Soc. A Math. Phys. Eng. Sci., 2014, doi: 10.1098/rspa.2014.0094. [39] M. A. Wadee and M. Farsi, “Local-global mode interaction in stringer-stiffened plates,” Thin-Walled Struct., 2014, doi: 10.1016/j.tws.2014.09.012. [40] M. A. Wadee and M. Farsi, “Imperfection sensitivity and geometric effects in stiffened 145 plates susceptible to cellular buckling,” Structures, 2015, doi: 10.1016/j.istruc.2015.04.004. [41] V. I. Weingarten, P. Seide, and J. P. Peterson, “Buckling of thin-walled circular cylinders,” NASA SP-8007, 1968. [42] B. Wang et al., “Experimental validation of cylindrical shells under axial compression for improved knockdown factors,” Int. J. Solids Struct., 2019, doi: 10.1016/j.ijsolstr.2019.01.001. [43] H. N. R. Wagner and C. Hühne, “Robust knockdown factors for the design of cylindrical shells under axial compression: potentials, practial application and reliability analysis,” Int. J. Mech. Sci., 2018, doi: 10.1016/j.ijmecsci.2017.11.020. [44] R. Wagner, “Robust design of buckling critical thin-walled shell structures,” DLR Dtsch. Zent. fur Luft- und Raumfahrt e.V. - Forschungsberichte, 2019. [45] B. O. Almroth, A. B. Burns, and E. V. Pittner, “Design criteria for axially loaded cylindrical shells,” J. Spacecr. Rockets, 1970, doi: 10.2514/3.30025. [46] D. A. S. DASt, DASt-Richtlinie 013 Juli 1980. Beulsicherheitsnachweise für Schalen. Stahlbau-Verlagsges., 1980. [47] M. W. Hilburger, “Developing the next generation shell buckling design factors and technologies,” 2012, doi: 10.2514/6.2012-1686. [48] N. W. Gardner, M. W. Hilburger, W. T. Haynie, M. C. Lindell, and W. A. Waters, “Digital image correlation data processing and analysis techniques to enhance test data assessment and improve structural simulations,” 2018, doi: 10.2514/6.2018-1698. [49] M. W. Hilburger, “On the development of shell buckling knockdown factors for stiffened metallic launch vehicle cylinders,” 2018, doi: 10.2514/6.2018-1990. [50] M. R. Schultz et al., “Test and analysis of a buckling-critical large-scale sandwich composite cylinder,” 2018, doi: 10.2514/6.2018-1693. [51] S. G. P. Castro, R. Zimmermann, M. A. Arbelo, R. Khakimova, M. W. Hilburger, and R. Degenhardt, “Geometric imperfections and lower-bound methods used to calculate knock- down factors for axially compressed composite cylindrical shells,” Thin-Walled Struct., vol. 74, pp. 118–132, 2014. [52] M. A. Arbelo, R. Degenhardt, S. G. P. Castro, and R. Zimmermann, “Numerical characterization of imperfection sensitive composite structures,” Compos. Struct., 2014, doi: 10.1016/j.compstruct.2013.09.041. [53] P. Hao, B. Wang, K. Tian, G. Li, Y. Sun, and C. Zhou, “Fast procedure for Non-uniform optimum design of stiffened shells under buckling constraint,” Struct. Multidiscip. Optim., vol. 55, no. 4, pp. 1503–1516, Apr. 2017, doi: 10.1007/S00158-016-1590-3. 146 [54] T. V. Kármán, “Festigkeitsprobleme im Maschinenbau,” in Mechanik, 1907. [55] E. Riks, “An incremental approach to the solution of snapping and buckling problems,” Int. J. Solids Struct., vol. 15, no. 7, pp. 529–551, 1979, doi: 10.1016/0020-7683(79)90081-7. [56] N. Vasios, “Nonlinear analysis of structures,” The Arc-Length method. Harvard, 2015. [57] B. A. Memon and X. Z. Su, “Arc-length technique for nonlinear finite element analysis,” J. Zhejiang Univ. Sci., vol. 5, no. 5, pp. 618–628, 2004, doi: 10.1631/jzus.2004.0618. [58] M. A. CRISFIELD, “A FAST INCREMENTAL/ITERATIVE SOLUTION PROCEDURE THAT HANDLES ‘SNAP-THROUGH,’” in Computational Methods in Nonlinear Structural and Solid Mechanics, 1981. [59] B. W. R. Forde and S. F. Stiemer, “Improved arc length orthogonality methods for nonlinear finite element analysis,” Comput. Struct., 1987, doi: 10.1016/0045-7949(87)90078-2. [60] J. F. Jullien and M. Araar, “Towards an optimal shape of cylindrical shell structures under external pressure,” Buckling Shell Struct. Land, Sea Air, pp. 21–32, 1991. [61] X. Ning and S. Pellegrino, “International Journal of Solids and Structures Imperfection- insensitive axially loaded thin cylindrical shells,” Int. J. Solids Struct., vol. 62, pp. 39–51, 2015, doi: 10.1016/j.ijsolstr.2014.12.030. [62] X. Ning and S. Pellegrino, “Experiments on imperfection insensitive axially loaded cylindrical shells,” Int. J. Solids Struct., vol. 115, pp. 73–86, 2017. [63] K. Kumar Yadav and S. Gerasimidis Assistant Professor, “Imperfection Insensitivity of Thin Wavy Cylindrical Shells Under Axial Compression or Bending,” 2020, doi: 10.1115/1.4045741. [64] K. K. Yadav and S. Gerasimidis, “Imperfection insensitive thin cylindrical shells for next generation wind turbine towers,” J. Constr. Steel Res., vol. 172, p. 106228, Sep. 2020, doi: 10.1016/J.JCSR.2020.106228. [65] B. S. Cox, R. Groh, and A. Pirrera, “Nudging axially compressed cylindrical panels towards imperfection insensitivity,” J. Appl. Mech., pp. 1–20, 2019. [66] R. L. Lincoln, P. M. Weaver, A. Pirrera, and R. M. J. Groh, “Imperfection-Insensitive Continuous Tow-Sheared Cylinders,” Compos. Struct., p. 113445, 2020. [67] N. Hu and R. Burgueño, “Tailoring the elastic postbuckling response of cylindrical shells: A route for exploiting instabilities in materials and mechanical systems,” Extrem. Mech. Lett., vol. 4, pp. 103–110, 2015, doi: 10.1016/j.eml.2015.05.003. [68] R. Burgueño, N. Hu, A. Heeringa, and N. Lajnef, “Tailoring the elastic postbuckling response of thin-walled cylindrical composite shells under axial compression,” Thin-Walled Struct., vol. 84, pp. 14–25, 2014. 147 [69] N. Hu and R. Burgueño, “Elastic postbuckling response of axially-loaded cylindrical shells with seeded geometric imperfection design,” Thin-Walled Struct., vol. 96, pp. 256–268, 2015. [70] N. Hu and R. Burgueño, “Harnessing seeded geometric imperfection to design cylindrical shells with tunable elastic postbuckling behavior,” J. Appl. Mech., vol. 84, no. 1, p. 11003, 2017. [71] N. Hu, “Tailoring the elastic postbuckling response of thin-walled axially compressed cylindrical shells,” Michigan State, 2015. [72] J. Guo, S. Liu, and R. Burgueño, “Tailoring the Elastic Postbuckling Response of Thin- Walled Cylindrical Shells for Applications in Mechanical Devices and Adaptive Structures,” in ASME 2017 Conference on Smart Materials, Adaptive Structures and Intelligent Systems, 2017, p. V002T03A037--V002T03A037. [73] P. Toan Thang, T. T. Nguyen, and J. Lee, “Nonlinear static analysis of thin curved panels with FG coatings under combined axial compression and external pressure,” Thin-Walled Struct., vol. 107, pp. 405–414, Oct. 2016, doi: 10.1016/J.TWS.2016.06.007. [74] S. Zhu, Z. Tong, J. Sun, Q. Li, Z. Zhou, and X. Xu, “Electro-thermo-mechanical post- buckling of piezoelectric functionally graded cylindrical shells,” Appl. Math. Model., vol. 98, pp. 309–322, Oct. 2021, doi: 10.1016/J.APM.2021.05.011. [75] L. W. Zhang, Z. X. Lei, K. M. Liew, and J. L. Yu, “Large deflection geometrically nonlinear analysis of carbon nanotube-reinforced functionally graded cylindrical panels,” Comput. Methods Appl. Mech. Eng., vol. 273, pp. 1–18, May 2014, doi: 10.1016/J.CMA.2014.01.024. [76] M. Bodaghi and M. Shakeri, “An analytical approach for free vibration and transient response of functionally graded piezoelectric cylindrical panels subjected to impulsive loads,” Compos. Struct., vol. 94, no. 5, pp. 1721–1735, Apr. 2012, doi: 10.1016/J.COMPSTRUCT.2012.01.009. [77] Z. X. Lei, L. W. Zhang, K. M. Liew, and J. L. Yu, “Dynamic stability analysis of carbon nanotube-reinforced functionally graded cylindrical panels using the element-free kp-Ritz method,” Compos. Struct., vol. 113, no. 1, pp. 328–338, Jul. 2014, doi: 10.1016/J.COMPSTRUCT.2014.03.035. [78] N. Hu, S. Liu, R. B.-I. and I. in Structural, and undefined 2016, “Controlled elastic instabilities in cylindrical shells for energy harvesting devices,” taylorfrancis.com, Accessed: Jun. 15, 2022. [79] M. Nazar, A. I. King-James Egbe, A. Matin Nazar, P. Jiao, A. H. Alavi, and K.-J. I. Egbe, “Harnessing postbuckling instability of piezoelectric cylinders with corrugation for energy harvesting,” https://doi.org/10.1117/12.2581669, vol. 11588, pp. 277–284, Mar. 2021, doi: 10.1117/12.2581669. 148 [80] N. Fazlalipour, H. Showkati, and S. Eyvazinejad Fioruzsalari, “Buckling behaviour of cylindrical shells with stepwise wall thickness subjected to combined axial compression and external pressure,” Thin-Walled Struct., vol. 167, p. 108195, Oct. 2021, doi: 10.1016/J.TWS.2021.108195. [81] F. Zhou, Z. Chen, H. Fan, and S. Huang, “Analytical study on the buckling of cylindrical shells with stepwise variable thickness subjected to uniform external pressure,” http://dx.doi.org/10.1080/15376494.2015.1068401, vol. 23, no. 10, pp. 1207–1215, Oct. 2016, doi: 10.1080/15376494.2015.1068401. [82] L. Chen, J. Michael Rotter, and C. Doerich, “Buckling of cylindrical shells with stepwise variable wall thickness under uniform external pressure,” Eng. Struct., vol. 33, no. 12, pp. 3570–3578, Dec. 2011, doi: 10.1016/J.ENGSTRUCT.2011.07.021. [83] H. Babaei, Y. Kiani, and M. R. Eslami, “Application of two-steps perturbation technique to geometrically nonlinear analysis of long FGM cylindrical panels on elastic foundation under thermal load,” https://doi.org/10.1080/01495739.2017.1421054, vol. 41, no. 7, pp. 847– 865, Jul. 2018, doi: 10.1080/01495739.2017.1421054. [84] H. S. Shen, “Postbuckling analysis of axially loaded piezolaminated cylindrical panels with temperature dependent properties,” Compos. Struct., 2007, doi: 10.1016/j.compstruct.2006.02.018. [85] A. Keller, “A note on single crystals in polymers : Evidence for a folded chain configuration,” Philos. Mag., vol. 2, no. 21, pp. 1171–1175, 1957. [86] A. I. Azad and R. Burgueño, “Semi-Analytical Model to Predict the Elastic Post-Buckling Response of Axially Compressed Cylindrical Shells With Tailored Distributed Stiffness,” no. 4, doi: 10.1115/1.4051093. [87] Simulia, “Abaqus User’s Manual version 2019,” Dassault Systèmes Simulia Corp.: Providence, RI, USA. 2019. [88] N. Hu, R. Burgueño, and N. Lajnef, “Structural Optimization and Form-Finding of Cylindrical Shells for Targeted Elastic Postbuckling Response,” in Proceedings of the ASME 2014 Conference on Smart Materials, Adaptive Structures and Intelligent Systems, 2014, pp. 1–8, doi: 10.1115/SMASIS20147446. [89] Y. Forterre, “Slow, fast and furious: understanding the physics of plant movements,” J. Exp. Bot., vol. 64, no. 15, pp. 4745–4760, Nov. 2013, doi: 10.1093/JXB/ERT230. [90] M. Gomez, D. E. Moulton, and D. Vella, “Critical slowing down in purely elastic ‘snap- through’ instabilities,” Nat. Phys. |, vol. 13, 2017, doi: 10.1038/NPHYS3915. [91] I. Wolfram Research, Mathematica, Version 12. Champaign, Illinois, 2020. [92] P. Jiao, Z. Chen, X. Tang, W. Su, and J. Wu, “Design of axially loaded isotropic cylindrical shells using multiple perturbation load approach – Simulation and validation,” Thin-Walled 149 Struct., vol. 133, pp. 1–16, Dec. 2018, doi: 10.1016/J.TWS.2018.09.028. [93] B. Wang et al., “Buckling of quasi-perfect cylindrical shell under axial compression: A combined experimental and numerical investigation,” Int. J. Solids Struct., vol. 130–131, pp. 232–247, Jan. 2018, doi: 10.1016/J.IJSOLSTR.2017.09.029. [94] W. Verduyn, I. E.-D. U. of Technology, D. of, and undefined 1982, “A testing machine for statistical analysis of small imperfect shells: Part I,” narcis.nl, Accessed: Feb. 10, 2022. [95] V. I. Weingarten, E. J. Morgan, and P. Seide, “Elastic stability of thin-walled cylindrical and conical shells under axial compression,” https://doi.org/10.2514/3.2893, vol. 3, no. 3, pp. 500–505, May 2012, doi: 10.2514/3.2893. [96] R. C. Tennyson, “Buckling modes of circular cylindrical shells under axial compression.,” https://doi.org/10.2514/3.5419, vol. 7, no. 8, pp. 1481–1487, May 2012, doi: 10.2514/3.5419. [97] W. Flügge, “Die Stabilität der Kreiszylinderschale,” Ingenieur-Archiv, vol. 3, no. 5, pp. 463–506, Dec. 1932, doi: 10.1007/BF02079822. [98] Esslinger, Maria, and Bodo Geier. Gerechnete Nachbeullasten als untere Grenze der experimentellen axialen Beullasten von Kreiszylindern. Deutsche Forschungs-und Versuchsanstalt für Luft-und Raumfahrt, 1972. [99] by Jobann Arbocx, C. D. Babcock, and J. Prepared hy, “Experimental investigation of the effect of general imperfections on the buckling of cylindrical shells,” 1968, Accessed: Feb. 10, 2022. [100] C. Babcock and E. Sechler, “The effect of initial imperfections on the buckling stress of cylindrical shells,” 1963, Accessed: Feb. 10, 2022. [101] B. ALMROTH, “INFLUENCE OF IMPERFECTIONS AND EDGE RESTRAINT ON THE BUCKLING OF AXIALLY COMPRESSED CYLINDERS,” Apr. 1966, doi: 10.2514/6.1966-1702. [102] D. P. Searson, “GPTIPS 2: An open-source software platform for symbolic data mining,” Handb. Genet. Program. Appl., pp. 551–573, Jan. 2015, doi: 10.1007/978-3-319-20883- 1_22/FIGURES/10. 150