Modal Superposition Using Implicit Boundary Finite Element Method and Application to Flapping Wing Design

MISSING IMAGE

Material Information

Title:
Modal Superposition Using Implicit Boundary Finite Element Method and Application to Flapping Wing Design
Physical Description:
1 online resource (96 p.)
Language:
english
Creator:
Zhang, Zhiyuan
Publisher:
University of Florida
Place of Publication:
Gainesville, Fla.
Publication Date:

Thesis/Dissertation Information

Degree:
Master's ( M.S.)
Degree Grantor:
University of Florida
Degree Disciplines:
Mechanical Engineering, Mechanical and Aerospace Engineering
Committee Chair:
Kumar, Ashok V
Committee Members:
Haftka, Raphael Tuvia
Ifju, Peter G

Subjects

Subjects / Keywords:
base -- boundary -- excitation -- implicit -- method -- vibration
Mechanical and Aerospace Engineering -- Dissertations, Academic -- UF
Genre:
Mechanical Engineering thesis, M.S.
bibliography   ( marcgt )
theses   ( marcgt )
government publication (state, provincial, terriorial, dependent)   ( marcgt )
born-digital   ( sobekcm )
Electronic Thesis or Dissertation

Notes

Abstract:
Modal superposition is a widely used technique in engineering design for structural dynamic analysis. It uses the natural frequencies and mode shapes to characterize the dynamic response of a linear system. It is more efficient than direct integration when doing repeated response analysis and it is more accurate because it avoids integration errors. This approach has been implemented in the finite element method to model forced vibration, free vibration and base excitation. Implicit Boundary Finite Element Method (IBFEM) is a mesh independent finite element method, which avoids mesh generation difficulties by using a background mesh for analysis. The background mesh is not dependent on the geometry and the geometry of the model is exactly represented by using equations of the boundary. The nodes of the mesh do not need to be on the boundary. This approach has been applied to various types of analysis such as solid mechanics, thermal conduction etc. In this thesis, the main goal is to extend the implicit boundary method to structural dynamics using modal superposition technique. In addition to evaluating mass and stiffness using a background mesh for modal superposition, an approach for imposing base excitation is developed that allows both translational and rotational input. Modal superposition technique for free vibration, forced vibration, frequency response and base excitation is introduced and its implementation using IBFEM is discussed. Several examples of these cases are simulated using the implemented algorithms and the results are compared with commercial finite element analysis software. The main application illustrated in this thesis is the design of flapping wing micro air vehicles, where the goal is to design the wing to flap in a desired fashion that mimics bird wings. The shape of the wing can be complicated and needs repeated analysis for optimization. Two different approaches for designing the wings are introduced in this thesis.
General Note:
In the series University of Florida Digital Collections.
General Note:
Includes vita.
Bibliography:
Includes bibliographical references.
Source of Description:
Description based on online resource; title from PDF title page.
Source of Description:
This bibliographic record is available under the Creative Commons CC0 public domain dedication. The University of Florida Libraries, as creator of this bibliographic record, has waived all rights to it worldwide under copyright law, including all related and neighboring rights, to the extent allowed by law.
Statement of Responsibility:
by Zhiyuan Zhang.
Thesis:
Thesis (M.S.)--University of Florida, 2013.
Local:
Adviser: Kumar, Ashok V.
Electronic Access:
RESTRICTED TO UF STUDENTS, STAFF, FACULTY, AND ON-CAMPUS USE UNTIL 2014-05-31

Record Information

Source Institution:
UFRGP
Rights Management:
Applicable rights reserved.
Classification:
lcc - LD1780 2013
System ID:
UFE0045485:00001


This item is only available as the following downloads:


Full Text

PAGE 1

1 MODAL SUPERPOSITION USING IMPLICIT BOUNDARY FINITE ELEMENT METHOD AND APPLICATION TO FLAPPING WING DESIGN By ZHIYUAN ZHANG A THESIS PRESENTED TO THE GRADUATE SCHOOL OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF MASTER OF SCIENCE UNIVERSITY OF FLORIDA 2013

PAGE 2

2 2013 Zhiyuan Zhang

PAGE 3

3 T o my parents

PAGE 4

4 ACKNOWLEDGMENTS I would like to express my sincere thanks to my advisor, Dr. Ashok V. Kumar, for his patien t guidance and encouragement during my research process Without his assistance, i t might have been impossible for me to finish this thesis His diligence and knowledgeab ility inspired me not only during this research process but also for my future life. I would like to thank the members of my supervisory committee, Dr. Raphael T. Haftka and Dr. Peter G. Ifju It is my honor to have them in m y committee and be guided for my thesis. I am grateful for their willingness to review this thesis and provide valuable suggestions. I thank Dr. Huihua Feng for his advice and understanding during my undergraduate study in Beijing Institute of Technology, China. It would not have been possible for me to purs u e for my M.S. degree here in University of Florida, US without his encouragement. Finally, I would like to show my hearty thanks to my parents for their constant love and support.

PAGE 5

5 TABLE OF CONTENTS page ACKNOWLEDGMENTS ................................ ................................ ................................ .. 4 LIST OF TABLES ................................ ................................ ................................ ............ 7 LIST OF FIGURES ................................ ................................ ................................ .......... 8 LIST OF ABBREVIATIONS ................................ ................................ ........................... 11 ABSTRACT ................................ ................................ ................................ ................... 12 CHAPTER 1 INTRODUCTION ................................ ................................ ................................ .... 14 1.1 Overview ................................ ................................ ................................ ........... 14 1.2 Goals and Objectives ................................ ................................ ........................ 17 1.2.1 Goal ................................ ................................ ................................ ......... 17 1.2.2 Objectives ................................ ................................ ................................ 17 1.3 Outline ................................ ................................ ................................ .............. 17 2 IMPLICIT BOUNDARY FINITE ELEMENT METHOD AND STRUCTURAL DYNAMICS ................................ ................................ ................................ ............. 19 2.1 Introduction to Implicit Boundary Finite Element Method ................................ .. 19 2.2 Traditional FEM on Structural Dynamics ................................ ........................... 20 2.3 Implicit Boundary Finite Element Method on Structural Dynamics .................... 21 2.3.1 Mass Matrix ................................ ................................ ............................. 23 2.3.2 Stiffness Matrix ................................ ................................ ........................ 23 3 THEORIES FOR MODAL SUPERPOSITION IN IBFEM ................................ ........ 25 3.1 Basic Principle ................................ ................................ ................................ .. 25 3.2 Free V ibration ................................ ................................ ................................ ... 26 3.2.1 Under D a mped S ystem ................................ ................................ ........... 27 3.2.2 Critically D amp ed S ystem ................................ ................................ ....... 28 3.2.3 Over D amp ed S ystem ................................ ................................ ............. 28 3.3 Forc ed V ibration ................................ ................................ ................................ 29 3.3.1 Harmonic Load ................................ ................................ ........................ 29 3.3.1.1 Under damped system ................................ ................................ ... 31 3.3. 1.2 Critically damped system ................................ ............................... 32 3.3.1.3 Over damped system ................................ ................................ ..... 34 3.3.2 Arbitrary I mpact ................................ ................................ ....................... 35 3.3.2.1 Un der damped system ................................ ................................ ... 38 3.3.2.2 Critically damp ed system ................................ ............................... 38

PAGE 6

6 3.3.2.3 Over damp ed system ................................ ................................ ..... 39 3.4 Ba se Excitation ................................ ................................ ................................ 40 3.4.1 Primary Base Excitation ................................ ................................ .......... 41 3.4.2 Secondary Base Excitation ................................ ................................ ...... 42 3.5 Frequency Response ................................ ................................ ........................ 42 4 RESULT AND DISCUSSION ................................ ................................ .................. 44 4.1 Free Vibration ................................ ................................ ................................ ... 44 4.1.1 Two Ends Clamped Beam ................................ ................................ ...... 44 4.1.2 Cantilever Plate ................................ ................................ ....................... 45 4.1.3 Cantilever 3D Beam ................................ ................................ ................ 45 4.2 Forced Vibration ................................ ................................ ................................ 46 4.2.1 Two En ds Clamped Beam ................................ ................................ ...... 46 4.2.1.1 Harmonic impact ................................ ................................ ............ 46 4.2.1.2 Arbitrary impact ................................ ................................ .............. 47 4.2.2 Ca ntilever Plate ................................ ................................ ....................... 47 4.2.2.1 Harmonic impact ................................ ................................ ............ 47 4.2.2.2 Arbitrary impact ................................ ................................ .............. 48 4.2. 3 Curved Beam ................................ ................................ .......................... 48 4.3 Base Excitation ................................ ................................ ................................ 48 4.3.1 Two Ends Clamped Beam ................................ ................................ ...... 49 4.3.2 Cantilever Beam ................................ ................................ ...................... 49 4.4 Frequency Response ................................ ................................ ........................ 50 5 APPLICAT ION IN FLAPPING WING DESIGN ................................ ....................... 70 5.1 Frequency response based design ................................ ................................ ... 70 5.1.1 Heaving Base Excitation ................................ ................................ .......... 71 5. 1.2 Heaving and Pitching Combined Base Excitation ................................ .... 72 5.2 Dynamics response based design ................................ ................................ .... 74 5.2.1 Stiff W ing Design ................................ ................................ ..................... 77 5.2.2 Composite Wing Design ................................ ................................ .......... 77 6 CO NCLUSION ................................ ................................ ................................ ........ 87 6.1 Summary ................................ ................................ ................................ .......... 87 6.2 Scope of Future Work ................................ ................................ ....................... 88 APPENDIX A INPUT FORMAT FOR MODAL SUPERPOSITION SOLVER ................................ 89 B DATA FOR THE DESIRED MOTION OF FLAPPING WING ................................ .. 91 LIST OF REFERENCES ................................ ................................ ............................... 94 BIOGRAPHICAL SKETCH ................................ ................................ ............................ 96

PAGE 7

7 LIST OF TABLES Table page B 1 Wing tip and wrist motion data in xyz coordinate ................................ ................ 92 B 2 Wing tip and wrist motion data in x'y'z' coordinate ................................ .............. 93

PAGE 8

8 LIST OF FIGURES Figure page 4 1 Two end clamped beam ................................ ................................ .................... 51 4 2 Under damped free vibration response of thin beam ................................ ......... 51 4 3 Cantilever plate geometrical model ................................ ................................ ... 52 4 4 Cantilever plate finite element modal ................................ ................................ 52 4 5 Cantilever plate free vibration response at t=0.025s ................................ .......... 52 4 6 Cantilever plate free vibration response at t=0.05s ................................ ............ 53 4 7 Under damped free vibration of cantilever plate ................................ ................ 53 4 8 Cantilever 3D beam geometry ................................ ................................ ............ 54 4 9 Cantilever 3D beam finite element modal ................................ .......................... 54 4 10 3D beam free vibration response at t=0.025s ................................ .................... 54 4 11 3D beam free vibration response at t=0.05s ................................ ...................... 55 4 12 Under damped free vibration of 3D beam ................................ .......................... 55 4 13 Two end clamped beam harmonic forced vibration response at t=0.025s ........ 56 4 14 Two end clamped beam harmonic forced vibration response at t=0.05s .......... 56 4 15 Under damped harmonic forced vibration of t wo end clamped beam ............... 57 4 16 Load curve for arbitrary forced vibration ................................ ............................ 57 4 17 Two end clamped beam arbitrary forced vibration response at t=0.025s .......... 58 4 18 Two end clamped beam arbitrary forced vibration response at t=0.05s ............ 58 4 19 Under damped arbitrary forced vibration of t wo end clamped beam ................. 59 4 20 Cantilever plate harmonic forced vibration response at t=0.025s ...................... 59 4 21 Cantilever plate harmonic forced vibration response at t=0.05s ........................ 60 4 22 Critical damping harmonic forced vibration of cantilever plate ........................... 60 4 23 Cantilever plate arbitrary forced vibration response at t=0.025s ........................ 61

PAGE 9

9 4 24 Cantilever plate arbitrary forced vibration response at t=0.05s .......................... 61 4 25 Critical damping arbitrary forced vibration of cantilever plate ............................ 62 4 26 Geometrical model for curved beam. ................................ ................................ .. 62 4 2 7 Curved beam finite element modal ................................ ................................ .... 63 4 2 8 Under damp ed harmonic forced vibration of curved beam. ................................ 63 4 2 9 Curved beam harmonic forced vibration response at t=0.0 42 5s ....................... 64 4 30 Geometrical model for two ends clamped beam base excitation ...................... 64 4 31 Two ends clamped beam base excitation deformation plot at t=0.195s ............ 64 4 32 Two ends clamped beam base excitation response ................................ .......... 65 4 33 Geometrical model for cantilever beam base excitation ................................ .... 65 4 34 Cantilever beam base excitation deformation plot at t=0.195s .......................... 66 4 3 5 Cantilever beam base excitation response ................................ ......................... 67 4 3 6 Frequency response example geometrical modal ................................ ............. 68 4 3 7 Cantilever plate frequency response deformation plot at fre=25Hz ................... 68 4 3 8 Cantilever plate frequency response deformation plot at fre=50Hz ................... 68 4 3 9 Cantilever plate frequency response deformation plot at fre=75Hz ................... 69 5 1 Flapping wing geometry for frequency response based design ......................... 79 5 2 Desired shape for frequency response flapping wing design ............................ 7 9 5 3 Normalized % RMS error versus the forcing frequency for graphite/epoxy laminated composite wing ................................ ......................... 80 5 4 Actual deformed shape of frequency response for a Graphite/Epoxy laminated composite wing flapping at 32 Hz ................................ ...................... 80 5 5 Normalized % RMS error of heaving and pitching versus the forcing frequency for graphite/epoxy laminated composite wing ........................ 81 5 6 Flapping wing shape for dynamics response based design .............................. 81 5 7 Coordinate systems for humming bird ................................ ............................... 81

PAGE 10

10 5 8 FFT of wing tip displacement ................................ ................................ ............. 82 5 9 FFT of wrist displacement ................................ ................................ .................. 82 5 10 The actual optimal result in z' direction when only minimize the error from wingtip ................................ ................................ ................................ ............... 83 5 11 The actual optimal result in z' direction when only minimize the error from wrist ................................ ................................ ................................ ................... 83 5 12 The actual optimal result in z' direction when only minimize the error from wingtip and wrist ................................ ................................ ................................ 84 5 13 The normalized error for different plies. ................................ .............................. 84 5 14 The actual optimal result in z' direction when only minimize the error from wingtip ................................ ................................ ................................ ............... 85 5 15 The actual optimal result in z' direction when only minimize the error from wrist ................................ ................................ ................................ ................... 85 5 16 The actual optimal result in z' direction when only minimize the error from wingtip and wrist ................................ ................................ ................................ 86 B 1 Global coordinate systems for humming bird ................................ ..................... 91 B 2 Coordinate systems for humming bird ................................ ............................... 92

PAGE 11

11 LIST OF ABBREVIATIONS CAD Computer Aided Design CFD Computational Fluid Dynamics DOFs Degree of Freedoms EBC Essential Boundary Condition FEM Finite Element Method FFT Fast Fourier Transform IBFEM Implicit Boundary Finite Element Method MAVs Micro Air Vehicles NUR BS Non Uniform Rational B Splines X FEM Extended Finite Element Method

PAGE 12

12 Abstract of Thesis Presented to the Graduate School of the University of Florida in Partial Fulfillment of the Requirements for the Degree of Master of Science MODAL SUPERPOSITION USING IMPLICIT BOUNDARY FINITE ELEMENT METHOD AND APPLICATION TO FLAPPING WING DESI GN By Zhiyuan Zhang May 2013 Chair: Ashok V. Kumar Major: Mechanical Engineering Modal superposition is a widely used technique in engineering design for structural dynamic analysis It u ses the natural frequencies and mode shapes to characterize the dynamic response of a linear system It is more efficient than direct integration when doing repeat ed response analysis and it is mo re accurate because it avoids integration errors This approach has been implemented in the f inite element method to model forced vibration, free vibration and base excitation. Implicit Boundary Finite Element Method (IBFEM) is a mesh independent finite element method, which avoid s mesh generation difficulties by using a background mesh for analysis. The background mesh is not dependent on the geometry and the geometry of the model is exactly represented by using equation s of the boundary. The nodes of th e mesh do not need to be on the boundary. This approach has been applied to various types of analysis such as solid mechanics, thermal conduction etc. In this the sis, the main goal is to extend the implicit boundary method to structur al dynamics using moda l superposition technique. In addition to evaluating mass and stiffness using a background mesh for modal superposition, an approach for imposing base excitation is developed that allows both trans lational and rotational input.

PAGE 13

13 Modal superposition techniqu e for free vibration, forced vibration frequency response and base excitation is introduced and its implement ation using IBFEM is discussed Several examples of these cases are simulated using the implemented algorithms and the results are compared with commercial finite element analysis software. The main application illustrated in this thesis is the design of flapping wing micro air vehicle s where the goal is to design the wing to flap in a desired f ashion that mimics bird wings. The shape of the wing can be complicate d and need s repeat ed analysis for optimization. Two different approaches for designing the wing s are introduced in this thesis.

PAGE 14

14 CHAPTER 1 INTRODUCTION 1.1 Overview In practical problems, most structures are under dynamic loads and it is necessary to perform dynamic analysis to simulate structural response and predict failure of structures. Generally, there are two techniques for simulating structural vibration, one i s modal superposition and the other is direct integration. For modal superposition, natural frequencies and modal shapes obtained from mod a l analysis are used to compute dynamic response. The other technique compute s the structural dynamic response by directly integrating through time domain, which highly depends on the accuracy of integration. For large scale structure, direct integration for dynamic analysis is computationally expensive Often it is more efficient to do mod a l analy sis once and use the result for modal superposition to find the dynamic response. Thus, the modal superposition method is popular in linear dynamic analysis. In order to solve engineering analysis problems, different numerical methods are developed and fi nite element method (FEM) is the most well known and widely used technique FEM discretiz e s the domain and boundar y of a problem by elements, and connect these elements by nodes. After being reported since 196 0s [1] [2] FEM ha s been well established and can be easily applied to many different areas of engineering analysis Structural vibration analysis using FEM is also well established However, traditional FEM has some shortcomings The first one is high c ost and difficult in creating a good quality mes h for the analysis Since traditional FEM need s high quality of mesh to get accura te results, the time spend on mesh generation is usually longer compared with that for computation.

PAGE 15

15 Even though meshes are created using accurate geometry from CAD software they do not fully represent the accurate geometry. Often the curves and surfaces of the boundary can only be approximated by the mesh. Some small fillet s or small hole s are deleted when generating the mesh, which can lead to stress infinity at the corner of some sharp shapes. Another disadvantage is that traditional FEM is not efficient in adaptive analysis which requires remesh ing the domain Since this process is difficult to fully automate and mesh must be recreated from the geometry in traditional FE M Due to these two main disadvantages of traditional FEM, three different approaches have been studies recently to improve the efficiency The first one that can overcome the two disadvantages of traditional FEM is meshfree (meshless) method. Liu etc. [ 5 ] gives the definition of meshfree method and Y. T. G u [6] summaried different meshfree method s and provides comparisons. Moving Least Square (MLS) method [7] Element Free Galerkin Method (EFGM) [8] and Meshless Local Petrov GalerKin Method [9] are examples of meshless method s There is no mesh generation process, so it does not have remeshing problem inefficiency This method was also implemented in structural dynamics [ 20][21]. However, the application of meshless method still has some unresol ved issues Firstly, it needs convergence and consistence proofs. C urrently meshless methods are not computationally efficient because its shape functions are expensive to evaluate due to lack of connectivity between nodes Besides, it is difficult to app ly precise boundary conditions, because most methods which are used to represent trial functions do not have Kronecker delta properties. I sogeometric analysis, which is first presented by T.J.R. Hughes [3] is the second method It could be regarded as a n extension f rom traditional FEM and it is an

PAGE 16

16 analysis method based on computer aided design (CAD) model Non uniform rational B Splines (NURBS) are used to represent geometry accurately in the mesh So the mesh could be easily recreated from the previous mesh which already has all the geometry information, rather than start from the geometry again This approach has also been applied to solve dynamic problems [4]. However, because NURBS always need high order to present the boundary and get accurate res ult, computation cost is usually higher and there are still some difficult ies for isogeometrical method [4] An other alt ernative method to eliminate th ose issues is mesh independent method, which means the mesh is independent of the geometr y. Thus, the mes h can be easily created automatically and the boundaries can be represented precisely. Belytschko et c [ 10 ] has presented extended finite element method (X FEM) which use a background mesh and implicit boundary representation to solve crack growth problem. The penalty method is used to describe the boundaries. The Implicit Boundary Finite Element Meth od (IBFEM) [11] is another method us ing a uniform background mesh to decreti z e the domain. Application to static problems such as elastic problem s and heat conduction problems have been well developed, but the feasibility and accuracy of applying IBFEM on structural dynamics have not been studied. Motivated by these, we explore the use of IBFEM for structural dynamics i n this thesis Dynamic analysis is often needed to find optimal structure s or input excitation to reduce vibration. A utomobile and aircraft structures which need reduce d vibration, are two well known areas in our daily life. In these examples the structures usually have complex geometry for which it is hard to generate mesh. For some other applications we need the structure to vibrate in a desired fashion An interesting application is wing

PAGE 17

17 design for flapping wing micro air vehicles From studying insects' wing structure, we know that the sh apes of designed flapping wings can be complex. In this thesis, we us e the mesh independent FEM approach to design the wing to have a desired frequency response when base excitation is applied. We also study the possibility of designing the wing to h ave non harmonic dynamic response by applying a periodic base excitation. 1.2 Goals and Objectives 1.2.1 Goal The main goal of this thesis is to implement m odal superposition dynamic analysis using the idea of impl icit boundary finite element method (IBFEM) This mesh independ ent analysis technique is then used to design the wings of a flapping wing micro air vehicle to have a specified dynamic response 1.2.2 Objectives The main objectives of the thesis are listed below: 1. Compute mass matrices and perform modal analysis for a variety of elements including 2D, shell like and 3D structures 2. U se model superposition approach to simulate the following dynamic response cases in IBFEM: free vibration forced vibration analysis, including harmon ic and arbitrary impact translational and rotational base excitation including harmonic and arbitrary base motion. 3. Extend IBFEM to compute frequency response. 4. Flapping wings design to achieve a desired steady state response or a specified dynamics response. 1.3 Outline The remaining chapters of this thesis are organized as follows: In Chapter 2, the theory of IBFEM is explained. S tructural dynamic analysis in traditional FEM is also summarized in this chapter After that the basic equation for implementing structural dynamics in IBFEM is introduced.

PAGE 18

18 In Chapter 3, the theory of modal superposition and its implementat ion using IBFEM is derived in detail. The formulas for free vibration, forced vibration, base excitation and frequency response are given in this chapter. Both harmonic and arbitrary impacts are considered In Chapter 4, we give some examples which include two ends clamped beams, cantilever plate and cantilever 3D beams to validate the implementation All examples are simulated in free vibration, forced vibration, base excitation and frequency response. In Chapter 5, th is new IBFEM implementation is applied on flapping wing design. Two different approaches are taken to find the optimized wings O ne is to obtain a desired steady sta te response and the other is based on dynamic response in time domain. In Chapter 6, the summary of the work and conclusions are provided. The future work prospect is also given in this chapter.

PAGE 19

19 CHAPTER 2 IMPLICIT BOUNDARY FINITE ELEMENT METHOD AND STRUCTURAL DYNAMICS Implicit Boundary Finite Element Meth od (IBFEM) is a method us ing a structured mesh to discretize the domain. P revious work [1 1 ] [1 5 ] ha s demonstrated application to structural analysis, electromagnetism and heat transfer analysis areas As one of the most widely used capability the implementation o f structural dynamics has great value In this chapter, we first briefly int roduce the theory of IBFEM and structural dynamics in traditional FEM Then we derive the basic equation of i mplementing structural dynamics in IBFEM. 2.1 Introduction to Implicit Boundary Finite Element Method As described in C hapter 1, IBFEM is a mesh independent method, which uses a structur ed background mesh to discretize the domain. The boundary of the domain is represented using equations. The solution structure of IBFEM could be written as following forms. (2 1) W here is an approximate step function of the bounda ry and is essential boundary condition (EBC) is a grid variable which is defined by piecewise interpolation or approximation. So Equation ( 2 1 ) could satisfy the condition on the boundary defined by Besides, the function value must be positive and non zero inside the domain. An approximate step function of the boundary can be defined as : (2 2 )

PAGE 20

20 W here, is the normal distance between point to the boundary lines and is the width defined for step function After substitut ing the step function to the integra l form o f governing equation, remaining process is identical to traditional FEM. In the following material, we will focus on the application o f structural dynamics. 2.2 Traditional FEM on Structural Dynamics Before we derive the equation for implementing structural dynamics in IBFEM, w e briefly recall the structural dynamics in traditional FEM. In general, the equilibrium equation for linear elastic problem could be written as: ( 2 3 ) W here is the stress, is body force, is velocity, is acceleration, is density and is damping ratio. Boundary condition: 1. Essential boundary condition: on 2. Natural boundary condition: p oint force t raction on b ody force on domain For dynamic problems, all the variables are time dependent. Thus we need initial condition s to solve this problem: and According to Galerkin method, the integral form for balance equation and boundary condition could be written as: (2 4 ) Rewriting this equation we can have:

PAGE 21

21 (2 5 ) So we can have: (2 6 ) W here and After applying appropriate boundary conditions and initial conditions we can solve the problem by solving Equation (2 6 ). 2.3 Implicit Boundary Finite Element Method on Structural Dynamics Usi ng the basic theory of structural dynamics with the solution structure defined earlier we can derive the equation for implement ing s tructural d ynamics in Implicit Boundary Finite Element Method T he weak form of Equation ( 2 4 ) could be written as: (2 7 ) Where is virtual strain, is Cauchy stress tensor, is traction vector, is virtual displacement and is body force. Unlike for static analysis, for structural dynamic analysis the displacement field is also a function of time. Especially, f or base excitation problems the displacement response usually consist s of two portions one is the static response because of the essential boundary load at a specific time t, the other portion is the dynamic response of the base excitation. Thus, in order to be more general, the grid variable need to be divided to two portions and is the static response because of the

PAGE 22

22 essential boundary load at a specific time t and is the dynamic respon se grid variable. The detail implementation of base excitation will be talked in C hapter 3 We can rewrite Equation (2 1) as: ( 2 8 ) Where is the boundary value is static grid variable and is dynamic grid variable T he stress and strain could also be decomposed as: (2 9 ) Equation ( 2 8 ) could be rewritten as: ( 2 1 0 ) By s implifying the above equation and we have: ( 2 1 1 )

PAGE 23

23 E liminat ing on both side of Equation (2 1 1 ), we can have: ( 2 1 2 ) F or different problems, the specific equation form for ( 2 1 2 ) could be different. For example, is for forced vibration problem, is for free vibration and is for pure base excitation problem. In the Chapter 3 we will discuss the detail s for different problem s separately The solution of Equation (2 1 2 ) is the dynamic grid variable which is noted as 2.3.1 Mass Matrix As mentioned in section 2.2 the mass matrix term could be written as: The mass matrix for one element could be written as: (2 1 3 ) 2.3.2 Stiffness Matrix The stiffness matrix is more complex than mass matrix. The homogeneous strain components could be written as: (2 1 4 ) W e can have the stiffness matrix for one element could be written as: (2 1 5 )

PAGE 24

2 4 W here, the strain displacement matrix is:

PAGE 25

25 CHAPTER 3 THEORIES FOR MODAL SUPERPOSITION IN IBFEM 3.1 Basic Principle For structural dynamic analysis, there are many differ ent method could be used. Modal superposition is a convenient and efficient method especially for models with larger number of nodes As shown in Chapter 2 t he dynamic equilibrium equation for an n degree of freedom linear mechanical system can be written as ( 3 1 ) In order to solve these equations, a transformation equation is implemented. is n by n nonsingular matrix and is called generalized displace ment. Substitute it to Equation ( 3 1 ), we have ( 3 2 ) Where Equation ( 3 2 ) is decoupled, if matrixes and are all diagonal matri c es. We need to find a matrix that can satisfy these conditions. Recall from the natural frequency and modal shape analysis, by neglecting damp term, we have (3 3 ) To solve this Equation ( 3 3 ), we assume the general solution is Then Equation (3 3 ) becomes an eigenvalue and eigenvector prob lem. We assume that for this n degree of freedoms (DOF s ) problem the eigenvalue is and the eigenvector is For the eigenvector it is easy to be proved that

PAGE 26

26 ( 3 4 ) By defining we can rewrite Equation ( 3 4 ) as It is obvious that is also a diagonal matrix. Specially, if we normalize eignvector according to mass matrix, we can get a new set of eignvector which has I f we set Equation ( 3 2 ) can be decoupled except for the damping term. If we assume proportional damping and express the damping matrix in Rayleigh damping form as then is a diagonal matrix. Then Equation ( 3 2 ) is a set of n uncoupled ordinary differential equations. ( 3 5 ) Thus, the dynamic force equilibrium equation is decomposed to n single degree of freedom equations. D ifferent solution s of Equation ( 3 5 ) could be obtained for different excit ation force S olving each single DOF dynamic force equilibrium equation we get solutions in the modal coordinate s system. After Equation (3 5 ) is solved and we get the response in modal coordinate system t he actual solution can be written as (3 6 ) Where is the actual response is the eigenvector matrix and is the response in modal coordinate system 3.2 Free V ibration For free vibration problem, the dynamic force equilibrium equation system in modal coordinate is

PAGE 27

27 ( 3 7 ) Rewrit ing the equation we have ( 3 8 ) W here, is natural frequency and is th e damping ratio. Assuming the solution in the following form and substituting it into Equation ( 3 8 ), we get ( 3 9 ) T he roots are Based on different there are three cases: 3.2.1 U nder D amped S ystem F or under damped system s the damping ratio satisfies After setting we have The complementary solution for Equation ( 3 8 ) is ( 3 10 ) ( 3 11 ) The constant A and B can be computed by the initial condition and By set ting we have and So,

PAGE 28

28 The free vibration response in modal coordinates is (3 1 2 ) 3.2.2 Critically D amp ed S ystem F or critically damped system, the damping ratio satisfies and The complementary solution for Equation ( 3 8 ) is ( 3 13 ) ( 3 1 4 ) The constant A and B can be computed by the initial condition and T hus, and W e have The free vibration response in modal coordinates is ( 3 1 5 ) 3.2.3 Over D amp ed S ystem F or over damped system, the damping ratio satisfies and w e have The complementary solution for Equation ( 3 8 ) is ( 3 16 ) ( 3 17 ) The constant A and B can be computed by the initial condition and

PAGE 29

29 T hus, and W e have The free vibration response in modal coordinates is ( 3 18 ) 3.3 Forced V ibration For forced vibration problem, the dynamic force equilibrium equation system in modal coordinate is ( 3 19 ) Here the term includes only the forces from as in Equation ( 2 12 ). There are two types of force vibration problems: harmonic analysis and modal time history analysis. For harmonic analysis the applied load is harmonic and for modal time history a nalysis the impact could be of an arbitrary type [22] 3.3.1 H armonic Load By assuming the load to be a sinusoidal function of frequency, the h armonic analysis force equilibrium equation becomes: ( 3 2 0 ) W here is the frequency and is initial phase

PAGE 30

30 Rewrite the Equation ( 3 2 0 ) we have ( 3 21 ) W here, is natural frequency, is damping ratio and The solution of this equation has two parts, complementary solut ion and particular solution. Th ese two parts could also be called the transient part and steady state part. The total solution will be the sum of the two parts. Complementary solution The complementary solution is the same as the free vibration response. Particular solution W e assume the particular solution Then we have and S ubstitut ing them into Equation (3 21) we can get Rewrite above equation we have: ( 3 2 2 ) By solving these equations we get ( 3 2 3 ) Hence the particular solution will be

PAGE 31

31 ( 3 2 4 ) ( 3 2 5 ) T here are three cases for complementary solution b ased on different : 3.3.1.1 U nder damped s ystem The complementary solution is ( 3 26 ) The total response of the system will be (3 27 ) ( 3 28 ) Appling i nitial condition and W e have (3 29 ) So,

PAGE 32

32 (3 3 0 ) (3 3 1 ) ( 3 3 2 ) 3.3.1.2 Critically damped s ystem The complementary solution is ( 3 3 3 ) The total response of the system will be

PAGE 33

33 (3 3 4 ) (3 3 5 ) Appling i nitial condition and W e have (3 36 ) So, (3 37 ) (3 38 )

PAGE 34

34 (3 39 ) 3.3.1.3 Over damped s ystem The complementary solution is ( 3 40 ) The total response of the system will be (3 41 ) ( 3 4 2 ) Appling i nitial condition and W e have (3 43 ) So,

PAGE 35

35 (3 44 ) (3 45 ) 3.3.2 A rbitrary I mpact For an arbi trary impact dynamic force equilibrium equation in modal coordinate becomes: (3 46 )

PAGE 36

36 Rewrite the Equation ( 3 46 ) we have ( 3 47 ) W here, is natural frequency, is damping ratio and The solution to this equation also has two parts, the complementary solution and particular solution [22] Complementary solution The complementary solution is the same as the free vibration response. Particular solution Th e part icular solution for an arbitrary force is usually generated by C onvolution integ ral, which is also called as Duhamel s integral. For above Equation ( 3 4 7 ), solution using Duhamel s integral could be expressed as: (3 4 8 ) W here In some practical structural dynamics problems, the excitation function is not known in an analytical expression form, but is given by a set of discrete values. Then be approximated as a piece wise linear function F or one interval as shown above, we assume the start time of this interval is and the end time interval is So we write the excitation force as (3 4 9 )

PAGE 37

37 We apply Duhamel s integral in this time interval, ( 3 50 ) By integra ting the expression, we have ( 3 5 1 ) W here, The complementary solution in the interval could be expressed as: U nder damped ( 3 5 2 ) C ritically damp ed ( 3 5 3 ) O ver damp ed ( 3 5 4 ) The total response of the system in the interval will be (3 5 5 ) Based on different damping case, we could have the solution be written in a general form as:

PAGE 38

38 (3 5 6 ) 3.3.2.1 U nder damped s ystem 3.3.2.2 Critically damp ed s ystem

PAGE 39

39 3.3.2.3 Over damp ed s ystem

PAGE 40

40 Equation ( 3 5 6 ) is the general iterati ve equation for solv ing an arbitrary excitation problem. By implementing appropriate interval the result could be accurate. However, this method could only be used on linear dynamics problem. 3.4 Base Excitation For base excitation, the input is from a constrained essential boundary condition. This analysis is widely used in seismology I n traditional FEA the total structure displacements split into quasi static and dynamic displacements which is a easy way to calculate base excitation [23] [24] Recalling Equation (2 9), in IBFEM we write the displacement in the following form. For base excitation problem s the dynamic force equilibrium equation system in modal coordinate is ( 3 5 7 ) Here we only consider the pure base excitation analysis, which means base motion is the only applied load The term only includes the forces from as in Equation ( 2 12 ). The modal superposition procedure is the same as for forced vibration,

PAGE 41

41 so we will not talk about it in details. Since term is the only difference compared with forced vibration, we will focus on the base excitation Recalling from ( 2 12 ) we have: When a problem has EBCs oscillat ing with respect to time, we call it a base excitation problem. Generally, there are two condition s for base excitation. One is that all the applied EBC s are identical functions of time. We will refer to this case as primary base excitation. When more than one EBC is applied on the structure and they are different functions of time, then the structure is statically deformed by the applied EBCs in addit ion to the deformation caused by the inertial loads. We will call this case secondary base excitation. It is obviously that, primary base excitation is simpler. We discuss both these cases below. 3.4.1 Primary B ase E xcitation P rimary base excitation is for models where all constrained EBCs' motion s are identical In this case, the structure has a rigid body motion in addi tion to structural deformation due to the inertial forces and is identical on the entire structure Then we can have and So,

PAGE 42

42 3.4.2 Secondary B ase E x citation For s econdary base excitation the motion for different constrained EBCs can be different At any given point in time, the structural is deformed statically by the applied EBCs. W e have: Where, is the static deformation field. In order to solve fo r we have: (3 5 8 ) S o from we can solve for T hen and can be derived. Therefore, base excitation problem could be solved based on Equation (3 5 7 ) using modal superposition method as mentioned in section 3.3 3.5 Frequency Response For computing the frequency response we make the following two a ssumption s: 1. N o damping 2. A ll the applied forces and base excitations are at the same frequency T he h armonic analysis dynamic force equilibrium equation in modal coordinate could be written as: (3 5 9 ) W e assume the particular solution so that and By substitut ing them into Equation (3 5 9 ) we have: (3 60 )

PAGE 43

43 R ewriting this equation, (3 6 1 ) By solving these equations we get (3 6 2 ) so (3 6 3 ) Thus, the frequency response is at any particular frequency is

PAGE 44

44 CHAPTER 4 RESULT AND DISCUSSION In order to validate the implementation of structural dynamics theories in IBFEM some simple examples are introduced in this chapter which include two ends clamped beams, cantilever plate and cantilever three dimensions ( 3D ) beams. The two end clamped beams cantilever plate and cantilever 3D beams examples are used for validating free vibration and forced vibration For base excitation, two different kinds of beam example s are provided. Cantilever plates examples are used for computing frequency response. For each example, the same mod e l was built in a commercial FEM software ( Abaqus ) to make a comparison. 4.1 Free Vibration In this section, we will present three examples about free vibration. Example 4.1.1 : Two Ends Clamped Beam A 1 D beam clamped a t bot h ends was analyzed to determine its free vibration response. The structure is as shown in the Figure 4 1 The material parameters were and This beam was modeled using plane stress element s. T he length of the beam is 1m and the cross section is 0.025mx0.025m. The initial force P=1000Pa. For the model in IBFEM, I2DBSpline16N elements are used and the mesh density is 10X3. In Abaqus, the CPS8 elements are used and the mesh densi ty is 100X3. Under damped condition is discussed here and damping ratio is 0.02 Three points which has the coordinate of A(0.25,0), B(0.5,0), C(0.75,0) are used to shown the deformation result verse time. The c omparison of IBFEM and Abaqus displacement

PAGE 45

45 results on y direction at point A B and C are shown in Figure 4 2 It is obvious that the result from IBFEM is the same as that from Abaqus. Example 4.1.2: Cantilever Plate A c ant ilever plate is as shown in Figure 4 3 The length is 0.0508m and the width i s 0.0254m. The material parameters are and This plate is modeled by us ing Bspline64N 3D shell element in IBFEM and S8R shell element in Abaqus. As shown in Figure 4 3, the initial force for free vibration is 15000 Pa. U nder damp ed condition is discussed and damping ratio is 0.02. The s tructure global deformation at t=0.025 t=0.5 are given. Three points which has the coordinate of A( 0.0254 ,0), B( 0.0 127,0), C(0,0) a re used to show the deformation result vers us time. The results are shown from Figure 4 5 to Figure 4 7 All these figures show the z component of displacement. From Figure 4 5 we find that the maximum displacement in z direction at t=0.025s in IBFEM is 3.552e 4m and it is 3.574e 4 in Abaqus. The relative difference between these two values is only about 0.6 %. From Figure 4 6 we can find that the relative difference between these two value at t=0.05s is only about 3.5 %. Example 4.1.3: Cantilever 3D Be am A three dimensional cantilever beam is as shown in the Figure 4 8 The dimension s of the beam are shown in the Figure 4 8 The material parameters are and This beam is modeled by using Quad4N 3D element in IBFEM In addition in Abaqus C3D20 elements are used to discretize this structure As shown in Figure 4 8 the initial force for free vibration is 1 0 000 Pa.

PAGE 46

46 U nder damp ed condition structure global deformation at t=0.025 t=0.5 a re given. Three points which ha ve the coordinate of A(3.0,0.5,0.25), B(2.0, 0.50,0.25), C(1.5,0.5,0.25) are used to shown the deformation result verse time. All the results are talking about displacement in y direction. The results are shown from Figure 4 10 to Figure 4 12 It indicates that the result from Abaqus and IBFEM are very close. From the above three examples, we can conclude that the free vibration capability works well for plane stress, shell 3D and 3D solid elements. As shown in Figure 4 12, the result of IBFEM and Abaqus do not match perfectly. That is because the model in Abaqus is using a coarse mesh. 4.2 Forced Vibration In this section, we will present three examples about forced vibration. In examples 4.2.1 and 4.2.2, we give harmonic forced vibration and arbitrary forced vibration result. In example 4.2.3, we only give the harmonic forced vibration result at under damped condition. Example 4.2.1: Two Ends Clamped Beam T he same structure as in Figure 4 1 is used for forced vibration The material and geometr ic constants are all the same Under damped condition is discussed here Harmonic and arbitrary loads are applied to the beam. 4.2.1.1 Harmonic i mpact T he harmonic load P is: P=1000 sin(200*pi*t)+2000 sin(1000*pi*t)+1500 sin(400*pi*t) (Pa) U nder damp ed condition is discussed and the structure global deformation at t=0.025 s and t=0.5 s are shown Three points which has the coordinate of A( 0.25 ,0),

PAGE 47

47 B( 0. 5 ,0), C(0 .75 ,0) are used to show the deformation result vers us time. The results are shown as in Figure 4 1 3 to Figure 4 1 5 4.2. 1 .2 Arbitrary i mpact The arbitrary load curve is as in Figure 4 1 6 The analysis is only from t=0 to t=0.05s. U nder damp ed condition is discussed and the structure global deformation at t=0.025 s and t=0.5 s are shown Three points which has the coordinate of A( 0.25 ,0), B( 0. 5 ,0), C(0 .75 ,0) are used to show the deformation result vers us time. The results are shown as in Figure 4 1 7 to Figure 4 19 All the se figures show the displacement component in the load direction. From the above result s of forced vibration for this beam example, we see that the result matches well for both harmonic and arbitrary load Example 4.2.2: Cantilever Plate T he same cantile ver plate used in Example 4.1.2 is subjected to forced vibration in this example The model is as in Figure 4 3 All the material and geometry constants are as described in Example 4.1.2 Critical damping condition is discussed here Harmonic and arbitrary loa ds are applied to the cantilever plate. 4.2.2.1 Harmonic i mpact T he harmonic load P is: P= 15000 sin(200*pi*t)+ 30000 sin(1000*pi*t)+ 22500 sin(400*pi*t) (Pa) G lobal deformation at t=0.025 and t=0.5 are shown Three points which has the coordinate of A( 0.0254 ,0), B( 0.0 127,0), C(0,0) are used to show the deformation result vers us time. The results are shown as in Figure 4 2 0 to Figure 4 2 2 All the figures show the displacement on z direction.

PAGE 48

48 4.2.2.2 Arbitrary i mpact The arbitrary forced vibration for ca ntilever plate is analyzed here. T he load curve is the same as shown above in Figure 4 1 6 Critical damping condition is discussed and global deformation at t=0.025 s and t=0.5 s are shown Three points which has the coordinate of A( 0.25 ,0), B( 0. 5 ,0), C(0 .75 ,0) are used to show the deformation result versus time. All the results in Figure 4 2 3 to Figure 4 2 5 show the displacement in the z direction. It is obviously that the results from IBFEM and Abaqus are close. From the results of cantilever plate forced vibration example, we see that the result matches well for both harmonic impact and arbitrary impact. Example 4.2. 3 : Curved Beam A curved beam is presented here for harmonic forced vibration. T he geometry is in Figure 4 26. The material parameters were E = 20GPa = 0.3 and density is 7800 kg/m3 The thickness is 0.01 mm Harmonic load is applied to the beam. T he harmonic load P is: P=1000 sin(200*pi*t) (Pa) For the model in IBFEM, I2DBSpline16N elements are u sed and the mesh density is 10X 10 In Abaqus, the CPS8 elements are used The mesh models are shown in Figure 4 27. Under damped condition is discussed here and damping ratio is 0.02 The c omparison of IBFEM and Abaqus displacement results on at point A(0. 1 ,0), B(0, 0 .1 ), C( 0. 1 ,0) are shown in Figure 4 2 8 Figure 4 29 gives the displacement comparison at t = 0.0425s. 4.3 Base Excitation In this section, we will present two examples about base excitation

PAGE 49

49 Example 4. 3 1 : Two Ends Clamped Beam A beam clamped at both sides is shown in Figure 4 30 The material parameters are and This beam is modeled by plane stress element s. The length of the beam is 1m and its cross section is 0.025mx0.025m. Rayleigh damp ing : alpha = 1, beta = 0.0002 The two boundaries are moving as U 1(t)=U2(t)=0.01*sin(2*pi*50*t) T he deformation plot at t=0.195 is in Figure 4 31 Three points are used to compare with the results in Abaqus at t ~ (0, 0.2) The coordinates of the three points are A(0.25,0), B(0.5,0), C(0.75,0) The results are as shown in Figure 4 3 2 All the results presented are displacement in the direction of the base excitation. From Figure 4 32 w e see that the initial response is a little different between IBFEM and Abaqus. This is because the initial condition setting in IBFEM and Abaqus is different. In Abaqus, base displacement, velocity and acceleration are all set to be zero at t=0. However, these variables in IBFEM are based on analytical derivation. For example, if the displacement is expressed as at the displacement is The velocity expression is and the value is at The initial condition only has influence on complementary solution so the steady state response is the same in IBFEM and Abaqus. Example 4. 3 2 : Cantilever Beam A cantilever beam is shown in Figure 4 33 The material parameters are and This beam is modeled by plane stress element s. The length of the beam is 1m an d its cross section is 0.025mx0.025m.

PAGE 50

50 Rayleigh damping : alpha = 1, beta = 0.02 The boundar y is moving as U 1(t)= 0.01*sin(2*pi*50*t) Five points on the beam (A, B, C, D and E) are used to compare with the results in Abaqus at t~(0, 0.2). The comparison is in Fi gure 4 3 5 T he deformation plot at t=0.195 is shown in Figure 4 3 4 All the results presented are displacement in the direction of the applied base excitation. The initial condition applied is different in IBFEM and Abaqus as mentioned before. Therefore, t he response in IBFEM and Abaqus do not match perfect ly in the beginning of the simulation 4.4 Frequency Response A c antilever plate shown in Figure 4 3 6 is used for frequency response analysis The length is 0.0508m and the width is 0.0254m. The material parameters are and This plate is modeled by using Bspline64N 3D shell element in IBFEM and S8R shell element in Abaqus as shown in Figure 4 4. D amping is not considered here. The cantilever plate is under base excitation and the amplitude is 0.05 mm For this cantilever plate, the response at three different frequencies (25Hz, 50Hz and 75Hz) is as shown from Figure 4 3 7 to Figure 4 3 9 The result s from the same model in Abaqus are also given for comparison.

PAGE 51

51 Figure 4 1. Two end clamped beam A B C Figure 4 2 Under damped free vibration response of thin beam A) Displacement at point A B) Displacement at point B C) Displacement at point C.

PAGE 52

52 Figure 4 3 Cantilever plate geometrical model A B Figure 4 4 Cantilever plate finite element modal A) E lement model in IBFEM B) E lement modal in Abaqus A B Figure 4 5 Cantilever plate free vibration response at t=0.025s A) R esult in IBFEM B) R esult in Abaqus A B C

PAGE 53

53 A B Figure 4 6 Cantilever plate free vibration response at t=0.05s A) R esult in IBFEM B) R esult in Abaqus A B C Figure 4 7 Under damped free vibration of cantilever plate A) Displacement at point A B) Displacement at point B C) Displacement at point C

PAGE 54

54 Figure 4 8 Cantilever 3D beam geometry A B Figure 4 9 Cantilever 3D beam finite element modal A) E lement model in IBFEM B) E lement modal in Abaqus A B Figure 4 10 3D beam free vibration response at t=0.025s A) R esult in IBFEM B) R esult in Abaqus

PAGE 55

55 A B Figure 4 11 3D beam free vibration response at t=0.05s A) R esult in IBFEM B) R esult in Abaqus A B C Figure 4 12 Under damped free vibration of 3D beam A) Displacement at point A B) Displacement at point B C) Displacement at point C

PAGE 56

56 A B Figure 4 13 Two end clamped beam harmonic forced vibration response at t=0.025s A) R esult in IBFEM B) R esult in Abaqus A B Figure 4 14 Two end clamped beam harmonic forced vibration response at t=0.05s A) R esult in IBFEM B) R esult in Abaqus

PAGE 57

57 A B C Figure 4 15 Under damped harmonic forced vibration of t wo end clamped beam A) Displacement at point A B) Displacement at point B C) Displacement at point C Figure 4 16 Load curve for arbitrary forced vibration

PAGE 58

58 A B Figure 4 17 Two end clamped beam arbitrary forced vibration response at t=0.025s A) R esult in IBFEM B) R esult in Abaqus A B Figure 4 18 Two end clamped beam arbitrary forced vibration response at t=0.05s A) R esult in IBFEM B) R esult in Abaqus

PAGE 59

59 A B C Figure 4 19 Under damped arbitrary forced vibration of t wo end clamped beam A) Displacement at point A B) Displacement at point B C) Displacement at point C A B Figure 4 20 Cantilever plate harmonic forced vibration response at t=0.025s A) R esult in IBFEM B) R esult in Abaqus

PAGE 60

60 A B Figure 4 21 Cantilever plate harmonic forced vibration response at t=0.05s A) R esult in IBFEM B) R esult in Abaqus A B C Figure 4 22 Critical damping harmonic forced vibration of cantilever plate A) Displacement at point A B) Displacement at point B C) Displacement at point C

PAGE 61

61 A B Figure 4 23 Cantilever plate arbitrary forced vibration response at t=0.025s A) R esult in IBFEM B) R esult in Abaqus A B Figure 4 24 Cantilever plate arbitrary forced vibration response at t=0.05s A) R esult in IBFEM B) R esult in Abaqus

PAGE 62

62 A B C Figure 4 25 Critical damping arbitrary forced vibration of cantilever plate A) Displacement at point A B) Displacement at point B C) Displacement at point C Figure 4 26 Geometrical model for curved beam

PAGE 63

63 A B Figure 4 2 7 Curved beam finite element modal A) E lement model in IBFEM B) E lement modal in Abaqus A B C Figure 4 2 8 Under damp ed harmonic forced vibration of curved beam. A) Displacement at point A B) Displacement at point B C) Displacement at point C

PAGE 64

64 A B Figure 4 2 9 Curved beam harmonic forced vibration response at t=0.0 42 5s A) R esult in IBFEM B) R esult in Abaqus Figure 4 30 Geometrical model for two ends clamped beam base excitation A B Figure 4 31 Two ends clamped beam base excitation deformation plot at t=0.195s A) R esult in IBFEM B) R esult in Abaqus

PAGE 65

65 A B C Figure 4 32 Two ends clamped beam base excitation response A) Displacement at point A B) Displacement at point B C) Displacement at point C Figure 4 33 Geometrical model for cantilever beam base excitation

PAGE 66

66 A B Figure 4 34 Cantilever beam base excitation deformation plot at t=0.195s A) R esult in IBFEM B) R esult in Abaqus

PAGE 67

67 A B C D E Figure 4 3 5 Cantilever beam base excitation response A) Point A B) Point B C) Point C D) Point D E) Point E

PAGE 68

68 Figure 4 3 6 Frequency response example geometrical modal A B Figure 4 3 7 Cantilever plate frequency response deformation plot at fre=25Hz A) R esult in IB FEM B) R esult in Abaqus A B Figure 4 3 8 Cantilever plate frequency response deformation plot at fre=50Hz A) R esult in IBFEM B) R esult in Abaqus

PAGE 69

69 A B Figure 4 3 9 Cantilever plate frequency response deformation plot at fre=75Hz A) R esult in IBFEM B) R esult in Abaqus

PAGE 70

70 CHAPTER 5 APPLICATION IN FLAPPING WING DESIGN Micro air vehicles ( MAVs ) are widely used in environmental monitoring, surve i llance and military. Different MAVs concept s have been developed, such as fixed wing, rotary wing flapping wing and etc They are usually defined as air vehicles with a maximal dimension of 15cm and 10m/s flight speed. Becau se MAVs oper ate in low Reynolds number regi on they have some technical difficult y such as low lift to drag ratio. Some birds, bats and insects share the same dimension and fight speed as MAVs which can give us an intuition and prototype in designing One common feature of the birds and insects in natural environment is they all have flapping wings. This indicates that fla pping wing is advantage ous in low Reynolds flights. Compared with fixed wing MAVs flapping wing MAVs can be more maneuverab le Most of the past research in flapping wing analysis involve aerodynamics and computational fluid dynamics ( CFD ) analysis [ 1 6 ] [ 17 ] For the wings of insects, the dynamic response du ring flapping is determined by its distribution of mass and stiffness Thus the analysis of dynamic response is necessary for flapping wing design. In this chapter, the flapping wing design is regarded as a structural design problem which involves base e xcitation to produce the dynamic response. The first section talks about the frequency response analysis based flapping wing design, in which only the maximum displacement of the wing is considered. The second section introduces a process for the wing design which tr ies to fit the dynamic response in the time domain 5.1 Frequency response based design In this section, the shape of wing is take n from the flapping wing stud y [ 18 ]. The root of wing is clamped and the material is g raphite/ e poxy laminated

PAGE 71

71 composite It has a width at root = 80 mm, width at tip = 20 mm, length of the wing = 200 mm and total thickness = 0.25 mm. Each lamina was transversely isotropic with elastic constants defined as: Pa, Pa, Pa, Pa and density The shape of the wing is shown in Figure 5.1. The flapping motion is s implified as a harmonic base excitation on z direction which has a displacement input at the root The displacement in the z direction is selected as the parameter to measure the difference between the desired deformed shape and actual frequency response shape. Two different steps are taken in this section to fit the desired shape The desired s hape is as shown in Figure 5 2. The first step only includes h eaving base excitation. The second step includes both heaving and pitching motion at the clamped root. In both steps, the laminate orientation and base excitation frequency are two design parameters. T h e sum of the squares of difference between desired shape and actual deformed shape is objective function, which is to be minimized. The optimized laminate orientation and base excitation frequency are found according to a grid points in the design space. 5.1.1 Heaving Base Excitation When only heaving base excitation applied to the wing, the optimization problem could be described as: Minimize (5 1) Where = Actual Displacement Vector = Desired Displacement Vector and = Scaling Factor

PAGE 72

72 I n order to find the optimal value of we set to zero the derivativ e of objective function with respect to We defined a root mean square (RMS) error normalized by the absolute maximum nodal z component of displacement in The normalized error was computed for various laminate orientations of the composite wing. For each orientation, this error was calculated at 2 0 equally spaced frequencies in the range 1 0 3 5 Hz. By the procedure described above the optimal laminate orientation and the frequency of the base excitation were found Figure 5 3 shows the normalized percentage RMS error calculated for the optimal ply orientation graphite/epoxy laminated composite wing for frequency form 10 Hz to 35Hz. As shown in the Fig. 5 3 for graphite/epoxy laminated composite wing the normalized RMS error is in the range of 4 30 % dependin g on the excitation frequency. The minimum normalized RMS error of 4 03 % at a forcing frequency of 32 Hz T he normalized % RMS error is n ormalized with respect to the m aximum absolute nodal z displacemen t of the desired deformed shape 5.1.2 Heaving and Pitching Combined Base Excitation The heaving based step is only based on harmonic base excitation, that is, harmonic displacement of the clamped end In order to fit the desired shape more

PAGE 73

73 closely, another excitation is introduced as an input and both heaving and pitching motion are considered as excitation The optimization p roblem can be stated as: Minimize (5 2) Where = Actual Displacement Vector of heaving = Actual Displacement Vector of pitching = Desired Displacement Vector = Scaling Factor for actual heaving displacement and = Scaling Factor for actual pitching displacement In order t o find the optimal value of and we set to zero the derivative of objective function with respect to and and The normalized error was computed for various laminate orientations of the composite wing. For each orientation, this error was calculated at 2 0 equally spaced frequencies in the range 1 0 3 5 Hz. By applying both heaving and pitching applied on the wing, the optimal ply orientation and the frequency of heaving and pitching were found in this step The forcing frequency should be in a reasonable region, so the optimal results were found a fter explor ing the design space The best results were found for which gave the

PAGE 74

74 minimum normalized RMS error of 3.37 % at a forcing frequency of 31 Hz as shown in Figure 5 5 5.2 Dynamics response based design S ection 5.1 introduced a wing designing procedure based on frequency response. However, the wing motions are usually unsymmetrical in practical problem, which cannot be represented by frequency response. Here we discuss another design method based on the dynamic resp onse of wings in time domain. Bret W. Tobalske et al [ 19 ] studied the kinematics of hummingbird flight and used synchronized high speed video cameras to measure the three dimensional wing kinematics A new wing is design ed to mimic the wing movement of hum mingbird according to his research T he length of the wing = 47 mm, average wing chord = 12 mm, s ingle wing area = 558mm 2 and total thickness = 0.25 mm. The wing shape was modeled as shown in Figure 5 6 Two steps will be used here for designing the wing In the first step, we assume the wing is made of high stiffness material and motion of the wing could be regarded as rigid body motion. A Fourier series input is found to make the dynamics response match the measured motion of hummingbird wing. In the sec ond step, we apply two layers graphite/epoxy laminated composite material The material will be different when the angle of two layers changes. The same procedure as in step one is taken for different material The optimal angle is selected to minimize the different between the dynamic response and measured motion of hummingbird wing. Here we first apply this two steps method to 12m/s forward flying. In order to simplify the problem, two points, tip and wrist, are compared instead of all the points on the wing. For f orward flying (12m/s) t he beat frequency is 41.633Hz

PAGE 75

75 and t he stroke plane angle relative to horizontal is ( ). Figure 5 7 shows the tip and wrist motion path of the win g. For the coordinate system, x direction is horizontal and z direction is vertical. Y direction is perpendicular to x and z direction. A nother coordinate system is created to simplify the process of modeling and simulation shown as the blu e one In this coordinate system, x' direction is perpendicular to s troke plane of the wing. Y direction is the same as y direction. Z' direction is perpenducular to the plane defined by the x' and y'. All t he motion data at tip and wrist according to these two coordin ate is in appendix C The motion of the wing can be decomposed into three portions, flapping motion, pitching motion and forward backward motion. Because we select two points on the wing only the flapping motion and forward backward motion are considered The wing tip and wrist displacement in the x and z' direction is shown in Table B 2 Figure 5 8 and 5 9 show Fast Fourier Transformation ( FFT ) of the data in Table B 2 In order to minimize the error of FFT, we assume the beat frequency of the wing is approximate ly 42 Hz, which is a little larger than the actual value 41.633 Hz. Because on x' direction the motion could be regarded as rigid body motion, the following analyzes are focusing on the flapping motion and only the z' displacement is considered. From Figure 5 8 and Figure 5 9 we see that the amplitudes at 42Hz, 84Hz and 126 Hz are much higher than that in other frequencies. Thus, these three frequencies are dominant. I t is obvious that the flapping motion of wing has significant influence on flig ht. Here we firstly discuss the flapping excitation and write the input flapping angle as:

PAGE 76

76 (5 3) Here we note the frequency response result under unit flapping angle in z' direction is for 0Hz, for 42Hz for 84Hz for 126Hz at wing tip and for 0Hz, for 42Hz for 84Hz for 126Hz at the wrist These eight values are constant for a specific wing. So we can write the dynamic response under in z' direction as: (5 4) (5 5) Because we are using a linear finite element model, the desired tip and wrist displacement are reduced to one third of the real measured result we note the desired motion in z' direction at the tip is and at the wrist is The differences between actual and desired motion represent the error of the whole period. T hen, we can write the objective function as: Minimize: (5 6) Where and is the w eight of error at tip and error at wrist = Actual displacement vector in z' direction of wing tip = Desired displacement vector in z'

PAGE 77

77 direction of wing tip = Actual displacement vector in z' direction of wrist and = Desired displacement vector in z' direction of wrist In order t o find the optimal value of and we set the derivative of objective function with respect to and to zero Then, the optimal input is calculated. 5.2.1 Stiff W ing D esign In this step we assume the material of wing is stiff enough that wing motion can be treated as rigid body motion. Isotropic material with E=2.0e12Pa is selected for this modal. The eight constant s and are computed based on frequency response. We can only minim ize the error at the tip by setting =1 and =0, and only minimize the error at the wrist by setting =0 and =1. W e can equal the weight of error at tip and wrist by set ting =0.25 and =1. Figure 5 10 shows the actual optimal displacement in z' direction at wing tip and wrist when we only minimize the error from tip. The displacement at wingtip matches per fectly between actual motion and desired motion but the y do not match well at the wrist Figure 5 11 shows the result when we only minimize the error from wrist. The actual and desired motions match well at the wrist, but do not match at the tip. Figure 5 12 shows displacement result in z' direction at wing tip and wrist when we minimize the error from both wing tip and wrist. It indicates that minimizing the error from both wing tip and wrist can get the best match comparing with the other two met hods. 5.2.2 Composite Wing D esign Inspired by the result of high stiffness wing that there are some dynamic deformation of the flapping wing, we seek a composite wing whose response matches

PAGE 78

78 the displacement at the tip and wrist. The same wing shape is use d here as in stiff wing design. The g raphite/ e poxy laminated composite material is applied. Each lamina was transversely isotropic with elastic constants defined as: Pa, Pa, Pa, Pa and density The wing is made of composite with two symmetric pl ies For each ply angle the same procedure introduced in stiff wing design is a pplied and a n error between actual optimal displacement in z' direction and the desired value is computed We defined the error as shown in Equation ( 5 7 ) which is a root mean square (RMS) error normalized by the absolute maximum displacement in a period. Figure 5 1 3 shows the normalized error value s for ply angle from 0 to 9 0 degree It indicates that when consider ing errors from both wingtip and wrist to optimize the input, we can get the smallest error. The optimal ply angle is 70 degree when we minimiz e the error from wingtip and wrist ( 5 7 ) Figure 5 1 4 shows the actual optimal displacement in z' direction at wing tip and wrist when we only minimize the error from tip and the ply angle is 50 degree Figure 5 1 5 shows the result when we only minimize the error from wrist and the ply angle is 90 degree Figure 5 1 6 shows displacement result in z' direction at wing tip and wrist when we minimize the error from both wing tip and wrist at ply angle =70 We can find t hat when we only minimize the error from tip, the motion will not match at the wrist. When we can make the wrist motion match, the difference at the wing tip is large. So we need

PAGE 79

79 to sacrifice the accuracy both at the wing tip and wrist and get compromise d result as in Figure 5 16. Figure 5 1 Flapping wing geometry for frequency response based design Figure 5 2 Desired shape for frequency response flapping wing design

PAGE 80

80 Figure 5 3 Normalized % RMS error versus the forcing frequency for graphite/epoxy laminated composite wing Figure 5 4 Actual deformed shape of frequency response for a Graphite/Epoxy laminated composite wing flapping at 32 Hz

PAGE 81

81 Figure 5 5 Normalized % RMS error of heaving and pitching versus the forcing frequency for graphite/epoxy laminated composite wing Figure 5 6 Flapping wing shape for dynamics response based design A B Fig ure 5 7 Coordinate systems for humming bird A ) D orsal view B ) L ateral view y' x z x S tart point z' x' x' y

PAGE 82

82 A B Fig ure 5 8 FFT of wing tip displacement A) x' directio n B) z' direction A B Fig ure 5 9 FFT of wrist displacement A) x' directio n B) z' direction

PAGE 83

83 A B Fig ure 5 10 The actual optimal result in z' direction when only minimize the error from wingtip A) Displacement at wingtip B) Displacement at wrist A B Fig ure 5 11 The actual optimal result in z' direction when only minimize the error from wrist A) Displacement at wingtip B) Displacement at wrist

PAGE 84

84 A B Fig ure 5 12 The actual optimal result in z' direction when only minimize the error from wingtip and wrist A) Displacement at wingtip B) Displacement at wrist Fig ure 5 13 The normalized error for different pl ies

PAGE 85

85 A B Fig ure 5 14 The actual optimal result in z' direction when only minimize the error from wingtip A) Displacement at wingtip B) Displacement at wrist A B Fig ure 5 15 The actual optimal result in z' direction when only minimize the error from wrist A) Displacement at wingtip B) Displacement at wrist

PAGE 86

86 A B Fig ure 5 16 The actual optimal result in z' direction when only minimize the error from wingtip and wrist A) Displacement at wingtip B) Displacement at wrist

PAGE 87

87 CHAPTER 6 CONCLUSION 6.1 Summary In this thesis, the modal superposition method was implemented using IBFEM. Initially, basic equations of modal superposition for traditional finite element and brief theory of IBFEM were introduced. Based on these the theory of modal superposition in IBFEM was presented. The main feature of this appro ach is that it combined the advantage of IBFEM and modal superposition together. A uniform structured background mesh was used for the analysis and this mesh can be generated automatically and easily Furthermore, based on uniform structured mesh integrati on is more accurate and modal superposition avoided the error of direct integration This modal superposition approach included t he analysis for free vibration, forced vibration, base excitation and frequency response It could deal with b oth harmonic and arbitrary impacts and the input could be acceleration, velocity or displacement for base excitation Some structures which include two ends clamped beams, cantilever plate and cantilever 3D beams were modeled in IBFEM. The results of these structures in f ree vibration, forced vibration, base excitation and frequency response were compared with commercial software packages. From the results, we can conclude that the implementation performs well in simulating dynamic response As an application of the dynami c analysis in IBFEM, wing design for flapping wing MAVs was introduced. Two approaches in flapping wing designing were presented. The first one used base excited frequency response analysis and assumed the wing was under a heaving motion and a pitching mot ion. The best combination of heaving and pitching were computed The other approach used total dynamic response

PAGE 88

88 in time domain and the input was regarded as a Fourier expansion The motion at the tip and wrist of the wing at flying speed of 12m/s was used as the desired motion. The optimal material for the wing and the amplitude of input Fourier expansion were generated. We can find that the motion of hummingbird's wi ng is mostly rigid body motion. In summary we have successfully implemented modal su perposition usi ng IBFEM and using this approach in flapping wing design. 6.2 Scope of Future Work We have given some examples to validate the implementation of modal superposition in IBFEM However, these examples only include plane stress, 3D stress, and 3D shell elements so they do not cover all element types. Some other element types such as beam element and Mindlin plate element also need to be validate d In the second approach for flapping wing design, the desired motion is described using only two points. It cannot present the pitching motion sufficiently Thus, the motion at another point needs to be generated as a desired motion to present the motion of the wing. Besides, the analysis in t his thesis only involve s the flight speed e qual to 12m/s. Some other flight speeds could also be analyzed.

PAGE 89

89 APPENDIX A INPUT FORMAT FOR MODAL SUPERPOSITION SOLVER # Solver Name FESmodalSuperposition # density numModes zetai prevSolverNum (density) (numofmodes) (damping ratio) (prevSolverNum)* # Harmonic load = # Number of harmonic load set ( nHL ) // numofharmonicloadset for (j=0; j< nHL ; j++) { # NumberOfFrequency=nF, BaseExcitationFlag(1=displ; 2=velocity; 3=accel) # B aseExcitationFlag not needed for forced vibration ( nF ) ( BaseExcitationFlag ) # LoadTypeNum ** (load type num 1) (load type num _m ) for (k=0; k< nF ; k++) { # Freq=Omega Amp=A Phase=phi (freq) (amplitude) < double> (phase) } } # Number of arbitrary load set (Load versus time curve) (numofarbitraryloadset)

PAGE 90

90 for (j=0; j< numofarbitraryloadset; j++) { # numOfSegmentsInLoadCurve, BaseExcitationFlag(1=displ; 2=velocity; 3=accel) # BaseExcitationFlag not needed for forced vibration ( numSegments ) ( BaseExcitationFlag ) # LoadTypeNums ** (load type num 1) (load type num _m ) for (k=0; k< numOfSegmentsInLoadCurve ;k++ ) { # NumOfPtsInSegment, Inte Flag(1=linear; 2=Cubic) ( numPtsInSegment ) (flag) for ( l =0; l < numPtsOnLdCurve ; l ++) { # Time Amp (time) (amplitude) } } } # Number of steady load (numofsteadyload) for (j=0; j< numofsteadyload; j++) { # LoadTypeNum ** load type }

PAGE 91

91 APPENDIX B DATA FOR THE DESIRED MOTION OF FLAPPING WING The data for the desired motion of flapping wing is from the work of Bret W. Tobalske E tc [ 19 ] Figure B 1 shows the tip and wrist motion path of the wing. For the coordinate system, x direction is horizontal and z direction is vertical. Y direction is perpendicular to x and z direction. We call this coordinate system global coordinate system an d t he coordinate system is as shown in following figure as the red coordinate system. The tip and wrist motion data is in Table B 1 Reverse method is us ed to get this data for Figure B 1 The maximum displacement is specified and the values on other point s are computed by scale from the plot. A B Figure B 1 Global coordinate systems for humming bird A ) D orsal v iew B ) L ateral v iew ( B 1) A nother coordinate system is created to simplify the process of modeling and simulation as the blue one shown in Figure B 2 In this coordinate system, x' direction is perpendicular to s troke plane of the wing. Y direction is the same as y direction. Z' direction is perpenducular to the plane which x' and y' included. Table B 2 shows the tip x y z x S tart point

PAGE 92

92 and wrist displacement in x'y'z' coordinate system. The transformation is presented as Equation B 1 and is t he stroke plane angle relative to horizontal direction. The following analysis is based on the motion data according to x'y'z' coordinate. A B Fig ure B 2 Coordinate systems for humming bird A ) D orsal view B ) L ateral view Tab le B 1 Wing tip and wrist motion data in xyz coordinate N Wing tip /mm Wing wrist /mm x y z x y z 1 14.2181 24.8226 35.428 1.45373 15.3318 16.3638 2 7.35045 34.1158 29.8209 3.96896 21.4613 15.0554 3 3.47911 45.3862 16.7376 5.61441 26.6022 7.20546 4 3.099185 49.3407 1.76595 8.805935 27.7885 2.13976 5 12.5779 43.4089 20.4564 14.1275 24.6249 11.8588 6 18.09295 34.1158 31.2968 14.22335 20.6704 14.4754 7 19.2538 27.3931 35.4087 15.9647 20.0772 15.2231 8 19.93185 29.3703 32.6051 16.6427 20.2749 15.2231 9 14.6103 33.3249 30.923 11.1277 22.8454 12.0457 10 3.77021 43.8044 18.5873 2.325265 26.4044 6.06474 11 23.0205 44.793 2.90666 5.41391 25.2181 4.96261 12 28.9225 35.3021 20.1019 9.09001 20.2749 13.7471 13 27.6658 26.7999 32.4375 2.80281 15.134 16.5507 14 14.2181 24.8226 35.428 1.45373 15.3318 16.3638 y y' x z x S tart point z' x' x' y

PAGE 93

93 Tab le B 2 Wing tip and wrist motion data in x'y'z' coordinate N Wing tip /mm Wing wrist /mm x y z x y z 1 0.581788 24.8226 38.17013 7.666156 15.3318 14.52988 2 4.748181 34.1158 30.34419 9.480111 21.4613 12.35092 3 3.261246 45.3862 16.78141 7.963264 26.6022 4.475118 4 2.175653 49.3407 2.82667 7.294312 27.7885 5.37737 5 3.692779 43.4089 23.7283 8.445328 24.6249 16.398 6 4.588779 34.1158 35.8579 7.522281 20.6704 18.8482 7 4.069939 27.3931 40.0989 8.83925 20.0772 20.2109 8 5.779015 29.3703 37.7753 9.464548 20.2749 20.473 9 1.521333 33.3249 34.1669 5.606461 22.8454 15.4108 10 10.6621 43.8044 15.6851 0.19981 26.4044 6.49215 11 20.1075 44.793 11.5793 3.07478 25.2181 6.669606 12 18.9039 35.3021 29.71934 3.06949 20.2749 16.19226 13 12.9766 26.7999 40.6103 3.81273 15.134 16.34761 14 0.581788 24.8226 38.17013 7.666156 15.3318 14.52988

PAGE 94

94 LIST OF REFERENCES [1] M J Turner R W Clough, H C Martin L C. Topp Stiffness and deflection analysis of Complex structures J. Aero. Sci. 23 ( 1956 ) 805 823 [2] J. H. Argyris D. W. Scharpf Finite elements in time and space Nucl. Eng. Des. 10 (1969) 456 464. [3] T.J.R. Hughes, J.A. Cottrell, Y. Bazilevs, Isogeometric analysis: CAD, finite elements, NURBS, exact geometry and mesh refinement, Comput. Methods Appl. Mech. Engrg. 194 (2005) 4135 419 5. [4] J.A. Cottrell, A. Reali, Y. Bazilevs, T.J.R. Hughes Isogeometric ana lysis of structural vibrations Comput. Methods Appl. Mech. Engrg. 195 (2006) 5257 5296. [5] G.R. Liu, Mesh free methods: moving beyond the finite element me thod CRC Press, Boca Raton 2003. [6] Y T Gu,. Meshfre e methods and their comparisons Int J Comput Methods 2 (2005) 477 515. [7] P. Lancaster K. Salkauskas Surfaces generated by moving least squares methods Math C omput 37 (1981) 141 158. [8] Y.Y. Lu, T. Belytschko L. Gu A new implementation of the element free Galerkin method Comput. Methods Appl. Mech. Engrg. 113 (1994) 397 414. [9] Z. D. Han S. N. Atluri The m eshless l ocal P etrov Galerkin (MLPG) a pproach for 3 d imensional e lasto dynamics CMC C omput. Mater. Con. 1 (2004) 129 140. [10] T Belytschko C Parimi N Mo es N. Sukumar S Usui Structured extended finite element methods for soli ds defined by implicit surfaces Int J Numer. Methods Engrg 56 (2003) 609 635. [11] A.V. Kumar, R Buria S Padmanabhan L X. Gu Finite element ana lysis using nonconforming mesh J. Comput Inf Sci Eng rg. 8 (2008) 031005 [12] A.V. Kumar S Padmanabhan, R K. Burla Implicit boundary method for finite element analysis using non conforming mesh or grid Int J Numer. Methods Engrg 74 (2008) 1421 1447 [13] R .K. Burla A V. Kumar Implicit boundary method for analysis using uniform B s pline basis and structured grid Int J Numer. Methods Engrg 76 (2008) 1993 2028. [14] R .K. Burla A V. Kumar, B V. Sankar. Implicit boundary method for determination of effective propertie s of composite microstructures Int J Solids Struct 46 (2009) 2514 2526.

PAGE 95

95 [15] S.U. Zhang, A V. Kumar Magnetostatic a nalysis u sing Implicit Boundary Finite Element Method I EEE Trans. Magn 46 (2010) 1159 1166. [16] W Shyy M Berg D. Ljungqvist Flapping and flexible wings for biological and micro air vehicles P rog. A erosp. S ci. 35 (1999) 455 505. [17] W Shyy, Y Lian, J Tang, D Viieru, H Liu Aerodynamic s of low Reynolds number flyers Cambridge University Press, New York, 2007. [18] B St anfor d P Beran M Kobayashi Aeroelastic o ptimization of o lapping w ing v enation: a c ellular d ivision a pproach in: 52nd AIAA/ASME/ASCE/AHS/ASC Structures, Structural Dynamics and Materials Conference, Denver, CO, 2011, AIAA 2011 2094 [19] B. W. Tobalske D. R. Warrick C J. Clark et D R. Powers T L. Hedrick G A. Hyder A A. Biewener al Three dimensional ki nematics of hummingbird flight J. Exp Biol 210 (2007) 2368 2382. [20] G.R. L iu, X.L. C hen, A mesh free method for static and free vibration analyses of thin plates of complicated shape J of S ound V ibrat 241 (2001) 839 855. [21] Y T Gu, G R Liu A meshless local Petrov Galerkin (MLPG) method for free and forced vibration analyses for solids Comput. M ech 27 (2001) 188 198. [22] R R Crai g A J Kurdil a Fundamentals of structural dynamic s John Wiley & Sons New York, 2006. [23] A K. Chopra, Dynamics of structures Prentice Hall, Englewo od Cliffs, NJ 1995. [24] M. Petyt Introduction to finite element vibration analysis Cambridge U niversity P ress, New York, 2010.

PAGE 96

96 BIOGRAPHICAL SKETCH Zhiyuan Zhang was born and brought up in a small city called Changzhi in Shanxi Province China. He entered Beijing Institute of Technology, China in 2006 and graduated with a Bachelor of Science in T hermal E nergy and D ynamic E ngineering in June 2010. After graduation, he joined the maste r 's program with honors in p ower m achinery and e ngineering t hermal p hysics at Beijing Institute of Technology in September 2010. He was specialized in noise, vibration and harshness an alysis of automobile. In 2011, he quit the study in China and enrolled in the master's program in Department of Mechanical and Aerospace Engineering at University of Florida. His areas of specialization include Finite Element Method, Boundary Element Method, numerical methods, optimization design and software development using o bject o rientation p rogramming t echniques.