Elastic wing response ’ s to an incoming gust

The behavior of thin elastic blade and wing subjected to a traveling disturbance is considered. The blade response to an incoming gust is predicted, then the pressure around the blade is coupled to the far field pressure in order to predict the intensity of acoustic radiation as well as the acoustic wave propagation in far field. The effect of the elasticity of the blade on the acoustic wave is predicted. The blade vibration induced by landing acoustic wave is investigated. The two dimensions inviscid flow aerodynamic theorem associated with the strip theorem are used to model the flow around the elastic thin wing. Bernoulli-Euler theorem are used in order to describe the wing motion. The fluid and the wing motions are coupled via the boundaries condition at the blade surface. The incoming gust considered here is a monochromatic wave traveling with a given speed. The problem formulation leads to a forced well known aeroelasticity Fung equation. The eigenvalue of the homogeneous part are computed and a formal solution of the forced equation is obtained


INTRODUCTION
The aerodynamic response of a rigid wing to an incoming gust dates back to the work of Possio (1938) and Sear (1941) where a transfer function is discovered establishing relationship between the amplitude of the incoming gust and the resulting force applied to the wing.Latter, among others, Goldstein and Atassi (1976), Atassi (1984) in a series of papers shown, to a second order terms accuracy in their expansion, that the unsteady lift caused by the gust can be constructed by linear superposition to the Sear lift of three independent components accounting separately for the effects of wing thickness, wing camber and the angle of attack.The deformation of the vortical part of the incoming gust by the wing generates an acoustic sources which has been modeled, among others, by Amiet (1976), Howe (1978) and experimentally explored by Arbey and Bataille (1983) and Roger and Moreau (2004).In practice applications, the wing is composed of an elastic material which may be strained under the action of the external forces.The induced deformation changes the incidence of the wing and consequently the lift.This induces a coupling between the wing motion and the external flow as well as the acoustic radiation inherent to this mechanism.The aim of the present report is to examine such an interaction between the flow, the acoustic waves and the deformation of the wing.
Interaction between the wing and the surrounded flow is, theoretically and experientially, investigated by Tang and Dowel (2001) where the wing is modeled by nonlinear beam theorem and fluid flow is described by a model driven from the flow around hard wing, It is found that the wing exhibit a limitcycle oscillation when the flow speed exceeds some critical value.Watanabe et al. (2002) resolved Navier stokes equation as well as potential flow to clarify the phenomenon of paper flutter, they conclude that the potential flow is very convenient for a parametric study of the paper flutter, Campost et al. (1998) study the correlation between the acoustic pressure and the the deformation of near by panel, they show that the acoustic pressure may cause a structure fatigue.Patil and Hodges (2001) studied the limit cycle oscillation induced by the coupling between the flow and the wing, Sun and Kaji (2002) examined the possibility of controlling the blade flutter by use of no rigid wall, Ballhaus and Goorjian (1978) Béeard et al. (2002), Brar et al. (1996), Rualli and Maute (2004) and Jiang and Wong (1998) used numerical method to study the flutter of an airfoil by indicial and numerical method.
A full numerical simulations of the aeroelaticity phenomenon in turbo machinery are undertaking by Jacquet-Richardet and Rieutord (1998), Willcox et al. (2002), Cinnella et al. (2004) and Copeland and Rey (2004), the list of course is not exhaustive.The bridges flutter phenomenon is studied, among others, by Scanlan and Jones (1999) and Salvatori and Spinelli (2006).The topic that we are developing in this report may be extended to be applied, among others fields, to turbo machinery and bridges stabilities.However, in this work we focus only on the wing flutter and the associated acoustic field and specially on the explicit derivation of the response function relating the the acoustic wave and the flutter modes to the incoming gust.To the best knowledges of the author such a transfer function does not exist in the literature.
The knowledge of the response of a streamlined body to an induced disturbance is a necessary step in the design of aeronautics and airspace machines.In order to predict the behavior of the system, constituted by the aerodynamic objects and the flow that it generates, subjected to an exterior disturbances, it is necessary to examine the interaction between the flow and the streamlined bodies under consideration.A careful examination of the behavior of the system requires solving the fluid equations and the solid equations while observing the continuity of the physical parameters, i.e. velocity and stress, at the interface fluid/solid.These procedure requires to update the domain of the fluid and the solid at each time step, if a numerical procedure is used for instant.Defining the fluid domain and the solid domain at each time step is essential when the deformation is large in comparison to the characteristic length of the system.Fortunately, in aeronautics and aerospace design often the deformation of the system has to be limited to small disturbances and small displacements to avoid failure.Therefore, the linearized equations lead often to a satisfactory results.
The set of equations suitable for the elasticity problem depend on flow Mach number, the rheological proprieties of the material used in the design and relative magnitude of the deformation of the solid under consideration.For large Mach number, the suitable equation for the description of the flow field are Boltzmann equation.In fact, the existence of shockwave and the dissociation of the gaze molecules associated with it creates a medium out of thermodynamic equilibrium and the Navier-Stokes equations fails to describe the fluid motion.For Mach number regime included in the range, to say 0.7 ≤ M ≤ 2, the compressible Navier-Stokes equations are convenient for the description of the fluid motion but a careful treatment of the shock-wave is required.For low Mach number, incompressible Navier-Stocks equations and the continuity equation approximate well the flow field generated by the aerodynamic tools.
The equation describing the solid deformation depend on the nature of the material under consideration and the magnitude of the deformation.Generally speaking, the equations of motion contain two sort of non linear terms.The first kind of non linear behavior is due to an eventual large displacement of the solid particles during their motion which is known as geometric nonlinearity.The second kind of nonlinearity is due to the constitutive equation describing stress-strain relationship.In addition, such a stress-strain relation might exhibit spatial dependence for composite material for instance.Fortunately, for the aeronautics and aerospace application, the linearized stress-strain relation and the small deformation constitute, often, a good basis for exploring and describing fluid flow and deformable solid structure interaction phenomenon.
To solve the coupled fluid and solid equations, while taking into account to the moving boundary conditions located at the interface formed by the contact between the fluid and the solid, is a hard task.In fact, the displacement of the domain demand a moving grid strategy.That is, the grid has to be updated at each time step, which requires data interpolation at each time step.Besides, the interpolation technique may lead to the instability of the numerical scheme.
In this work, we limit our self to potential flow with singularities distributed continuously on the blade camber to describe the fluid motion.The solid motion is described using Bernoulli-Euler approximation.Thus, we consider the blade as cantilever beam having two degree of freedom.These approximations allow us to find a solution of the six order partial differential forced equation in closed analytical form.The obtained solution gives detailed information about the blade response to an incoming gust.Specially, it possible to clarify how the whole range of the proper modes of the coupled blade-flow system are excited by the incoming gust and quantify the amplification rate or attenuation rate of each frequency after being excited.The obtained solution allows us to compute the pressure around the wing.Then, the computed pressure is used as a boundary condition for the far field pressure equation.Solving the far field pressure equation, one can see the effect of the elasticity of the wing on the acoustic wave and vice-versa.

STEADY AERODYNAMIC LOAD
We summarize in this section the very known elementary results concerning the aerodynamics of loaded thin wing of infinite span wise length as it has been drawn up by Katz and Plotkin.The frame of reference and the thin wing are sketched in figure (1).In this theorem the velocity potential and the velocity around the wing are The x and z components of the velocity and the pressure at the surface of the wing read (1) where γ stands for a density of vorticity distributed along the cord.The lift dL acting onto the length dx per unity of length in y direction of the wing is  1) and recall the integral of Glauert to show that at the surface of wing (7) Equation ( 4) then implies that (8) The orthogonality of the trigonometric functions leads to (9) The force applied by the fluid onto the wing, i.e. lift of the wing, and the aerodynamic moment at the leading edge (point o in figure (1)) are ( 10) (11) The point of application of the aerodynamic forces is given by the equation M(x cp ) = M(0) + x cp L = 0 Using Eq. ( 10)-( 11) we eliminate L and M(0), knowing that M(x cp ) = 0, gives (12) The moment relative to the elastic axis is (13) Figure 1 A sketch of the blade in (x, z)-plane

UNSTEADY AERODYNAMIC LOAD
In this section the wing is allowed to have two degrees of freedom, namely h and α, a translation and rotation around the elastic axis.In a frame attached to the wing, figure (1), h measured in the opposed direction of z and the rotation measured by α is about the elastic axis, α is oriented by y axis.h and α depend on y.As the system of differential equations is linear, it is enough that we consider only one Fourier component of the incoming gust.Thus, (14) A is the ratio of the amplitude of the incoming gust to the main flow, ω is the frequency of the incoming gust.The fluid velocity at the wing is now (15) as the problem is linear, we can subtracts the steady part of the solution, considering by that the undisturbed equilibrium dynamic position as a state of reference, and limit our self to (16) At the surface of the wing x = c 2 (1 -cos θ), therefore, the precedent equation becomes (17) where k is the reduced frequency, k = ωc 2U .Taking account to the well known following relation (18) the velocity at the wing becomes (19) Therefore, the identification of the coefficients of the trigonometric functions in Eq.( 7) and in Eq. ( 19), the coefficients A n could be found which in unsteady case read

LIFT AND MOMENT IN THE UNSTEADY CASE
The lift of the wing and the aerodynamic moment at the leading edge (point o) in the unsteady case can be obtained in a similar way that has been done in the steady case.Thus We consider here a cantilever.This is a beam clamped in one side and free to move in the other side.Thus, the boundary conditions are (31) (32)

DIMENSIONLESS EQUATIONS
Let the distance be scaled by c, the time by c U , the mass by ρc 2 ×1m, then EI and GJ will be scaled by 1 and so on.Then eliminating the star for brevity, the dynamic equations become where we have supposed that the cantilever wing is at undisturbed position at t = 0. Therefore h and β and their derivatives are null at t = 0. Performing Laplace transform, the equations of motion become (36) (37)

PARTICULAR NON HOMOGENEOUS SOLUTION
The forcing terms in Eq. ( 36)-(37) do not depend on y.Therefore, a particular non homogeneous solution, to say ( hp' β ~p) , may be sought as a function independent of y, thus (38) (39) Let ∆ p be the determinant of the system which reads (40) A careful inspection of the precedent equation shows that ∆p is a polynomial of degree four in s.Therefore, it posses four roots each of them is a pole for the general solution in Laplace plane.It is obvious that those roots do not depend on the elasticity modulus, rather, they depend on the mass and the inertia of the wing and the location of its elastic axis and inertia axis.For the compactness of writing of the particular non homogeneous solution, we note using the precedent notation, the particular solution could be cost in the form where σ is a real number larger than any real part of the eigenvalues s. h and β are given by Eq. ( 61)-(62), respectively.From Eq. ( 61)-(62), we can see that the poles of the function h and β are the forced reduced frequency 2ik, the zeros of ∆ p and the zeros of Det(M).We suppose in the following that there is no resonance phenomena which requires a spacial treatment which we delay for on other paper.In other words, the zero of Det(M) and the zeros of ∆ p differ from the forced frequency.Thus, the solution is composed of two parts, the first part is due to the fact that the system is directly effected by the incoming gust (Forcing term).This part of the solution represented by the pole located at s = 2ik.The other part of the solution is due to the eigenvalues of the system excited by the incoming gust which is represented by the zeros of Det(M) and the zeros of ∆ p .So, the solution can be written in the form (64) (65) In the precedent equations s = 2ik in the terms between two brackets and in the factor ∆ -1 p .So, Det(B j ), Det(M), a 1 , c 1 , a 2 , c 2 and k j have to be computed taking into account to the fact that s = 2ik.The sum in precedent equations extend over all the eigenvalues of the system except the forced frequency s = 2ik.s m , in the precedent equations, stand for the eigenvalues of the homogeneous system and the zeros of ∆ p .r m is the radius of a small circle enclosing the eigenvalue s m .The radius of the circle, r m , must be small enough to ensure the convergence of Laurent series and to not enclose in the circle more then one eigenvalue.The result of the integral must be independent of the value of r m as long as r m is smaller than the radius of the convergence of Laurent series and the circle does not enclose more than one eigenvalue.

PRESSURE FLUCTUATION
The pressure at the surface of the wing is described by Eq. (1) and Eq. ( 6).Thus, the difference in dimensionless pressure between the two surfaces of the wing is (66) where, A 0 , A 1 , ... , A n are given by Eq. ( 20)- (22).Applying Bernoulli equation and taking account to the fact that u << U, the pressure around the wing is where H n (1) , H (2)  n are Hankel functions.Applying Bernoulli equation, the pressure around the wing is given by Eq.( 68) and the pressure in far field is given by Eq. ( 71).Matching the two equations at the circumference of a circle of radius rc gives the unknown coefficients in Eq. (71).where is the temporal eigenvalue of the wave equation.By identification, has to match the argument of the exponential terms in Eq.( 64)-( 65) .The general solution is the sum of all the contribution coming from all the eigenvalues.Therefore, where, 0 = 2k, j = -is j .k is the reduced frequency of the incoming gust and s j are the eigenvalue of the homogeneous system and the zeros of ∆ p .The coefficients (A j 0 , A j 1 , A j 2 ) are the coefficients of the exponential terms having the time as an arguments in the expansion of (A 0 , A 1 , A 2 ) given by Eq. ( 20)-( 22) after elimination of h and β by Eq. ( 64)-( 65).
The asymptotic analysis of the functions H (1)  1 , H (1)  2 and H (2)  1 , H (2)  2 shows that H (1)  1 and H 2 (1) represent an incoming acoustic wave while H 1 (2) and H 2 (2) represent a radiating out acoustic wave.If there is no out acoustic source, the terms containing H 1 (1) and H 2 (1) must be dropped and the solution becomes (77) In some circumstance, for instance one it is to know if the acoustic wave may cause a material fatigue, it suitable to quantify the effect of an acoustic wave due to a source located at some distance from the blade.To do so, we neglect the effect of the blade radiation, then the solution becomes (78) In this case, φ(r, ϕ, t) is known and the coefficients A 0 j , A 1 j , A 2 j have to be computed in order to estimate the force applied by the acoustic wave on the blade.

RESULTS
The unsteady two dimension aerodynamic theorem associated with strip theorem are used to describe the load on an elastic blade or wing.Partial differential equations, describing the interaction between the flow and the elastic wing subjected to an incoming gust are established, Eq. ( 36)-(37).The flow is a potential flow with singularities distributed continuously on the camber of the wing.The wing is considered as a cantilever obeying to Euler-Bernoulli approximation.The incoming gust is a traveling monochromatic wave with a given frequency.Solving mathematically the dynamic equations, Eq.( 36)-(37), the response of the wing to and incoming gust is quantified, Eq. ( 64)-(65).Then the pressure near the wing is coupled to the pressure in far field in order to examine the effect of the elasticity of the wing on the far field acoustic wave.
A numerical procedure is written using FORTRAN language under Linux to find the eigenvalues of the homogeneous dynamic equations describing the interaction between the wing and the surrounded fluid, the curve are plotted using "G N U P L O T" graphic procedure under Linux.In order to find the eigenvalues, a procedure is written to find the roots of equation 48.This procedure establishes functional relationship between the frequency and the temporal amplification/dampening rate s and the wave number noted k i , i = 1, ..., 6. Newton-Raphson method is used to fin the root s of the characteristic equation 55, namely, Det(M) = 0.At each Newton-Raphson iteration, equation 48 are solved to find the new values of k i (s).The procedure is repeated until convergence.In order to compute the response functions, equations 64 and 65 give h and β for a given amplitude A and given reduced frequency k of the incoming gust.Only the forced part of the solution is taken hare (The part brackets).Thus we suppose the other eigenvectors are not exited.As a result, the displacement h and the angle β are obtained as function of the of A and k.Remember that in this part of the solution s = 2ik because this is the contribution coming from the pole s = 2ik.Then the response functions are defined as a transfer function relating the input/output data i.e. (h/A, β/A) = (R(k);M(k)).The acoustic waves are computed as follow: For an incoming gust of amplitude A and reduced frequency k, we solve equations equations 61 and 62 for ( h, β), then equations 20-22 give A 0 , A 1 , A 2 , then from equation 78 we deduce the pressure of the acoustic wave for a given s, k, r and ϕ.
In all the numerical experiment undertaking here the dimensionless Young modulus and shear modulus are taking such that E = G, the dimension less wingspan l = 10, x 0 = 1/4 and x α = 0.1.Figure 4 shows the imaginary part of the response function versus the real part when the reduced frequency, k of the incoming gust varies.In this figure s = ik.The response function presented here take account to the forced part of the solution only.Thus, we consider only the term between brackets in Eq. ( 64)-(65).We disregarded the remaining terms because for large t they go to zero.However, the disregarded terms are very important in transitional state when the incoming gust lands on the wing, i.e. for short time.The evaluation of those terms is in hand.The response functions introduced here are the ratios of the force and moment applied to the blade, to the amplitude of the incoming gust.The loop appeared in the response functions, in figure (4) , is due to the fact that the frequency of the incoming gust in this point is close to the proper frequency of the wing-flow, in this case s = -0.0314+ 3.5158i.Figure (4) proves that the elasticity of the wing may not be ignored if the frequency of the incoming gust is close to the proper frequency of the wing-flow.The response function related to the moment applied to the wing is shown in fgure (5).The loop appeared in the fgure occurs too when the frequency of the incoming gust is close to the proper frequency of the wing-flow.In figures ( 4) and ( 5) , s = -0.0314+3.5158i.

CONCLUSION
In this paper, we consider the behavior of an elastic wing under the action of an incoming gust.The two dimensional non viscous flow theorem associated with strip theorem are used, to describe the wing response to an exterior solicitation.The Euler-Bernoulli approximation are used to describe the wing motion.The flow and the wing motion are coupled via the boundary condition at the surface of the wing.The resulting equation are solved analytically and the velocity of the fluid, the displacement of the wing and the pressure in the fluid are computed.The response functions defined as the ratio of the force and moment applied by the fluid to the wing to the amplitude of the incoming gust is computed for a range of the frequency of the incoming gust.The pressure in far field are computed, and the effect of the elasticity of the wing on the acoustic wave is predicted.It is found that the effect of the elasticity becomes very important when the frequency of the incoming gust is close to the proper frequency of coupled system, i.e. wing-flow.This effect appeared as a loop in the response function when the frequency of the incoming gust is close to the proper frequency of the coupled wing-flow system.Further, for a small Mach number, we find that the intensity of the acoustic wave, radiated out by the wing as a response to the incoming gust, is directly related to the unsteady force and the unsteady moment applied by the incoming gust onto the wing.The force and the moment applied by an acoustic wave due a source located at some distance from the wing could be estimated by our solution.As expected, the fatigue caused by such acoustic wave appears to be important if the acoustic wave has a frequency close to the frequency of the wing-flow system.

( 2 )
Thus, the lift applied to a strip of unity length of the wing reads (3) At the surface of the wing the impermeability condition implies that the velocity component in z direction is (4) Y = Y (x) being the equation of the surface of the wing.As usual in classical aerodynamics theorem we introduce a new variable θ.For (5) we expand γ as follows (6) Following Katz and Plotkin (2001), we substitute the last equation in the second equation of the set of Eq. ( The equation of motion of a cantilever 2 wing are (27) (28) In our equations, the term L and M contain the forcing terms represented by the Bessel functions.These terms are absent in the equation formulated by Fung (1993).Eliminating L and M leads to an explicit form of the cantilever wing equations, namely (29) Int.Jnl. of Multiphysics Volume 1 • Number 2 • 2007 (30) introduce Laplace transform of h and β (35) Int. Jnl. of Multiphysics Volume 1 • Number 2 • 20075 HOMOGENEOUS SOLUTIONDropping the forcing terms from the dynamic equations, i.e.Eq. (36)-(37), leads to a homogeneous system of differential equations.Let the solution of the homogeneous system of equations be ( hh , βh ) which satisfies the following equations Looking for a solution in the form of a normal mode, eliminate h from the two precedent equations to obtain the dispersion relation (48)5.2EIGENVALUES AND EIGENVECTORS OF THE SYSTEMIn the last section, we found that the dispersion equation has 6 roots, to say, (49) hence, we can solve for k as function of s.Therefore, the homogeneous solution is () is the determinant of the matrix M which is given by Eq. (55).The matrices M j are obtained from the matrix M by replacing the column j by the vector ( hp , 0, βp , 0, 0, 0) T .The superscript T stands for the vector transpose.Using the explicit form of the particular non homogeneous solution, the general solution becomes (60) In order to reveal explicitly the singular points of the functions h and β , we replace the particular solution by its value.Hence (61) (62) where the matrices B j is obtained from the matrix M by replacing the column j by the vector (b 1 c 2 -b 2 c 1 , 0, a 2 c 1 -a 1 c 2 , 0, 0, 0) T .It should be noted that in the last equation the functions Det(B j ), Det(M) and ∆ p depend on s.The roots of the functions Det(M) and ∆ p are the poles of the functions h and β.To express the displacement and the torsion as function of time, one has to invert Laplace transform.The inverse formula are (63) where r 2 > ζ 2 .Int. Jnl. of Multiphysics Volume 1 • Number 2 • 2007 7.1 FAR FIELD PRESSURE FLUCTUATION It is well known, see Fung (1993), that the far field pressure fluctuation in a frame attached to the wing satisfies the following wave equation, namely (69) Let replace the dimensional variables (t, x, y) by the dimensionless variables (tU/u, x/c, y/c), then the wave equation becomes (70) M = U a is the Mach number and a is the speed of sound.We introduce the new coordinate then the solution in cylindrical coordinate ( r, φ ) reads (71)

7. 2
LOW MACH NUMBER APPROXIMATION At low Mach number, φ ≈ ϕ, r ≈ r a , δ ≈ 0. The pressure on the surface of a circle of radius r c which include the wing is Let expand the pressure in Fourier series, namely Thus, a 0 = 0 and for j ≠ 0, the Fourier coefficient are such that 7.3 FIRST ORDER APPROXIMATION For large r c one can approximate the integrand (ξ 2 + r 2 c -2r c ξ sin(ϕ)) -1 using Taylor expansion, Thus, all a j = 0, but (72) and all the b j = 0, but (73) and the pressure at the circumference of the circle is (74) Matching the coefficients of the trigonometric functions in Eq.(74) and (71), and taking into account that ϕ ≈ ∼ ϕ, the far field pressure has to satisfies (75)

Figure 2
Figure2shows the real part versus the imaginary part of the eigenvalues of the system.For the configuration taking into account here, all the modes are dampened (negative real part of the eigenvalues s)

Figure ( 3 )
Figure(3) shows the variation of the displacement h of elastic axis and the angle β versus time at different station at the wing.As expected the displacement and the angle are in phase and an increase in the angle β leads to an increase in the lift, and consequently to the elevation of the wing.Note that the amplitude of the displacement h and the angle β are a monotone function of the distance y.However, h and β might exhibit an oscillation in y direction if the eigenvalues of the system of high frequency are exited by the incoming gust.

Figure 4 1 .
Figure 4 The imaginary part, R i , versus real part, R r , of the response functions at some station localized at various blade span wise distance.The response function in this figure is the ratio of the force applied to the wing to the amplitude of the incoming gust.The curve is obtained by varying the s.(a): y l = 0, (b): y l = 1 2 , (c): y l = 1.The loop appeared in the response function occurs when the frequency of the incoming gust is close to the proper frequency of the wing-ow.C = 1, EI = GJ = 10 6 , I α = 0.1, l = 10, s = 0.1, x c = 0.25, x 0 = 0.1, ρ = 100, Only the forced part of the solution is considered here, the part between brackets in equations 64-65, in which s = 2ik.

Figure 5
Figure 5 The imaginary part, M i , versus real part, M r , of the response functions at some station localized at various blade span wise distance.The response function in this figure is the ratio of the moment applied to the wing to the amplitude of the incoming gust.The curve is obtained by varying s. (a): y l = 0, (b): y l = 1 2 , (c): y l = 1.The loop appeared in the response function occurs when the frequency of the incoming gust is close to the proper frequency of the wing.C = 1, EI = GJ = 10 6 , I α = 0.1, l = 10, S = 0.1, x c = 0.25, x 0 = 0.1, ρ = 100, Only the forced part of the solution is considered here, the part between brackets in equations 64-65, in which s = 2ik..

Figure ( 6 )
Figure(6) shows the acoustic wave in far field versus the radial distance and the frequency of the incoming gust.It appears that, the acoustic wave is significantly affected by the elasticity of the wing if the frequency of the incoming gust is close to the proper frequency of the wing-flow.

Figure 6
Figure 6 The acoustic wave in far field versus the radial distance and the reduced frequency f = k of the incoming gust.One can see the effect of the elasticity of the wing when the incoming frequency f = k is near the proper frequency of the wing which is in this case s = -0.0314+ 3.5158i.C = 1, EI = GJ = 10 6 , I α = 0.1, l = 10, s = 0.1, x c = 0.25, x 0 = 0.1, ρ = 100.