cylindrical shells conveying subsonic compressible fluid flows with general boundary conditions

Author’s Accepted Manuscript Dynamics and stability of conical/cylindrical shells conveying subsonic compressible fluid flows with general boundary conditions M. Rahmanian, R.D. Firouz-Abadi, E. Cigeroglu www.elsevier.com/locate/ijmecsci

PII: DOI: Reference:

S0020-7403(16)30589-6 http://dx.doi.org/10.1016/j.ijmecsci.2016.10.037 MS3484

To appear in: International Journal of Mechanical Sciences Received date: 19 June 2016 Revised date: 20 September 2016 Accepted date: 26 October 2016 Cite this article as: M. Rahmanian, R.D. Firouz-Abadi and E. Cigeroglu, Dynamics and stability of conical/cylindrical shells conveying subsonic compressible fluid flows with general boundary conditions, International Journal of Mechanical Sciences, http://dx.doi.org/10.1016/j.ijmecsci.2016.10.037 This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting galley proof before it is published in its final citable form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

Dynamics and Stability of Conical/Cylindrical Shells Conveying Subsonic Compressible Fluid Flows with General Boundary Conditions M. Rahmaniana,b , R.D. Firouz-Abadia,∗, E. Cigeroglub a Department

of Aerospace Engineering, Sharif university of Technology, Tehran, P.O.Box 11155-8639, Iran of Mechanical Engineering, Middle East Technical University, 06800 Ankara, Turkey

b Department

Abstract A fast and efficient reduced order formulation is presented for the first time to study dynamics and stability of conical/cylindrical shells with internal fluid flows. The structural and fluid formulations are developed based on general assumptions to avoid any deficiency due to modeling. Their respective solutions and the final solution to the coupled field problem are also developed in a way to be capable of capturing any desirable set of boundary conditions. In addition to the flexibility provided by the solution methodology and generalization provided by the formulation, current solution proposes an additional advantage over others which is the minimal computational cost due to the special reduced order model proposed. Therefore, stability margins of the problem at hand can be obtained both efficiently and accurately. Proposed formulation is verified by comparing the results of the present study with the results available in literature for cylindrical/conical shells at different boundary conditions. Comprehensive parameter studies are performed in order to draw general insights over the effects of boundary conditions, semi-vertex angle and compressibility on the dynamics and stability margins of conical shells with internal fluid flows. Keywords: Conical shells, Fluid-Structure Interaction, Reduced Order Model, Compressibility Effects, Static Condensation.

1. Introduction There are a variety of engineering problems where a structural component is in contact with some sort of stationary or flowing fluid. Such problems fall into the category of fluid-structure interactions where the coupled system is inherently different from any of the segregated fluid or structural components. Many practical applications of such problems can occur in shell-type structures. Among all the available shell geometries, cylindrical shells have attracted the most attentions and a very detailed insight is obtained during the past decades over the problems of cylindrical shells interacting with fluids. Conical shells interacting with fluid flows, on the other hand, are not fully addressed in literature. The few available studies on conical ∗ Email:

fi[email protected] ; Telephone: +982166164606

Preprint submitted to International Journal of Mechanical Sciences

November 10, 2016

shells loaded with fluid have also arisen fundamental negotiations. Moreover, cylindrical shells and annular plates are special cases of conical structures; hence, a detailed understanding of such problems can provide further insights over other practical geometries as well. Some practical and sensitive applications of conical shells conveying fluid flows include wind and water tunnels which are mostly used in the aerospace and marine engineering applications to simulate a special condition, civil and military aircraft engines and space shuttle nozzles. Pipelines conveying pressurized and high-velocity refinery products such as gas or fuels are another major examples of shells with internal compressible fluid flows. In the mentioned applications, the internal medium is usually in the high-subsonic regime where compressibility effects are significant; Hence, compressible fluid models are advised for subsonic regimes in such applications to obtain more accurate results, specially, if light fluids are used as the flowing medium. We will first mention the most reputable studies on cylindrical shells and then confine ourselves to the case of conical shells. The first serious studies on vibrations and stability of cylindrical shells dates back to 1950s. First studies were mainly devoted to the development of governing equations and then providing some solution methodologies. Flutter of thin cylindrical shells conveying fluid by Paidoussis et. al. [1] was among the first reputable studies of that time. They used Flugge assumptions for the structural description and the potential flow to model the fluid inside the shell. Their results showed that increasing the fluid velocity leads to a flutter instability for a clamped-free cylindrical shell. If the boundary condition is altered to clampedclamped, the instability will also change qualitatively to a divergent type followed by a coupled-mode flutter. Their study was followed by a similar investigation by Weaver et. al. [2]. They used the Flugge-Kempner equations and potential flow assumptions to study dynamic stability of fluid conveying pipes. Paidoussis then tried to resolve a paradox arose in his previous studies [3]. He discussed that both shell and beam models predict the same behaviors qualitatively and as the length of the shell increases, both beam and shell models predict identical results. Moreover, existence of coupled-mode flutter in beam modes of a simply supported pipe conveying fluid was shown to be possible for the first time. Actually, it was believed for several years that pipes conveying fluid flows can not have coupled mode flutter instabilities. Another interesting study was conducted by Everstine [4] proposing a new symmetric formulation in the solution of the potential flow inside the shell. This can significantly reduce computational costs by using the advantages of symmetric formulation and also remove the fictitious modes generated in non-symmetric formulations. Selmane et. al. [5] have also considered the problem of fluid loaded cylindrical shells submerged and conveying fluid at the same time. They have also proposed a novel method of solution, known as the hybrid finite element method, to solve the coupled fluid-structure system of equations. Detailed stability analysis of both isotropic and composite cylindrical shells conveying internal fluid flows are also provided by Kochupillai et. al. [6, 7]. Amabili has deeply studied the problem of cylindrical shells in contact with fluid flows for linear [8, 9, 10] and nonlinear [11, 12, 13, 14] cases. Finally, a great and comprehensive review on the dynamics of pipes and cylindrical shells with fluid flows is provided by Paidoussis in a two-volume book series [15, 16]. The 2

two books have summarized approximately all the literature and findings until their publication date. Dynamic behavior of conical shells in vacuo has been the subject of many studies [17, 18, 19, 20, 21, 22, 23, 24, 25] while the first serious study on vibration and stability of conical shells with air loadings was conducted by NASA in 1970 [26]. Donnell equations of motion along with inviscid two dimensional quasisteady aerodynamic model were used to describe the problem and an approximate solution was obtained. They have verified their results with similar results in cylindrical shells. There is a large gap in literature on stability analysis of fluid (air/water) loaded conical shells between the 1970s and 1990s and only a few publications can be found. During these two decades, Ueda [27] investigated theoretically and experimentally the supersonic flutter of a conical shell. Bismarck-Nasr et. al. [28] presented a finite element solution procedure for supersonic flutter predictions of conical shells with external flows. Sunder et. al. [29] have also performed a detailed analysis to determine the optimum cone angles in an aeroelastic sense. They employed the finite element solution method along with the static condensation procedure and specifically studied the clamped-clamped conical shells under supersonic external flows. They utilized the aerodynamic piston model to capture the fluid behavior. Ultimately, they proposed an empirical relation determining the flutter boundary pressure and consequently obtaining the optimum vertex angle with the most resistance to flutter. In a latter study, during the 1980s Mason et. al. [30] performed an interesting study regarding the aeroelastic instabilities in rocket nozzles using the finite element approach. They modeled the nozzle by means of thin, truncated conical shells subjected to internal supersonic flow. The supersonic flow on the other hand is modeled by isentropic gas relations. Several test cases were provided to validate different aspects of the study. Dynamic behavior of anisotropic conical shells with internal fluid was studied by Lakis et. al. [31]. A novel finite element approach was introduced where conical finite elements were used to discretize the domain. The velocity potential formulation and Bernoulli equation were used to describe fluid inside the shell. Conical shells with internal supersonic fluid flows were stuided by Aleksandov et. al. [32]. The conical shell in their study was clamped-free with the larger end being free. A supersonic gas flow entering the shell from the smaller and clamped edge is considered and its dynamic loading is obtained by the aerodynamic piston theory. In another research by Jeong et.al. [33], a theoretical formulation was suggested on free vibrations of thin conical shells filled with an ideal fluid with both simply supported and clamped boundary conditions via a finite element approach. Jhung et. al. [34] also studied the fluid-structure interaction of bridged conical shells filled with ideal fluid. A velocity potential formulation was used to model the internal fluid and a finite element solution was provided for the coupled system of equations. Caresta et. al. [35] have developed an analytical model to study the modal vibrations of fluid loaded conical shells in the low frequency range. Solution to the shell equation was presented in a power series format. The fluid loading was also obtained by dividing the shell into a finite number of narrow strips which were considered to be locally cylindrical with a specific mean radius. However, it is identified that such a model is feasible only in the low frequency range. Sabri et. al. [36] have studied the supersonic aeroelastic analysis of conical shells 3

by means of a novel hybrid finite element solution procedure. The shell is governed by the Sanders equations of motion and the fluid flow is described by the linear first order piston theory with curvature corrections. They have investigated empty and partially filled conical shells. Mahmoudkhani et. al. [37] also studied supersonic flutter of functionally graded concial shells. They included the aerodynamic heating effects and obtained the aero-thermoelastic stability margins. Love’s assumptions were utilized to derive the structural equations and the fluid flow was modeled by the first order aerodynamic piston theory. Recently, dynamics of conical shells with subsonic internal flowing fluid has attracted attentions. In a first study, S. Kumar et. al. [38] proposed a semi-analytical finite element formulation to study vibrations and stability of conical shells conveying fluids. They used Bernoulli equation to determine pressure acting on the walls. In addition to free vibration characteristics of the coupled-field problem, they obtained the stability margins and critical fluid velocities. In a later study by Kerboua et. al. [39], a hybrid finite element method was proposed to study vibrations of conical shells with internal flowing fluid. In the proposed methodology, the displacement functions were derived from exact solutions of Sanders equilibrium equations. The fluid flow was modeled by the velocity potential equation and pressure acting on the walls was estimated by Bernoulli equation. Hybrid finite element solution was eventually shown to be consistent with the semi-analytical solution provided by [38] except for the critical velocities. The third major contribution in vibrations of conical shells with flowing subsonic internal fluid is the recent study by Bochkarev et. al.[40]. Authors, used finite element method to solve both structural and fluid equations. Stability analysis of both cylindrical and conical shells are as well provided and effects of different boundary conditions of the velocity potential equation was investigated. Both studies by Kerboua et. al. [39] and Bochkarev et. al.[40] are compared to the first study on this subject by S. Kumar et. al. [38]. Comparisons show that major differences in prediction of the critical fluid velocity exists, especially at lower and higher circumferential wave numbers. It’s also noteworthy that all the studies have neglected the density variations due to compressibility effects and have obtained the mean flow velocity by a simple approximation of the conservation of mass. Moreover, in both studies by S. Kumar et. al. [38] and Kerboua et. al. [39], boundary conditions on the fluid flow is not discussed which is shown to be significantly effective on the final outcomes of the analysis [40]. A final remark is that even for the smallest fluid velocity, the system will be non-conservative due to the gyroscopic effects appearing in the coupling conditions. In previous studies [38, 39], the final non-conservative system of equation is treated as a conservative system which means that the frequency variations are the only criteria to determine the stability boundaries. This is not true for such non-conservative systems; consequently, the instabilities reported are both quantitatively and qualitatively subject to change when taking variations of damping into account. In the present study, a conical shell with internal subsonic flowing fluid (i.e, water or air) is considered. The fluid is assumed to be compressible, irrotational and inviscid and the governing equations of motion are obtained in their most general linear representation. The structural equations of motion are also obtained 4

using the variational principle and taking both rotary inertia and shear deformation effects into account. An extra modification is also introduced in the structural equations of motion to remove the inconsistency of rigid body rotations in the solution. This modification is usually known as the Sander’s modification. Actually, the structural and fluid governing equations are obtained in their most general forms to avoid any accuracy loss during the analysis and justify the discrepancies in literature more accurately. The main novelties of the current study can be divided into two categories. The first one being the solution methodology. The structural equations of motion are solved by means of a very flexible and rapid-convergent semi-analytical Galerkin type method and the fluid equation of motion is solved by the finite element method. Although, the required degrees of freedom in the structural and fluid solutions are minimized by selecting the most proper trial/shape functions, a reduced order model is also developed to further minimize the degrees of freedom in the fluid solution. The final coupled field solution is significantly faster than other conventional solution algorithms. The second novelty can be attributed to the investigation of different boundary conditions of the structural and fluid domains. Verification studies are then conducted for both fluid loaded conical and cylindrical shells. Furthermore, stability analysis of conical shells at different geometries and boundary conditions are provided to illustrate general insights over the stability analysis of conical shell structures.

2. Governing Equations 2.1. Structural Equations of Motion In general, the three dimensional (3D) elasticity relations should be used to obtain the equations of motion for a general structure. Since such equations are very complex to solve, some simplifying assumptions have to be introduced to reduce the complexities of the problem at hand. Shells are among the most practical structures where several simplifying theories have been introduced to provide such reductions. All shell theories have been proposed to simplify a specific category of problems with the minimum loss of accuracy. For example, in the case of very thin shells, including effects of rotary inertia and shear deformation can at most alter the final results up to one percent; hence, these two effects are negligible in nearly all thin shell theories. Actually, the cost of getting more accurate results should be justified. The degree of accuracy must be tailored by selecting a proper shell theory. The basic idea behind all shell theories is to reduce the 3D problem to an equivalent 2D representation by approximating the displacements across the thickness and then solving for the displacements on the mid-surface. Discussions on the specific assumptions made in each theory and their validity is beyond the scope of this study and interested readers are referred to the monograph by Leissa [41]. Different classifications are proposed for shell theories. One basic classification is based on the degrees of freedom included and order of approximation which are namely as: Classical Shell Theories (CST), First order Shear Deformable Theories (FSDT) and Higher order Shear Deformable Theories (HSDT). 5

Classical shell theories are fundamentally developed for thin shells according to the Kirchhoff-Love kinematic assumption. Therefore, the effects of shear deformation and rotary inertia are neglected in this category of theories. FSDTs are one step ahead of CSTs. In the FSDT category of theories, the effects of shear deformation and rotary inertia are included. Although FSDTs are much more accurate compared to CSTs, they still deliver inaccurate results if the shear correction factor is not justified accurately for the problem at hand, or if the thickness of the shell is significantly large. To remove these inconsistencies, the most accurate and computationally expensive category, i.e. HSDT, are introduced. Even though, HSDTs are the most accurate, one should trade off between the degree of accuracy needed and the computational cost of solution. Another important point is that for metallic structures or thin/moderately thick shells, FSDTs are the most reasonable choice, since predictions of HSDTs and FSDTs are closely identical for such problems [42] where no stretching/contraction occurs in the thickness direction. In the current study, it is as well intended to have an algorithm that is capable of modeling small or thick conical shells where shear deformation effect plays an important role and hence, accurate results can be obtained. Moreover, in order to make comparisons with the available results in literature for short conical shells, a shear deformable theory is more reasonable and reliable. In conclusion, FSDTs have proved themselves to be sufficiently genuine in practical applications. They provide reasonable compromise between the solution accuracy and the computational effort required. In the current study, the Flugge FSDT shell theory is chosen to derive the governing equations. One further refinement, which is the Sander’s modification, is also appended to the final formulations to remove the rigid body rotation inconsistency [43]. This inconsistency is present in the original form of most shell theories and hence, leads to rigid body rotations in cases where a torsional loading is applied to the structure. Figure 1 demonstrates a schematic representation of a conical shell of length l, semi-vertex angle α and thickness h. Variables u, v and w represent the longitudinal, circumferential and radial displacements, respectively. The origin of the cylindrical coordinate (s, θ, R) is located on the mid-surface and at the input edge of the cone. The coordinate variables s, θ and R represent the longitudinal, circumferential and radial directions, respectively. A set of three translational and two rotational springs (corresponding to each degree of freedom) with the respective spring stiffnesses of ku , kv , kw , kφ1 , kφ2 are attached to the mid-surface of the cone at both ends. Proper justification of the spring stiffnesses can give highly accurate results for the corresponding boundary condition. For the sake of simplicity, only one set of spring attachments are shown in Fig. 1. The input fluid velocity is denoted by U and the radius of curvature of the cone is r, which is perpendicular to the middle surface of the cone. For conical shells the following definitions are used to obtain r and R values. R

=

r

=

R1 + s sin α R1 + s tan α cos α 6

(1) (2)

where R1 is the input radius of the cone.

⎧ ⎪ ⎪ U (s, θ, z, t) ⎪ ⎨ 1





Based on assumptions in FSDT, the displacement field U = ⎫ ⎪ ⎪ ⎪ ⎬

⎧ ⎪ ⎪ u(s, θ, t) ⎪ ⎨

⎫ ⎪ ⎪ ⎪ ⎬

U1

U2

U3

⎧ ⎪ ⎪ φ (s, θ, t) ⎪ ⎨ 1

takes the following form.

⎫ ⎪ ⎪ ⎪ ⎬

= +η v(s, θ, t) U2 (s, θ, z, t) φ2 (s, θ, t) ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ U (s, θ, z, t) ⎭ ⎩ w(s, θ, t) ⎭ ⎭ ⎩ 0 3

(3)

where φ1 and φ2 are the rotations of the normals to the mid-surface around θ and s coordinates, respectively; and η is the thickness variable being defined over the interval [− h2 , h2 ]. To obtain the governing equations of motion for the conical shell, Hamiltonian of the structure is constructed as follows  δH =

t2 t1

δ(T − V + Wnc )dt = 0,

(4)

where T and V are the kinetic and potential energies of the shell and Wnc is the work done by nonconservative external forces. Kinetic energy of the shell is defined as 1 T = 2

θ s  0

0

h 2

−h 2

η ρs (U˙ 12 + U˙ 22 + U˙ 32 )R(1 + ) dη ds dθ, r

(5)

where over dot symbol represents differentiation with respect to the temporal coordinate and ρs is the density of the structural material of the shell. Substituting the displacement field of Eq.(3) into the kinetic energy definition given by Eq.(5) and performing some manipulations gives the following variational representation of the kinetic energy θ s δT = − 0

0

⎡ ⎣

ρ(¨ uδu + v¨δv + wδw) ¨ + I2 (φ¨1 δφ1 + φ¨2 δφ2 ) +I1 (φ¨1 δu + φ¨2 δv + u ¨δφ1 + v¨δφ2 )

⎤ ⎦ Rds dθ,

(6)

where the following definitions are employed T

[ρ, I1 , I2 ] =



h/2

 T η ρs 1, η, η 2 (1 + )dη. r −h/2

(7)

In the next step, a general definition is obtained for the work done by external forces on the structure. The primary external loading on the shell originates from the internal fluid. The amount of pressure exerted on

7

the surface of shell is denoted by P . Therefore, the following variational formulation is determined 

θ



s

δWnc =

P δwR dsdθ. 0

(8)

0

The last expression needed to construct Hamiltonian of the system is the potential energy expression which is the summation of both of the elastic strain energy and the potential energy stored in the boundary springs. The final potential expression in variational representation is defined as δVtotal

=

δVelastic

=

δVsprings

=

δVelastic + δVsprings ⎡ ⎤ σ (δε + ηδκ ) + σ (δε + ηδκ ) 1 22 2,0 2  θ  s  h/2 ⎢ 11 1,0 ⎥ η ⎢ ⎥ 0 1 0 1 ⎢ +τ12 (δω1 + ηδw1 ) + τ21 (δω2 + ηδw2 ) ⎥ R(1 + )dηdsdθ r ⎣ ⎦ 0 0 −h/2 +τ23 δγ23,0 + τ13 δγ13,0 ⎡ ⎤  θ s s1 (kus1 uδu + kvs1 vδv + kw wδw + kφs11 φ1 δφ1 + kφs12 φ2 δφ2 )s=s 1 ⎣ ⎦ Rdsdθ, s2 0 0 +(kus2 uδu + kvs2 vδv + kw wδw + kφs21 φ1 δφ1 + kφs22 φ2 δφ2 )s=s

(9)

2

where kij denotes the spring stiffness corresponding to the ith degree of freedom at the jth position. The strain components are also defined as follows ⎧ ⎪ ⎪ ε1,0 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ κ1 ⎪ ⎨ ω10 ⎪ ⎪ ⎪ ⎪ ⎪ ω11 ⎪ ⎪ ⎪ ⎪ ⎩ γ 13,0 ⎧ ⎪ ⎪ ε2,0 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ κ2 ⎪ ⎨ ω20 ⎪ ⎪ ⎪ ⎪ ⎪ ω21 ⎪ ⎪ ⎪ ⎪ ⎩ γ 23,0

⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎬ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎭ ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎬

=

⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎭

=

⎧ ⎫ ⎧ ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ u 0 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ φ1 ⎪ 0 ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎬ ⎨ ⎬ ∂ + v 0 ⎪ ⎪ ⎪ ∂s ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ φ2 ⎪ ⎪ 0 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ w ⎭ ⎩ φ ⎪ ⎭ 1 ⎧ ⎧ ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ v ⎪ u ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ φ1 ⎪ φ2 ⎪ ⎪ ⎨ ⎬ sin α ⎪ η −1 ∂ ⎨ + (1 + ) ( u −v ⎪ r R∂θ ⎪ R ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ 1 ⎪ ⎪ ⎪ ⎪ ⎪ −φ2 [1 + ( r − ⎪ φ1 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ w ⎭ ⎩ 0

(10)

⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎬ ⎪ ⎪ ⎪ R ⎪ ⎪ )η] 2 ⎪ r cos α ⎪ ⎪ ⎪ ⎭

+

⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ 1⎨

w 0

0 r⎪ ⎪ ⎪ ⎪ ⎪ 0 ⎪ ⎪ ⎪ ⎪ ⎩ rφ − v 2

⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎬ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎭

)(11)

The following stress-strain relation holds for a homogeneous orthotropic material ⎧ ⎪ ⎪ σ11 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ σ ⎪ ⎨ 22 τ12 ⎪ ⎪ ⎪ ⎪ ⎪ τ13 ⎪ ⎪ ⎪ ⎪ ⎩ τ 23

⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎬

⎡ c11

⎢ ⎢ ⎢ c12 ⎢ ⎢ =⎢ 0 ⎢ ⎪ ⎪ ⎪ ⎢ ⎪ ⎪ ⎢ 0 ⎪ ⎪ ⎣ ⎪ ⎪ ⎭ 0

c12

0

0

c22

0

0

0

c66

0

0

0

c66

0

0

0 8

⎤⎧ ⎪ ⎪ ⎥⎪ ⎪ ⎥⎪ ⎪ ⎪ 0 ⎥⎪ ⎥⎪ ⎨ ⎥ 0 ⎥ ⎥⎪ ⎥⎪ ⎪ ⎪ 0 ⎥⎪ ⎪ ⎦⎪ ⎪ ⎪ c66 ⎩ 0

⎫ ⎪ ε11 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ε22 ⎪ ⎪ ⎬ , γ12 ⎪ ⎪ ⎪ ⎪ γ13 ⎪ ⎪ ⎪ ⎪ ⎪ γ23 ⎭

(12)

where the constant coefficients cij are defined in terms of Young’s modulus (E) and Poisson’s ratio (ν) as c11

=

c12

=

c66

=

E 1 − ν2 Eν c21 = 1 − ν2 E . 2(1 + ν) c22 =

(13)

Implementation of all the aforementioned expressions result in the governing equations of motion with the rigid body rotation inconsistency. Based on this inconsistency, a rigid body rotation applied to the structure gives a nonvanishing torsion except for flat plates and spherical shells [43]. In order to remove this inconsistency, the equilibrium equation about the normal to the differential element of the shell must be satisfied, which can be obtained by the modification of following strain components in Eqs.(10-11) as [43] ω ˜ 10

=

ω10 − φ3

(14)

ω ˜ 20

=

(15)

ω ˜ 21

=

ω20 + φ3 φ3 , ω21 + r

(16)

where φ3 =

η ∂u 1 ∂v (1 + )−1 (R + v sin α − ). 2R r ∂s ∂θ

(17)

Employing force and moment resultant definitions for the shear correction factor Ks = 5/6 as ⎧ ⎪ ⎪ N11 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ N ⎪ ⎨ 12 Q13 ⎪ ⎪ ⎪ ⎪ ⎪ M11 ⎪ ⎪ ⎪ ⎪ ⎩ M 12 ⎧ ⎪ ⎪ N22 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ N ⎪ ⎨ 21 Q23 ⎪ ⎪ ⎪ ⎪ ⎪ M22 ⎪ ⎪ ⎪ ⎪ ⎩ M21

⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎬



⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎭ ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎬

=

⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎭

=

h/2

⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨

σ11 τ12

⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎬

η (1 + ) dη Ks τ13 ⎪ ⎪ r −h/2 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ σ η 11 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ τ η ⎪ ⎭ 12 ⎧ ⎫ ⎪ ⎪ ⎪ ⎪ σ 22 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ τ 21 ⎪ ⎪  h/2 ⎨ ⎬ dη, Ks τ23 ⎪ −h/2 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ σ22 η ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ τ21 η ⎭

(18)

(19)

and using Eqs.(6,8,9) in the Hamiltonian relation of Eq.(4) and utilizing the Flugge theory assumption, i.e. keeping the thickness variables up to the 3rd order after applying the Taylor expansion for (1 + ηr )−1 , leads 9

to the following integral equation of motion. ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ 2π l ⎢ ⎢ ⎢ ⎢ ⎢ 0 0 ⎢ ⎢ ⎢ ⎢ ⎢ ⎣

ˆ 12 ∂δu N22 M ∂δu R sin α)Rδu + RN11 ∂s + (N21 − 2r ) ∂θ ˆ 21 ∂δv M +(ρ¨ v + I1 φ¨2 − NR21 sin α − Qr23 )Rδv + N22 ∂δv ∂θ + R(N12 + 2r ) ∂s ∂δw +(ρw ¨ − P + Nr22 )Rδw + RQ13 ∂δw ∂s + Q23 ∂θ ∂δφ1 1 +(I2 φ¨1 + I1 u ¨ + Q13 + MR22 sin α)Rδφ1 + RM11 ∂δφ ∂s + M21 ∂θ ∂δφ2 2 +(I2 φ¨2 + I1 v¨ + Q23 − MR21 sin α)Rδφ2 + M22 ∂δφ 12 ∂s ∂θ + RM   s1 + (kus1 uδu + kvs1 vδv + kw wδw + kφs11 φ1 δφ1 + kφs12 φ2 δφ2 ) s=s1  s2 + (kus2 uδu + kvs2 vδv + kw wδw + kφs21 φ1 δφ1 + kφs22 φ2 δφ2 ) s=s2

(ρ¨ u + I1 φ¨1 +

ˆ 12 = M ˆ 21 = Here, M

 h/2 −h/2

⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ds dθ = 0. ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦

(20)

ητ12 dη and δX , X ∈ {u, v, w, φ1 , φ2 } is an arbitrary function which satisfies only

the geometric boundary conditions of the problem. As previously mentioned, the current solution procedure is capable of modeling all combinations of boundary conditions including elastic and classical ones. In this study, we confine ourselves to the following classical definitions of boundary conditions Free condition (F): u = v = w = φ1 = φ2 = 0. Clamped condition (C): u = v = w = φ1 = φ2 = 0. Simply-Supported condition (SS): v = w = 0. All combinations of the above mentioned definitions, including symmetric and asymmetric ones are carefully studied. The symmetric boundary conditions involve C − C and S − S while asymmetric conditions are C − F, C − S, S − F, S − C, F − C andF − S. 2.2. Fluid Equations of Motion In order to derive a mathematical representation of the fluid inside the shell, a different treatment (compared to the previous section) is followed. Different mathematical models are proposed to simulate stationary/flowing fluid dynamics. Although, the most accurate model is the nonlinear Navier-Stokes equation which includes all the effects of viscosity, thermal dissipations, and body forces, dealing with such a complex model in the problems of fluid-structure interaction (FSI) is not a feasible choice. Therefore, it’s customary in FSI problems to neglect the viscosity effects and body forces. Also, the fluid is usually considered to be irrotational and isentropic. In the current study, the fluid perturbation amplitude is considered to be in the order of the thickness of the shell and a linear fluid vibration assumption is totally feasible. In addition, the pressure is perpendicular to the surface of shell and no cavitation occurs in the flow [38]. Starting from the Navier-Stokes equation [44], neglecting body forces and viscosity terms, the Euler ˆ equation is obtained. Since the flow is irrotational, based on Kelvin’s theorem a velocity potential Φ can be obtained so that a mean flow velocity U is be defined in terms of the velocity potential as U = ˆ Substituting the velocity potential definition into the Euler equation and performing some further ∇Φ. 10

manipulations by considering a small perturbation assumption, the following wave/acoustic equation is obtained, ˆ = ∇2 Φ

ˆ ∂ 2 ˆ ∂M 2 ∂ Φ 1 ∂ + U ) . ( Φ + c2∞ ∂t ∂s ∂s ∂s

(21)

where c∞ is the speed of sound in the fluid domain and the Mach number is defined as M = U /c∞ . Variations of the Mach number appears in Eq.(21) due to the fact that the fluid velocity, and hence Mach number, varies in the longitudinal direction in a conical conduit. Due to the compressibility assumptions made in the current study, the mean flow velocity (U ) can not be obtained by a simple relation of conservation of mass (A ×U = cte., where A is the area of the corresponding section). In all previous studies of conical shells conveying fluid flows [38, 39, 40], such simple relations are used to determine the mean flow velocity. In the current study, the isentropic flow relations are utilized to determine the mean flow velocity. The main advantage of the Mach-area relations in isentropic flows is that variations of density due to the compressibility effects are also included and the relation ρf × A × U = cte. is satisfied at each cross section of the cone. Isentropic Mach-area relation is defined as [45]

γ+1   2(γ−1) γ−1 A2 M1 1 + ( 2 )M22 = 2 A1 M2 1 + ( γ−1 2 )M1

(22)

where the specific heat ratio γ is equal to 1.4 for air while for water it has to be equal to unity in an ideal case, but in order to avoid singularity issues, a value of γ = 1.0001 is taken for water containments. Subscripts 1 and 2 in Eq.(22) correspond to parameters of the first and the second cross-sections, respectively. Several physically different boundary conditions can be attributed to the fluid equation at either of input or output boundaries, which are namely zero flux ˆ n=0 ∇Φ.¯

(23)

ˆ = 0), where n and zero perturbed pressure (Φ ¯ denotes the normal vector to the undeformed configuration. The zero perturbed pressure boundary condition corresponds to open ends while the zero flux condition represents a closed end condition. As explicitly expressed by Shayo et. al. [46] and Paidoussis et. al. [47], a proper justification of the fluid end conditions can significantly alter the dynamics of the whole system. Shayo et. al. [46] and Paidoussis et. al. [47] suggested different models to damp out the perturbed velocity ˆ in a finite length after the output of the shell. They also found that the upstream fluid potential (Φ) condition has minimal effects on the dynamics of the system while the downstream conditions are vital in determining the stability margins both qualitatively and quantitatively. In the current study, the output 11

models proposed by Shayo and Paidoussis [46, 47] are not utilized. Instead, the effect of different fluid boundary conditions on the stability and dynamic behavior is studied by imposing different combinations of open and closed conditions. To address each set of boundary condition, open ends (zero pressure) are denoted by number one and closed ends (zero flux) by number zero; hence a fluid boundary condition of 10 (”fbc10”) denotes a shell with open input and closed output conditions. Similarly, ”fbc11” denotes a shell with both ends open. As stated by Shayo et. al. [46] and Paidoussis et. al. [47], and verified in the this study, the fluid boundary conditions at the input (upstream) has nearly no effect on the final results for different structural boundary conditions. This means that each of the following pairs (fbc10,fbc00) and (fbc11,fbc01) are nearly identical. Consequently, only ”fbc10” and ”fbc11” are considered as representatives of the two sets of fluid boundary conditions . In other studies by Senthil Kumar et. al. [38] and Kerboua et. al.[39], a zero flux condition is taken into account at both input and output edges (”fbc00”). 2.3. Coupling Condition Internal fluid dynamics is directly dependent upon the structural motions and boundary conditions. Conical shells - without asymmetric/localized imperfections - are axisymmetric structures therefore only two distinct boundary conditions have to be specified on the fluid-structure interface and axis of symmetry other than the input and output boundaries. As mentioned previously, at the input and output edges, a zero perturbed pressure or zero flux condition can be applied. For the axis of symmetry of the shell, a zero flux condition has to be taken into account. Such condition on the axis of symmetry implies that no fluid flows are permitted normal/oblique to the axis of symmetry and only tangential fluid flows with respect to the axis of symmetry are allowed. The only remaining condition is the one on the fluid-structure interface where the fluid and structural motions are coupled. A reasonable assumption would be the impermeability condition where the fluid and structural particles are assumed to be perfectly attached (with no cavitation or diffusion) and move with an identical velocity in the same direction. One can obtain the impermeability condition as [44]

 ˆ n ∇Φ.¯

wall

=q=

∂w ∂w +U . ∂t ∂s

(24)

Substituting Eq.(24) in the weak form representation of the fluid equation, Eq.(21), couples the fluid velocity potential to the structural states. As shown in the structural equations of motion, the pressure loading (P ) due to the internal fluid dynamics is included in the formulation in Eq. (20). This pressure can be obtained by the Bernoulli relation

12

[44] as P = −ρf (

ˆ ˆ ∂Φ ∂Φ +U ), ∂t ∂s

(25)

where ρf is the fluid density. Substitution of Eq.(25) into Eq.(20) gives the structural equations in terms of fluid states. Explicit expressions for the coupling conditions are given in the following section and a detailed discussion is also provided for solving the coupled system of equations.

3. Solution Procedure One of the main objectives of the current study is to obtain a flexible and fast convergent solution for the coupled system of equations with an acceptable accuracy. The solution procedure would be even much more fascinating if the computational efforts are minimized. To this aim, extended Galerkin approach (energy based procedure) is selected as the solution method of the structural equations and the finite element method is chosen for the fluid domain. For the structural solution, specific trial functions are proposed which are capable of capturing any set of boundary conditions including symmetric or mixed boundary conditions. The trial functions have an outstanding inherent characteristic leading to incredibly fast convergence of the results even with very small number of expansion terms. For the fluid part, quadrilateral 9-noded element is used to discretize the domain and accelerate convergence of the results (minimize the number of elements for a reasonable convergence). Finally a reduced order model is developed along with the static condensation technique to minimize the number of required degrees of freedom in the fluid domain. The final reduced order system of equations with the static condensation has shown to be much faster than the original system of equation. In the following subsections, details of the structural, fluid and coupled field analysis are presented. 3.1. Solution of the Structural Equations In the original Galerkin approach, a differential form of equations has to be used as well as the corresponding boundary conditions. The solution must be a summation of a set of independent comparison functions (satisfying all the available boundary conditions). If the functions included in the solution are not independent, the final system of equation is indeterminate (rank-deficient) with the degree of indeterminacy equal to the number of dependent functions. Suggested comparison functions must only satisfy all the available boundary conditions thus the only remaining task would be justification of the expansion coefficients so that the partial differential equations are as well satisfied. Therefore, substituting the suggested solution into the governing equation results in a finite residual term, i.e. R. A specific weighting function, say ψi (x), is introduced and the following weighted integral is satisfied. Note that this weighting function is identical 13

to the comparison functions used in the solution expansion. ∫ ψi (x)R(x, cj )dΩ = 0 ; i, j = 1, 2, ..., N.

(26)

Ω

If the temporal coordinate is identically satisfied via a harmonic motion assumption, the above relation actually introduces N algebraic equations. Satisfaction of the obtained system of equations is achieved by a proper justification of the unknown coefficients cj which are then referred as eigenvectors. The problem arises when a proper set of comparison functions are looked for, even in the simplest shell problems. This is mainly due to the complexities present in satisfaction of natural boundary conditions. One proven trick to bypass such complexities is to transform the differential equations into their integral weak form representation by utilizing integration by parts or utilizing the energy form of equations as shown in Eq.(20). In the weak form representation, all the natural boundary conditions are automatically satisfied in the integral equations; therefore, the approximate functions are only required to satisfy the geometric boundary conditions. Such functions are called admissible functions. It should be noted that in the weak form representation, the maximum order in the deliberative terms is halved; hence, a smaller order of continuity needs to be satisfied by the trial functions used. Moreover, in the integral representation, both Ritz and Galerkin methods are the same if the system under consideration is self-adjoint. The comparison/admissible functions are usually expressed in terms of trigonometric functions. Other expressions in terms of power polynomials and Chebyshev orthogonal polynomials are also reported in literature [48]. The power polynomial expansion is very flexible in satisfying the boundary conditions. The basic problem with the power expansion method is that the lower order expansions do not form a complete set and hence cannot provide an accurate and converged result. On the other hand, higher order power expansions lead to a numerically unstable solution due to the ill-conditioning of the final matrices. Such ill-conditionings are due to the round-off errors present in the integrations. Chebyshev polynomials are better than simple power polynomial functions. Since they constitute a set of orthogonal functions and as a result a well-posed system of equations will be obtained. It’s shown by Kurylov et. al. [49] that the convergence of Chebyshev polynomials for some boundary conditions are even better than trigonometric expansions. Effects of a variety of expansions (i.e. Chebyshev polynomials of the 1st and 2nd kind, Hermite polynomials and Legendre polynomials) on the vibration behavior of conical/cylindrical shells have been also investigated and similar conclusions were reached by Qu. et. al.[48]. In the current study a general approximate trigonometric solution is utilized for all sets of boundary conditions which were previously introduced by Su et. al. [50] for free vibration analysis of shells of revolution. These approximate functions are very powerful to capture the behavior of the structure with minimal degrees of freedom. In general, to solve the integral equations of motion , Eq.(20), the following separation of variables are

14

applicable u(s, θ, t)

=

u ˜(s, θ)eωmn t

v(s, θ, t)

=

v˜(s, θ)eωmn t

w(s, θ, t)

=

w(s, ˜ θ)eωmn t

φ1 (s, θ, t)

=

φ˜1 (s, θ)eωmn t

φ2 (s, θ, t)

=

φ˜2 (s, θ)eωmn t ,

(27)

where ωmn is the general eigenvalues of the system corresponding to the mode shape with n circumferential and m longitudinal half waves. In general ωmn can take a complex value while for undamped and conservative systems it can take a purely imaginary value corresponding to the natural frequencies of the system. The symmetric spatial trial functions are also defined as

u ˜(s, θ)

=

N Nc a   n=0

v˜(s, θ)

=

n=0

w(s, ˜ θ)

=

φ˜1 (s, θ)

=

n=0

φ˜2 (s, θ)

=

where λm =

mπ l

 bkn ξk (s) sin(nθ)

k=1

w ˆmn cos(λm s) +

2 

 ckn ξk (s) cos(nθ)

k=1

φˆ1mn cos(λm s) +

m=0

N Nc a  

n=0

vˆmn cos(λm s) +

2 

m=0

N Nc a  

akn ξk (s) cos(nθ)

k=1

m=0

N Nc a  

n=0

u ˆmn cos(λm s) +

m=0

N Nc a  



2 

2  k=1

φˆ2mn cos(λm s) +

m=0

2 

(28)

 dkn ξk (s) cos(nθ)  ekn ξk (s) sin(nθ),

k=1

and the auxiliary functions ξk (s) are any set of sufficiently smooth functions defined over

the shell’s length. In the present study they are defined as follows ξ1 (s)

=

ξ2 (s)

=

s s( − 1)2 l s2 s ( − 1), l l

(29) (30)

and u ˆmn , vˆmn , w ˆmn , φˆ1mn , φˆ2mn , akn , bkn , ckn , dkn , ekn are real constant Fourier coefficients to be determined later. In Eq. (28) Nc and Na are the total number of expansion terms in the circumferential and longitudinal directions, respectively. Substituting Eqs.(28-30) into Eq.(20), and following the Galerkin approach, the

15

circumferential Fourier expansion is identically satisfied using the following relations, ⎧ ⎨ 0 n=0 sin2 nθ dθ = ⎩ π n≥1 0 ⎧  2π ⎨ 2π n = 0 cos2 nθ dθ = ⎩ π n≥1 0 

2π

(31)

(32)

In other words, the governing equations hold for any value of n and hence a preallocation of the n value is required. Neglecting the external loading P for a free vibration analysis, an eigenvalue problem is obtained as

2 Ms )ζ = 0, (Ks + ωmn

(33)

where Ks and Ms are the stiffness and mass matrices, respectively. The two matrices are not necessarily symmetric for the current formulation. The vector of generalized coordinates is also denoted by ζ and is defined as ζT =

 u ˆmn

akn

vˆmn

bkn

w ˆmn

ckn

φˆ1mn

dkn

φˆ2mn

 ekn

.

(34)

Natural frequencies as well as their corresponding mode shapes are determined by solving the above eigenvalue problem with a generalized eigenvalue solver. 3.2. Solution of the Fluid Equation In the fluid domain, the finite element method is selected to be the solution method, which is actually piece-wise application of the Ritz method where the solution on each element is approximated to be linear/quadratic or even cubic depending on the shape functions used. The common property in both Ritz and finite element methods is to use the weak form of the governing equations. Weak form representation of the fluid equation, Eq.(21), can be obtained by using the Laplacian definition in the cylindrical coordinate and multiplying both sides of the fluid equation by a trial function, say Ψ, and integrating by parts in which Mach number differentiations cancel out identically. Hence, in the weak form formulation, all derivatives distribute evenly over the trial function and the field variable as,   Ω

  ˆ ˆ ˆ ˆ ˆ ˆ ∂Ψ ∂ Φ Ψ ∂ Φ 1 ∂ Φ ∂Ψ ∂ Φ 2M ∂ 2 Φ ∂Ψ Ψ ∂2Φ 2 + (1 − M ) − + + dΩ = Ψˆ + Ψ q dΓ, c2∞ ∂t2 c∞ ∂s∂t ∂s ∂s R ∂R R2 ∂θ ∂θ ∂R ∂R Γ

16

(35)

where Ω stands for the fluid domain and Γ represents the boundaries of the fluid domain. Moreover, qˆ describes the input/output flux on the boundaries of the shell. Employing a Fourier expansion in circumˆ trial function (Ψ) and flux (ˆ ferential direction for all the velocity potential (Φ), q ) variables as Ψ

=

NT cos(nθ)

(36)

ˆ Φ

=

NΦ cos(nθ)

(37)

q cos(nθ)

(38)

qˆ =

where q is the vector of fluid flux defined in Eq.(24) and N is the vector of shape functions of a quadrilateral  T with individual shape functions defined by Zienkiewicz 9-noded element being N = N 1 · · · N9 [51]. The final axisymmetric system of equations are obtained as ¨ + Cf Φ ˙ + Kf Φ = D, Mf Φ

(39)

where Mf , Cf and Kf are the mass, damping and stiffness matrices, respectively. Furthermore, D is the vector of external loadings. The final mass, damping and stiffness matrices as well as the vector of external loading are constructed by assembling the individual elemental sub-matrices after determining the following integrals over each element. Mf

(e)

Cf (e)

=

=

1 c2∞ 2 c∞



NT N dΩ(e)



Ω(e)

M NT

∂N dΩ(e) ∂s

(41)

 n 2 ∂NT ∂N NT ∂N ∂NT ∂N − + ( ) NT N + dΩ(e) ∂s ∂s R ∂R R ∂R ∂R

(42)

  Kf (e)

(1 − M 2 )

= Ω(e)



D(e)

=

(40)

Ω(e)

NT q dΓ(e) .

(43)

Γ(e)

It is clear that the size of stiffness, damping and mass matrices depend on the total number of degrees of freedom in the fluid domain. 3.3. Solution of the Coupled System of Equations In the previous section, it is explained that the fluid and structural dynamics are coupled via the flexible wall movements, which is associated as the impermeability condition. Employing the impermeability ˆ to condition of Eq.(24) in the right hand side of the fluid equation, Eq.(39), couples the fluid states (Φ) the structural states (ζ). In the next step, utilizing Bernoulli relation, Eq.(25), in the structural equations, 17

Eq.(20), gives a coupled formulation for the fluid and structural states. Both coupled fluid and structural equations have to be solved simultaneously to obtain the fluid/structural states at any instant of time. The final coupled system of equations is eventually as follows ⎫ ⎡ ⎤⎧ ⎨ Φ ¨ ⎬ Cf ⎦ +⎣ −CwΦ Ms ⎩ ζ¨ ⎭

⎡ ⎣

Mf

−CΦw

0

0

0

⎫ ⎡ ⎤⎧ ⎨ Φ ˙ ⎬ Kf ⎦ +⎣ ⎩ ζ˙ ⎭ −KwΦ

−KΦw Ks

⎫ ⎤⎧ ⎨ Φ ⎬ ⎦ = 0, ⎩ ζ ⎭

(44)

where the explicit definition of the coupling sub-matrices are as follows Φw

Ne  

=

C

=

C

ρf

Ne  

Ne  

=

K

e=1

KwΦ

=

ρf

and w ¯ ∈

 1

cos λ1 s

···

w ¯ T NR ds

(e)

(e)

s2

U NTf

(e) s1

Ne  

(46)

∂w ¯ ds ∂s

(e)

s2

(e) s1

e=1

(45)

(e)

s2 s1

e=1

Φw

NTf wds ¯

(e)

s1

e=1 wΦ

(e)

s2

Uw ¯T

(47)

∂N R ds, ∂s

(48)

 cos λM s

ξ1

. The variable Ne defines the number of elements on the

ξ2

(e)

fluid-shell interface. Since integrations are carried out over the boundary elements, variables s1

(e)

and s2

stand for the positions of the first and the last nodes in the corresponding element. In order to transform the above system of equations into a generalized eigenvalue problem format, state space transformation have to be utilized. This is due to the presence of damping matrix. Hence, the final coupled system of equations in the state space representation is as ⎡ ⎣

where y T =

 Φ

ζ

T

0

M∗

M∗

C∗

⎧ ⎫ ⎡ ⎨ ⎬ y ˙ −M∗ ⎦ d +⎣ dt ⎩ y ⎭ 0

⎫ ⎤⎧ ⎨ y˙ ⎬ ⎦ = 0, K∗ ⎩ y ⎭

⎤

0

(49)

and ⎤

⎡ M∗

=

⎣ ⎡

C∗

=

⎣ ⎡

K∗

=

⎣

Mf

0

0

Ms

⎦

Cf

−CΦw

−CwΦ

0

Kf

−KΦw

−KwΦ

Ks

18

⎤ ⎦

(50)

⎤ ⎦.

(51)

The above state space system is usually ill-conditioned due to considerable differences in the order of the fluid and structural sub-matrices. Because of the significant ill-conditionings present in the problem, it is not practical to use a double precision eigenvalue solvers. Actually ordinary solvers fail to deliver reliable results in cases where the condition number of a matrix is greater than 1e16. For such cases, a quadruple precision eigenvalue solver provided by Advanpix [52] is used to obtain the eigenvalues with no accuracy loss due to ill-conditioning issues in the MATLAB environment. 3.4. Reduced Order Modeling In order to increase the efficiency of the solution procedure, a reduced order model is developed along with the static condensation technique. In fact, for dense fluid containments such as water, the difference between eigenvalues of the fluid and structure are usually so large (with larger natural frequencies attributed to the fluid system). Hence; it can be implied that fluid participates in the whole response in a static manner. For other light fluid flows such as air, a specific number of modes (which are the first dynamic modes with smaller natural frequencies) contribute to the response in a dynamic way and higher modes can still be considered as static modes. Thus, the following expansion is proposed as the solution of the fluid system and comprises both static and dynamic modes, Φ = μ + Rϑ(t),

(52)

where R is the matrix of the dominant right eigenvectors of the fluid system and ϑ(t) is the generalized modal coordinate. The total static response of the fluid is denoted by μ and defined as follows, Φw ˙ μ = K−1 ζ + KΦw ζ). f (C

(53)

This relation shows that μ is a quasi-static quantity, since it depends on the time-derivatives of the structural states. Using Eq.(53) in Eq.(52), substituting the resultant relation into Eq.(44) and premultiplying the fluid equation by the left eigenvectors ,LT , of the fluid system, the following reduced order model is constructed ⎡ ⎣

ˆf Ms − M

0

¯ D

¯f M

⎫ ⎡ ⎤⎧ ⎨ ζ¨ ⎬ ˆf C ⎦ +⎣ ⎩ ϑ ¨ ⎭ ¯ E

¯ wΦ −C ¯f C

⎫ ⎡ ⎤⎧ ⎨ ζ˙ ⎬ ˆf Ks − K ⎦ +⎣ ⎩ ϑ˙ ⎭ 0

19

¯ wΦ −K ¯f K

⎫ ⎤⎧ ⎨ ζ ⎬ ⎦ =0 ⎩ ϑ ⎭

(54)

where ˆf M

=

Φw CwΦ K−1 f C

(55)

ˆf C

=

Φw Φw −(CwΦ K−1 + KwΦ K−1 ) f K f C

(56)

ˆf K

=

Φw KwΦ K−1 f K

(57)

V

=

Φw Φw −1 Mf K−1 (Ms − CwΦ K−1 ) f C f C

(58)

¯f M

=

LT (Mf R + VCwΦ R)

(59)

¯f C

=

LT (Cf R + VKwΦ R)

(60)

¯f K

=

LT Kf R

(61)

¯ wΦ C

=

CwΦ R

(62)

¯ wΦ K

=

KwΦ R   ¯ = LT V(CwΦ K−1 KΦw + KwΦ K−1 CΦw ) + Mf K−1 KΦw + Cf K−1 CΦw D f f f f   −1 −1 ¯ = LT Cf K KΦw − V(Ks − KwΦ K KΦw ) . E f f

(63) (64) (65)

A closer look at the reduced order system of equations given by Eq.(54) reveals that the number of degrees of freedom in the fluid domain is reduced to the number of modes considered in the reduced order formulation. A similar state space representation as mentioned in Eq.(49) is utilized to transfer the coupled reduced order system of equations given in Eq.(54) to a standard eigenvalue problem format. For incompressible and dense fluids inside the shell, the fluid is expected not to get involved in the total dynamic response of the coupled system, and hence, a minimum number of dynamic mode shapes of the fluid can be used. The special case, where the number of dynamic fluid mode shapes is taken to be zero in the reduced order model, corresponds to the case where the fluid equation in the coupled system (54) is nullified. Even in the structural equation of motion (first equation in the coupled system (54)) the inertia and damping coupling effects are nullified and the fluid affects the final equations in a static manner. This is equivalent to solving the structural equation along with the ordinary Laplacian equation (∇2 Φ = 0). For light and compressible fluids inside the shell, a dynamic response of the fluid is expected, and hence, the number of fluid mode shapes required for a reasonable convergence and reliable result would increase but it is still much less than the total degrees of freedom in the original system of equations without reduction. In what follows, the effect of number of dynamic mode shapes included in the reduced order model, on the degree of accuracy of the final results is evaluated.

20

4. Results and Discussion In order to obtain results and implement all previous discussions, a computer code is developed. The code is written in the MATLAB programming environment. In the first step, the code is verified with the available literature for different geometries and boundary conditions. As mentioned earlier, different boundary conditions are addressed by a proper justification of the spring stiffnesses at the extremes of the shell. A simple study shows that the following spring stiffnesses given in Table 1 are appropriate to simulate each set of boundary conditions [50]. Table 1: Boundary condition definitions and the corresponding proper stiffness values, k∞ = 1e15.

Boundary Condition Clamped (C) Simply Suported (SS) Free (F)

Definition u = v = w = φ1 = φ2 = 0 v=w=0 No Constraint

ku k∞ 0 0

kv k∞ k∞ 0

kw k∞ k∞ 0

kφ1 k∞ 0 0

kφ2 k∞ 0 0

An important task prior to obtaining general solutions to the problem and parameter studies, is to gain a comprehensive insight over the convergence rates. Some parameters such as the number of expansion terms in the Galerkin solution, the degrees of freedom (number of elements) in the fluid domain and ultimately the minimum number of dynamic mode shapes required in the reduced order model have to be evaluated. Convergence rates have to be determined for different boundary conditions, semi-vertex angles and circumferential wave numbers to provide a solid confidence over the convergence of the results in the following sections. Having studied the convergence rates and validity of the results, we’re then allowed to investigate the effects of boundary conditions, semi-vertex angles and fluid velocities on the natural frequencies and stability margins of the fluid-loaded shells. 4.1. Convergence Studies In the first place, a convergence study of the required number of expansion terms (Na ) in the approximate solution, for different boundary conditions is provided to demonstrate the capability of the current solution method. Natural frequencies (Hz) of conical/cylindrical shells are given in Table 2 for different boundary conditions, semi-vertex angles, fluid velocities and number of expansion terms (Na ). The bold frequencies are those where the second digit does not change anymore as the number of expansion terms (Na ) increases. Considering the stationary fluid containments at all presented boundary conditions, the minimum Na value increases as the semi-vertex angle increases. For the flowing fluids, the trend is not the same as the stationary fluids. In such cases as the fluid velocity increases, the minimum number of expansion terms (Na ) also increases culminating to U = 100m/s where an approximate solution with 40 expansion terms is sufficient for a reasonable convergence up to two digits. Other velocities greater than this value have been also checked as special cases, where their critical velocity is above 100m/s and similar results are obtained. In the following simulations of current study, the number of expansion terms is taken as a fixed value of Na = 40. 21

In the second convergence study, we will determine the minimum degrees of freedom (DOF) required for the fluid domain solutions. As mentioned earlier, the axisymmetric fluid domain is discretized by the quadrilateral 9 - noded elements. Since the meshing algorithm is supposed to be a structured one, we have to define the number of elements on any of the input and wall boundaries of the shell. The number of degrees of freedom is then obtained for any generic mesh. Table 3 is provided to demonstrate the sensitivity of the natural frequencies (Hz) of a conical shell with α = 30 to the number of DOF, fluid velocity and boundary conditions. Results show that for most of the cases, even a mesh of 900 elements, being equivalent to 3731 DOFs, is sufficient for a reasonable convergence. However, it is possible to have some other cases where the second digits of the natural frequencies are not fixed yet. Therefore the minimum number of elements is determined to be 1200 (4961 DOFs), which ensures the convergence of the coupled natural frequencies up to two digits. As discussed, the minimum number of expansion terms in the structural equations and DOFs in the fluid equations are now determined to be Na = 40 and DOF = 4961, respectively. According to the expansions in Eq. (28), the total degrees of freedom in the structural part is determined by using the following relation, 5 × (Na + 3). Also, the total number of degrees of freedom of the coupled system in the state space form becomes 10352 which is high. In order to minimize total degrees of freedom, the so called reduced order formulation is proposed. Results provided in Fig. (2) are aimed to show the convergence of the natural frequencies (Hz) of the coupled system at different fluid velocities (stationary and flowing) for different boundary conditions. As depicted in Fig. (2), reduced model including only the first five dynamic modes (eigenfunctions) of the fluid equations, gives accurate results. Here, since the fluid containment is an incompressible fluid (i.e.,water), using a reduced order model even with zero modes will alter the results by less than 0.5%, which is totally negligible. For the current study, the number of reduced order modes is taken to be ROM = 5. This concludes a final system composed of 440 equations in the state space form with a dimension of 440 × 440 compared to the original system composed of 10352 equations. From a memory management perspective, the dense square matrix in the original representation, requires memory allocation of about 817.6MB in double precision floating points, while for the reduced order system (with much smaller matrices), a memory allocation of about 1.48MB is required. The ratio of the memory allocation in the original system to the reduced system is 553, while the accuracy of the two models are almost identical. In quadruple precision, which is utilized in the current study, the amount of memory allocation is about 31 times more than double precision in the MATLAB environment. This means that for a reduced order system with a 440 × 440 dense square matrix, a memory of 45.9MB is required while for the original system a memory of about 24.7GB is needed.

22

4.2. Verifications of the Results Since most of the available literature on shells containing fluid flows are the ones discussing cylindrical shells, we will first compare the results of the present study with those of cylindrical shells. As reported in [9, 53], a steel cylindrical shell with

L R

=2,

h R

= 0.01 , E = 206 GPa , ρ = 7850 kg/m3 , ν = 0.3 is considered.

The internal fluid is assumed to be water with a density of ρf = 1000 kg/m3 . Four sets of boundary conditions ¯ which are namely C-C, S-S, C-S and S-C conditions are investigated. The dimensionless frequencies (Λ) ¯ ) are clearly reported in Figs. (3), where the dimensionless frequency versus the dimensionless velocities (U ¯ = U/ω0 l and Λ ¯ = Λ/ω0 , respectively and the denominator of the dimensionless and velocity are defined as U  2 2 Eh relations are defined as ω0 = πl2 12ρ(1−ν 2 ) . It can be observed from the results given in Figs. (3), that results of the present study agree with the results available in literature [9, 53]. Slight differences obtained in the case of C-C boundary condition have been also reported in previous studies such as [53]. It’s claimed in the study by Ugurlu et. al. [53] that such differences observed in other studies are mainly due to the lack of convergence at higher fluid velocities. In this study, the effect of fluid velocity on the convergence of the results is taken into account. In the second stage of verification, conical and cylindrical shells made of steel are considered in a way that the mean radii of conical and cylindrical shells are identical and is equal to a = 12 (Rmin + Rmax ) = 0.876m. Hence the input and output radii of conical shells at specific semi-vertex angles are varied in a way that the mean radius remains intact. Other geometrical and mechanical parameters are l cos α = 0.9144m, h = 1.5mm , ρ = 7850 kg/m3 , E = 210GPa , ν = 0.3. We will use these physical dimensions for a variety of semi-vertex angles ranging from zero to 45 degrees throughout the current study. Such geometrical dimensions are basically selected for comparison purposes with the available literature. For this part of  study, a dimensionless parameter is defined as Ω = ωω0 × 100 where ω0 = (1/R1 ) E/ρs . Moreover it is worth noting that, throughout this study, the critical velocities are reported in the middle section of the cone and not at the input. In this way, it is possible to compare the critical velocity variations at a variety of semi-vertex angles with those of cylindrical shells. Figures (4,5,6) present a through comparison of variations of the dimensionless frequencies with those of Bochkarev et. al. [40] and Kumar et. al. [38] at three different boundary conditions of C-C, S-S and C-F. For each case, the circumferential wave numbers with available data in literature are selected. The corresponding variations of the real part of the eigenvalues (damping) are also plotted for each circumferential wave number. The instability conditions are those corresponding to the first appearance of a positive damping value. Both divergence-type and flutter-type instabilities are also depicted in Figs. (4,5,6). As can be seen from the figures, current results are in agreement with those of Bochkarev et. al. [40] while considerable differences are present between the results of the current study and those of Kumar et. al. [38], especially for the case of C-F conical shells. It should be noted that, similar discrepancies are present in the results reported by Bochkarev et. al.[40]. Such differences are found to be due to: 23

1. The boundary conditions imposed on the fluid velocity potential equation is the zero flux condition (fbc00) in the study by Kumar et. al. [38] while in the present study and the study by Bochkarev et. al. [40], a zero perturbed pressure (opening condition, fbc11) is imposed. 2. Special theory used by Kumar et. al. [38] may have some shortcomings in determining the critical fluid velocity at lower and higher circumferential wave numbers. 3. The differences between present results and those of Kumar et. al. [38] augments as the fluid velocity increases. There is a possibility that their results are not converged for higher velocities. This is a common error as indicated by Ugurlu et. al. [53]. 4. The data presented in Figs. (4,5,6) are obtained via a digitization process of the available figures in [38, 40] and small digitization errors are inevitable. 4.3. Stability of Conical/Cylindrical Shells Internal fluid flows (i.e, air, water, ...) become substantially important in cases where the fluid loadings exerted on the structure are large enough (compared to the overall stiffness of the structure) to cause instabilities. If the internal fluid is a light one, such as air, there would be often no concern regarding the stability issues in the subsonic flow regime. The only exception can be shells made of elastomeric material where the stiffness of the structure is negligible and a small aerodynamic loading can raise stability issues. In the current study, we will confine ourselves to the case of water filled metal shells for the stability investigations, where the hydrodynamic loading is sufficient to make ordinary shell structures unstable. The geometrical configuration and mechanical specifications are identical to the previous section which are also reported by Kumar et. al. and Bochkarev at. al. [38, 40]. In order to detect the instability occurrences, one has to investigate variations of the real and imaginary components of the eigenvalues, as the fluid velocity increases. According to the solution expansion expressed in Eq.(27), real and imaginary components of the eigenvalues are, respectively, attributed to the damping and frequency of the coupled system. In systems with complex eigenvalues, simultaneous investigation of both real and imaginary components is required. The point where a real component of the eigenvalue becomes positive signals an instability condition. If the corresponding imaginary component is zero, the instability is regarded as static or divergence; on the other hand, if the corresponding imaginary component is positive, it represents a dynamic type or flutter instability. In the present study, the post-divergence and post-flutter dynamics are disregarded and we aim to determine only the stability onset. In the study of Kumar et. al. [38], only variations of the frequencies with respect to the fluid velocity are considered and the velocity corresponding to the first occurrence of zero frequency is associated with the critical velocity. According to previous discussions, results reported in Kumar et. al. [38] cannot be accurate for most cases, since flutter type instabilities are totally missed. Even for the divergent-type instabilities, the point corresponding to immediate nullification of frequency at a specific velocity can not be attributed to divergence. This is due 24

to the fact that for many cases, frequency continues to take a zero value for larger fluid velocities and the real component of the eigenvalue gets positive afterwards. Following the stability analysis scheme, a comprehensive study is conducted to determine the stability margins at eight structural boundary conditions as well as two fluid boundary conditions of fbc11 and fbc10, for both cylindrical and conical shells. Results are summarized in Tables (4-11). The detailed analysis are not reported here for the sake of brevity. Actually, for each boundary condition and at a specific circumferential wave number and semi-vertex angle, variations of the frequency and damping versus the fluid velocity are obtained with small-enough velocity increments. Critical fluid velocities are then searched among the first 10 longitudinal wave numbers m = 1, 2, ..., 10. The minimum fluid velocity where at least one damping is positive is then taken as the critical velocity for that specific circumferential wave number and semi-vertex angle. Note that damping variations occur at significant rates and in case that velocity increments are not small enough, some of the instability modes may be missed. Results obtained agrees with those reported by Bochkarev et. al. [40], while agreements with the values reported by Kumar et. al. [38] are less, especially, at lower and higher circumferential wave numbers. This is mainly due to the fact that in Kumar. et. al. [38] only the frequency variations are taken into account to determine the instability conditions which introduces both quantitative and qualitative errors in the analysis. Moreover, the theory implemented by Kumar et. al. [38] may be subject to over/under estimations of the critical velocities. The same disagreements are also reported by Bochkarev et. al. [40]. The fluid boundary condition is clearly affecting the instability margins at all circumferential wave numbers both qualitatively and quantitatively. Although predictions of the two fluid boundary conditions are different for each circumferential wave number, the qualification of instability is not different for the lowest critical velocities at each semi-vertex angle and structural boundary condition. The bold and underlined values in Tables (4-11) denote the lowest of the critical velocities which can be called the overall onset of instability for the corresponding configuration. A brief investigation of the critical velocities of cylindrical shells reported in Tables (4-11) shows that, the corresponding critical velocities for fbc11 and at the following pairs of boundary conditions are identical: C-F, F-C ; C-S,S-C and S-F,F-S. This is an expected result due to the symmetry present in cylindrical shells. Furthermore, except for the case of S-F and F-S, the critical conditions occur in the form of single or coupled mode flutter at least for the first few circumferential wave numbers. This is not the case for the fbc10 boundary conditions due to the asymmetric fluid boundary condition affecting the final outcomes of the analysis. According to Table 9 for S-F boundary condition, the overall stability behavior of conical/cylindrical shells with fbc10 are more reasonable due to the more qualitatively acceptable results. For other structural boundary conditions the quantification differences in the fbc10 and fbc11 conditions are clear at all circumferential wave numbers. The common characteristic between the overall onset of instability of fbc10 and fbc11 is that they are all static (divergent) type instabilities. Although, fbc10 is a more reasonable choice according 25

to the reasonings provided in the previous sections (based on [46, 3]), we have found that this condition is incompatible with short shells at C-F structural boundary conditions because no instability condition can be obtained even with very small velocity increments. If longer shells are considered at such a structural and fluid boundary condition, the instabilities are detected with no specific difficulty. As the semi-vertex angle increases the shell becomes asymmetric (since the input and output radii are not the same anymore), and hence, no identical results are shown up at pairs of mixed-boundary conditions for fbc11. On the other hand, as the semi-vertex angle increases, the corresponding critical velocity at a specific boundary condition usually experiences a decrease. This implies that regardless of the boundary conditions imposed, conical shells with larger vertex angles are much more prone to stability issues at lower velocities. Furthermore, provided results reveal that C-C boundary condition delivers the highest critical velocities for both fbc11 and fbc10. This can be also inferred from the fact that required energy to trigger an instability situation in the cylindrical/conical shells with C-C boundary conditions are the highest compared to any other set of boundary conditions. In fact, there is a special order for the critical fluid velocities at various boundary conditions for each of fbc11 and fbc10 at all semi-vertex angles. The critical fluid velocities are cr. cr. cr. cr. cr. cr. in the following order for fbc11: UC−C > UC−S > US−C > US−S > UC−F > UFcr.−C > UFcr.−S > US−F while cr. cr. cr. cr. cr. the order is different for fbc10 as UC−C > US−F > US−C > UC−S > US−S > UFcr.−C > UFcr.−S . Note that for

the case of cylindrical shells with fbc11, pairs of mixed boundary conditions deliver identical results. For most sets of boundary conditions, semi-vertex angles and circumferential wave numbers, instabilities occur in the form of divergence of the first mode. C-C boundary condition, in contrary, suggests the most susceptible condition to flutter instabilities compared to other cases. The reason may be attributed to the higher stiffness of the overall structure with C-C boundary conditions, which is equivalent to the instability conditions due to higher fluid velocities. In fact, for higher velocities real part of the eigenvalues may grow larger, which helps to physically avoid the occurrence of divergence-type instabilities. 4.4. Semi-vertex angle effect The effect of semi-vertex angle is implicitly studied in the results provided in Tables (4-11). Based on those results, as the semi-vertex angle increases, the critical fluid velocity experiences a decrease regardless of the boundary conditions. It may be judged that concluding only based on the tabular data provided for three sets of semi-vertex angles is not acceptable. Therefore, results provided in Fig. 7 are given to demonstrate the direct and actual effects of semi-vertex angles on the critical fluid velocities at some boundary conditions and circumferential wave numbers. Note that the actual critical fluid velocity of a configuration (i.e, a specific semi-vertex angle and boundary condition) is determined by the circumferential wave number with the minimal critical fluid velocity. As shown, variations of the critical fluid velocities with respect to semi-vertex angles are very unpredictable for most circumferential wave numbers at fbc11, however, it can be concluded that as the semi-vertex angle increases, the actual critical velocity (the minimal 26

critical velocity) experiences a decrease. For fbc10 condition, variations of the critical fluid velocity versus the semi-vertex angles are more reasonable and show an overall reduction in the critical velocity as the semi-vertex angle increases. According to the presented results in Fig. 7, the overall stiffness (critical velocity) can change significantly as the semi-vertex angle increases and one can optimize a design in an aeroelastic/hydroelastic sense to get the largest critical velocities. 4.5. Compressibility effects It is claimed that the present formulation is much more accurate than other solutions provided in literature to study stability problems of conical shells with internal compressible flows. The main advantage is that the mean flow velocity is not predicted by a simple approximation of the conservation of mass, instead an isentropic flow relation is incorporated to satisfy both conservation of mass and density variations due to compressibility effects. In deed, compressibility effects are expected not to get involved in the analysis of shells conveying dense (incompressible or nearly incompressible) fluids such as water. The main discrepancy between the compressible and incompressible models is revealed when light fluid containments, such as air, are taken into account in the subsonic flow regime. Actually, shells containing light fluids, are not much susceptible to stability issues due to negligible aeroelastic loadings exerted on the structure. In order to investigate the effects of compressibility on the dynamics of air loaded shells, several results are provided in Tables (12-17) . Tabular data are provided in order to illustrate the discrepancies more precisely. Based on the results, compressibility effects have significant impacts on the air-loaded natural frequencies at the lower circumferential wave numbers. As the circumferential wave number increases, both compressible and incompressible models present identical results. This holds both in cylindrical and conical shells. The rate of convergence of compressible and incompressible models depend upon the boundary condition. This finding can have significant effects on situations where a low-modulus conical shell, say an elastomeric shell, is prone to instability issues in its lower circumferential wave numbers. The analyst/designer has to study for such effects prior to any decision. As expected, natural frequencies of the coupled system decreases with increasing the air velocity. Changes are negligible due to the low air loadings in the subsonic regime (M < 0.85) where current model is valid.

5. Conclusion In the current study, dynamics and stability of conical and cylindrical shells conveying compressible/incompressible fluid mediums are studied. The structural equations of motion are obtained in a general form using the Flugge shell theory. A flexible and fast-convergent solution method is proposed to obtain free/forced vibration response of conical/cylindrical shells at different structural and fluid boundary conditions. Moreover, a novel reduced order formulation is proposed to study the dynamic behavior of shells with 27

internal fluid flows. Convergence studies are conducted to determine the minimal degrees of freedom in the coupled system of equations. Verifications of results with those of cylindrical and conical shells conveying internal fluid flows are then performed and reasonable agreements are obtained and the few disagreements are justified by logical reasonings. Regarding the stability margins, a comprehensive study is then performed to determine the critical fluid velocity at eight combinations of structural boundary conditions and two fluid boundary conditions for different semi-vertex angles and circumferential wave numbers. Such comprehensive parameter study is reported for the first time for conical/cylindrical shells with internal fluid flows. A variety of interesting results are then reported and the fbc10 fluid boundary condition is found to be more realistic. It is also shown that semi-vertex angles can play an important role in stabilizing/destabilizing the system of coupled fluid-structure. In fact, as the semi-vertex angle increases, stability margins are usually diminished and the critical fluid velocity is reduced. Effects of compressibility are also studied to show both the capabilities of the current solution methodology and to outline the discrepancies between compressible and incompressible models when light fluid containments are employed.

Acknowledgement Authors would like to gratefully acknowledge the support provided by the TUBITAK (The Scientific and Technological Research Council of Turkey), Grant no: 2216 Research Scholarship Program for International Researchers.

28

References [1] M. Paidoussis, J.-P. Denise, Flutter of thin cylindrical shells conveying fluid, Journal of Sound and Vibration 20 (1) (1972) 9 – 26. [2] D. S. Weaver, T. E. Unny, On the Dynamic Stability of Fluid-Conveying Pipes, Journal of Applied Mechanics 40 (1) (1973) 48–52. [3] M. P. Paidoussis, Flutter of Conservative Systems of Pipes Conveying Incompressible Fluid, Journal of Mechanical Engineering Science 17 (1) (1975) 19–25. [4] G. C. Everstine, A symmetric potential formulation for fluid-structure interaction, Journal of Sound Vibration 79 (1981) 157–160. [5] A. Selmane, A. Lakis, Vibration Analysis of Anisotropic Open Cylindrical Shells Subjected To a Flowing Fluid, Journal of Fluids and Structures 11 (1) (1997) 111 – 134. [6] J. Kochupillai, N. Ganesan, C. Padmanabhan, A semi-analytical coupled finite element formulation for shells conveying fluids, Computers and Structures 80 (34) (2002) 271 – 286. [7] J. Kochupillai, N. Ganesan, C. Padmanabhan, A Semi-Analytical Coupled Finite Element Formulation for Composite Shells Conveying Fluids, Journal of Sound and Vibration 258 (2) (2002) 287 – 307. [8] M. Amabili, R. Garziera, Vibrations of Circular Cylindrical Shells with Nonuniform Constraints, Elastic Bed and Added Mass; Part I: Empty and Fluid-Filled Shells, Journal of Fluids and Structures 14 (5) (2000) 669 – 690. [9] M. Amabili, R. Garziera, Vibrations of Circular Cylindrical Shells with Nonuniform Constraints, Elastic Bed and Added Mass. Part II: Shells Containing or Immersed In Axial Flow, Journal of Fluids and Structures 16 (1) (2002) 31 – 51. [10] M. Amabili, R. Garziera, Vibrations of Circular Cylindrical Shells with Nonuniform Constraints, Elastic Bed and Added Mass. Part III: Steady Viscous Effects on Shells Conveying Fluid, Journal of Fluids and Structures 16 (6) (2002) 795 – 809. [11] M. Amabili, F. Pellicano, M. Paidoussis, Non-Linear Dynamics and Stability of Circular Cylindrical Shells Containing Flowing Fluid. Part I: Stability, Journal of Sound and Vibration 225 (4) (1999) 655 – 699. [12] M. Amabili, F. Pellicano, M. Paidoussis, Non-Linear Dynamics and Stability of Circular Cylindrical Shells Containing Flowing Fluid, Part II: Large-Amplitude Vibrations without Flow, Journal of Sound and Vibration 228 (5) (1999) 1103 – 1124. [13] M. Amabili, F. Pellicano, M. Paidoussis, Non-Linear Dynamics and Stability of Circular Cylindrical Shells Containing Flowing Fluid. Part III: Truncation Effect without Flow and Experiments, Journal of Sound and Vibration 237 (4) (2000) 617 – 640. [14] M. Amabili, F. Pellicano, M. Paidoussis, Non-Linear Dynamics and Stability of Circular Cylindrical Shells Containing Flowing Fluid. Part IV: Large-Amplitude Vibrations with Flow, Journal of Sound and Vibration 237 (4) (2000) 641 – 666. [15] M. P. Paidoussis, Slender Structures and Axial Flow, Volume 1, Fluid-Structure Interactions, Academic Press, ISBN 978-0-12-544360-9, URL http://www.sciencedirect.com/science/bookseries/18745652/1, 1998. [16] M. P. Paidoussis, Slender Structures and Axial Flow, Volume 2, Fluid-Structure Interactions, Academic Press, ISBN 978-0-12-544361-6, URL http://www.sciencedirect.com/science/bookseries/18745652/2, 2003. [17] H. Garnet, J. Kempner, Axisymmetric Free Vibrations of Conical Shells, Journal of Applied Mechanics 31 (1964) 458–466. [18] T. Irie, G. Yamada, Y. Kaneko, Free vibration of a conical shell with variable thickness, Journal of Sound and Vibration 82 (1) (1982) 83 – 94. [19] L. Tong, Free vibration of orthotropic conical shells, International Journal of Engineering Science 31 (5) (1993) 719 – 733. [20] X. Hu, T. Sakiyama, H. Matsuda, C. Morita, Vibration analysis of twisted conical shells with tapered thickness, International Journal of Engineering Science 40 (14) (2002) 1579 – 1598.

29

[21] C. Shu, An efficient approach for free vibration analysis of conical shells, International Journal of Mechanical Sciences 38 (89) (1996) 935 – 949. [22] K. Liew, T. Ng, X. Zhao, Free vibration analysis of conical shells via the element-free kp-Ritz method, Journal of Sound and Vibration 281 (35) (2005) 627 – 645. [23] O. Civalek, Vibration analysis of laminated composite conical shells by the method of discrete singular convolution based on the shear deformation theory, Composites Part B: Engineering 45 (1) (2013) 1001 – 1009. [24] L. Tong, Free vibration of composite laminated conical shells, International Journal of Mechanical Sciences 35 (1) (1993) 47–61. [25] G. Jin, Z. Su, T. Ye, X. Jia, Three-dimensional vibration analysis of isotropic and orthotropic conical shells with elastic boundary restraints, International Journal of Mechanical Sciences 89 (2014) 207 – 221. [26] S. C. Dixon, M. L. Hudson, Flutter, Vibration, and Buckling of Truncated Orthotropic Conical Shells with Generalized Elastic Edge Restraint, Report, 1970. [27] K. S. Ueda, T., M. Kihira, Supersonic Flutter of Truncated Conical Shells, Transactions of the Japan Society for Aeronautical and Space Sciences 20 (1977) 13–30. [28] M. N. Bismarck-Nasr, H. R. C. Savio, Finite-Element Solution of the Supersonic Flutter of Conical Shells, AIAA Journal 17 (10) (1979) 1148–1150. [29] P. Sunder, C. Ramakrishnan, S. Sengupta, Optimum cone angles in aeroelastic flutter, Computers and Structures 17 (1) (1983) 25 – 29. [30] D. R. Mason, P. T. Blotter, Finite-element application to rocket nozzle aeroelasticity, Journal of Propulsion and Power 2 (6) (1986) 499–507. [31] A. Lakis, P. V. Dyke, H. Ouriche, Dynamic analysis of anisotropic fluid-filled conical shells, Journal of Fluids and Structures 6 (1992) 135–162. [32] V. Aleksandrov, S. Grishin, Dynamics of a conical shell containing supersonic gas flow, Journal of Applied Mathematics and Mechanics 58 (4) (1994) 703–712. [33] K.-H. Jeong, T.-W. Kim, K.-B. Park, M.-H. Chang, Hydroelastic vibration of fluid filled conical shell, 1999. [34] M. J. Jhung, J. C. Jo, K. H. Jeong, Modal Analysis of Conical shell filled with fluid, Journal of Mechanical Science and Technology 20 (11) (2006) 1848–1862. [35] M. Caresta, N. J. Kessissoglou, Vibration of fluid loaded conical shells, Journal of acoustical society of America 124 (4) (2008) 2068 – 2077. [36] F. Sabri, A. A. Lakis, Hybrid finite element method applied to supersonic flutter of an empty or partially liquid-filled truncated conical shell, Journal of Sound and Vibration 329 (3) (2010) 302–316. [37] S. Mahmoudkhani, H. Haddadpour, H. Navazi, Supersonic flutter prediction of functionally graded conical shells, Composite Structures 92 (2010) 377386. [38] D. S. Kumar, N. Ganesan, Dynamic analysis of conical shells conveying fluid, Journal of Sound and Vibration 310 (12) (2008) 38 – 57. [39] Y. Kerboua, A. Lakis, M. Hmila, Vibration analysis of truncated conical shells subjected to flowing fluid, Applied Mathematical Modelling 34 (2010) 791–809. [40] S. Bochkarev, V. Matveenko, Natural vibrations and stability of shells of revolution interacting with an internal fluid flow, Journal of Sound and Vibration 330 (13) (2011) 3084 – 3101. [41] A. Leissa, Vibration of Shells, NASA, Washington, United States, 1973. [42] R. D. Firouz-Abadi, M. Rahmanian, M. Amabili, Free Vibration of Moderately Thick Conical Shells Using a Higher Order Shear Deformable Theory, Journal of Vibration and Acoustics 136 (5) (2014) 051001. [43] J. Reddy, Theory and Analysis of Elastic Plates and Shells, Second Edition, CRC Press, Taylor and Francis Group, 2007.

30

[44] R. Firouzabadi, M. Noorian, H. Haddadpour, A fluid-structure interaction model for stability analysis of shells conveying fluid, Journal of fluids and structures 26 (2010) 747–763. [45] J. Anderson, Modern Compressible Flow: With Historical Perspective, third edition, McGraw-Hill Global Education, University of Maryland-College Park, 2003. [46] L. Shayo, C. Ellen, Theoretical studies of internal flow-induced instabilities of cantilever pipes, Journal of Sound and Vibration 56 (4) (1978) 463 – 474. [47] M. Padoussis, V. Nguyen, A. Misra, A theoretical study of the stability of cantilevered coaxial cylindrical shells conveying fluid, Journal of Fluids and Structures 5 (2) (1991) 127 – 164. [48] Y. Qu, X. Long, G. Yuan, G. Meng, A unified formulation for vibration analysis of functionally graded shells of revolution with arbitrary boundary conditions, Composites Part B: Engineering 50 (2013) 381–402. [49] Y. Kurylov;, M.Amabili, Polynomial versus trigonometric expansions for nonlinear vibrations of circular cylindrical shells with different boundary conditions, Journal of Sound and Vibration 329 (2010) 1435 – 1449. [50] Z. Su, G. Jin, S. Shi, T. Ye, X. Jia, A unified solution for vibration analysis of functionally graded cylindrical, conical shells and annular plates with general boundary conditions, International Journal of Mechanical Sciences 80 (2014) 62–80. [51] O. C. Zienkiewicz, R. L. Taylor, J. Zhu, The Finite Element Method: Its Basis and Fundamentals, Seventh Edition, Elsevier Ltd., ISBN 978-1-85617-633-0, URL http://www.sciencedirect.com/science/book/9781856176330, 2013. [52] L. Advanpix, Multiprecision Computing Toolbox for MATLAB V.3.9.1.10015, URL http://www.advanpix.com, 2015. [53] B. Ugurlu, A. Ergin, A hydroelasticity method for vibrating structures containing and/or submerged in flowing fluid, Journal of Sound and Vibration 290 (35) (2006) 572 – 596.

31

v

w u

h x z s

k G1,k G2 r

y R U

B R

k u ,k v ,k w

l

Figure 1: Schematic of a conical shell and the corresponding annotations

Table 2: Convergence of the natural frequencies (Hz) of conical shells with internal water at different geometries and boundary conditions, (n, m) = (1, 1), l cos α = 0.9144m; a = 0.876m; h = 1.5mm, Q9 elements, Mesh = 4961DOF, fbc11. α

U (m/s)

Na = 10

20

30

40

50

0

0 25 50 75 100

123.02 122.08 119.22 114.24 105.69

122.93 121.99 119.14 114.16 98.78

122.91 121.98 119.13 114.14 95.84

122.91 121.98 119.12 114.14 95.05

122.91 121.98 119.12 114.14 95.05

30

0 25 50 75 100

71.58 71.43 70.99 70.23 69.17

71.42 71.27 70.83 70.09 69.04

71.39 71.24 70.80 70.07 69.02

71.37 71.23 70.80 70.06 69.01

71.37 71.23 70.80 70.06 69.01

0

0 25 50 75 100

158.59 156.90 151.72 142.68 119.60

158.59 156.90 151.70 142.65 106.84

158.59 156.90 151.70 142.64 106.78

158.59 156.90 151.70 142.64 106.76

158.59 156.90 151.70 142.64 106.76

30

0 25 50 75 100

114.27 114.02 113.26 111.98 110.13

114.22 113.97 113.21 111.93 110.08

114.22 113.97 113.21 111.93 110.08

114.21 113.96 113.21 111.93 110.07

114.21 113.96 113.21 111.93 110.07

0

0 25 50 75 100

160.58 158.86 153.61 144.39 123.36

160.51 158.76 153.39 144.00 107.11

160.51 158.75 153.35 143.94 107.01

160.51 158.74 153.34 143.93 106.99

160.51 158.74 153.34 143.93 106.99

30

0 25 50 75 100

131.21 130.81 129.61 127.58 124.69

131.15 130.74 129.52 127.44 124.46

131.14 130.73 129.50 127.41 124.41

131.14 130.73 129.49 127.40 124.40

131.14 130.73 129.49 127.40 124.40

C-F

S-S

C-C

32

Table 3: Convergence of the natural frequencies (Hz) of conical shells with internal water at different boundary conditions, (n, m) = (1, 1), l cos α = 0.9144m; a = 0.876m; h = 1.5mm, Q9 elements, Na = 40, α = 30, fbc11. U (m/s)

Mesh = 1891

2821

3731

4961

7381

C-C

0 25 50 75 100

131.14 130.73 129.50 127.40 124.40

131.14 130.73 129.50 127.40 124.40

131.14 130.73 129.49 127.40 124.40

131.14 130.73 129.49 127.40 124.40

131.14 130.73 129.49 127.40 124.40

S-S

0 25 50 75 100

114.22 113.97 113.22 111.93 110.08

114.22 113.97 113.21 111.93 110.08

114.22 113.97 113.21 111.93 110.08

114.22 113.97 113.21 111.93 110.07

114.22 113.97 113.21 111.93 110.07

C-F

0 25 50 75 100

71.38 71.23 70.80 70.06 69.01

71.38 71.23 70.80 70.06 69.01

71.38 71.23 70.80 70.06 69.01

71.38 71.23 70.80 70.06 69.01

71.38 71.23 70.80 70.06 69.01

(a) U = 0

(b) U = 100m/s

Figure 2: Natural frequencies (Hz) of a conical shell with internal water versus the required mode shapes in the ROM ; (n, m) = (1, 1), l cos α = 0.9144m, a = 0.876m , h = 1.5mm , Na = 40, Q9 elements , Mesh = 4961DOF, α = 30◦ , fbc11.

33

(a) S-S

(b) S-S

(c) C-C

(d) C-C

(e) C-S

(f) S-C

¯ variations versus the dimensionless fluid velocity (U ¯ ) for a cylindrical shell with Figure 3: Dimensionless frequency (Λ) 2;

h R

= 0.01 ; E = 206GPa ; ρ = 7850

kg ; m3

ρf = 1000

kg ; m3

ν = 0.3 at different boundary conditions.

34

L R

=

(a) n = 12 , α = 10◦

(b) n = 11 , α = 30◦

Figure 4: Comparison of the dimensionless frequencies (Ω) versus the fluid velocity for Clamped - Clamped conical shells l cos α = 0.9144m, a = 0.876m , h = 1.5mm , Na = 40, Q9 elements , Mesh = 4961DOF, fbc11.

35

(a) n = 4 , α = 10◦

(b) n = 7 , α = 30◦

Figure 5: Comparison of the dimensionless frequencies (Ω) versus the fluid velocity for Simply supported - Simply supported conical shells l cos α = 0.9144m, a = 0.876m , h = 1.5mm , Na = 40, Q9 elements , Mesh = 4961DOF, fbc11.

36

(a) n = 8 , α = 10◦

(b) n = 15 , α = 30◦

Figure 6: Comparison of the dimensionless frequencies (Ω) versus the fluid velocity for Clamped - Free conical shells l cos α = 0.9144m, a = 0.876m , h = 1.5mm , Na = 40, Q9 elements , Mesh = 4961DOF, fbc11.

37

38

130.84

130.27 129.77

98.25

75.46 70.71

79.59

147.67

146.84 145.60

87.11

68.86 65.54

82.55

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15

Ref. [40]

Ref. [38]

n 120.51(F) 120.45(F) 120.30(F) 120.06(F) 120.57(D) 115.11(D) 97.74(D) 84.33(D) 75.45(D) 71.01(D) 70.74(D) 74.19(D) 80.49(D) 80.49(D) 81.06(D)

fbc11

α=0

120.25(F) 119.25(D) 119.25(D) 119.00(D) 118.75(D) 109.50(D) 92.00(D) 79.50(D) 71.25(D) 67.50(D) 67.50(D) 71.25(D) 77.00(D) 78.00(D) 79.00(D)

fbc10

133.15 131.90 119.54 97.62 81.30 70.43 63.87 61.39 62.21 67.10 73.74 76.24 78.81

133.98

Ref. [38]

77.52

72.73 68.59

94.05

122.27 121.60

122.95

Ref. [40] 111.25(D) 111.18(D) 111.00(D) 110.73(D) 110.29(D) 108.00(D) 93.48(D) 80.92(D) 72.72(D) 68.89(D) 69.17(D) 72.94(D) 77.93(D) 78.01(D) 78.97(D)

fbc11

α = 10

111.56(D) 111.56(D) 111.36(D) 111.15(D) 110.53(D) 104.58(D) 88.76(D) 76.83(D) 69.22(D) 65.92(D) 66.75(D) 71.07(D) 83.62(F) 76.83(D) 77.86(D)

fbc10

89.57 88.33 73.09 57.39 47.29 40.96 38.15 38.98 42.01 46.63 50.28 53.93 56.41

89.59

Ref. [38]

62.73

54.14 54.54

64.85

88.84 87.52

102.04

Ref. [40] 78.45(D) 78.39(D) 78.21(D) 77.83(D) 76.92(D) 73.07(D) 64.10(D) 57.11(D) 54.03(D) 54.64(D) 57.58(D) 59.45(D) 60.21(D) 61.66(D) 63.61(D)

fbc11

α = 30

79.21(D) 79.09(D) 78.97(D) 78.49(D) 77.40(D) 71.86(D) 62.19(D) 55.40(D) 52.62(D) 53.71(D) 57.59(D) 59.65(D) 60.13(D) 61.70(D) 63.76(D)

fbc10

Table 4: Critical fluid velocities (m/s) at the middle of Simply Supported - Simply Supported conical shells with different vertex angles and their corresponding unstable longitudinal half-wave (m) as well as instability type (F = Flutter or D = Divergence).

39

129.17

128.17

96.21

87.65

86.06

124.03

96.65

91.26

92.92

130.37

147.26

145.60

130.93

148.50

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15

Ref. [40]

Ref. [38]

n 120.69(F) 120.60(F) 120.42(F) 120.18(F) 120.96(D) 120.48(D) 118.11(D) 114.93(D) 104.16(D) 96.06(D) 90.72(D) 88.02(D) 87.48(D) 88.05(D) 87.15(D)

fbc11

α=0

121.00(F) 120.75(F) 125.25(D) 120.50(F) 121.75(F) 135.75(F) 120.00(F) 113.00(D) 103.00(D) 95.25(D) 90.25(D) 88.00(D) 114.00(F) 110.75(F) 95.00(F)

fbc10 138.96 137.80 137.30 136.57 134.40 129.78 116.15 104.94 96.27 89.59 86.51 85.45 87.15 87.93 88.77

Ref. [38]

84.26

85.16

92.68

122.84

143.17

137.66

139.66

Ref. [40] 119.62(F) 119.50(F) 121.22(F) 125.01(F) 125.26(F) 124.89(F) 120.68(D) 109.92(F) 99.96(D) 92.54(D) 87.75(D) 85.56(D) 86.00(D) 108.64(F) 111.08(F)

fbc11

α = 10

115.67(F) 131.66(F) 112.38(D) 120.00(F) 120.59(F) 119.77(F) 120.39(D) 108.69(D) 99.24(D) 92.05(D) 87.52(D) 85.67(D) 95.13(F) 116.49(D) 118.13(F)

fbc10 96.65 96.07 96.65 93.77 91.25 79.72 68.85 61.42 56.39 53.92 53.18 54.75 57.91 60.76 63.45

Ref. [38]

88.60

80.44

68.34

90.26

114.83

115.41

107.90

Ref. [40]

80.84(D) 92.97(F) 95.05(F) 97.79(F) 97.66(F) 92.02(D) 80.82(D) 71.13(D) 81.49(F) 83.64(F) 85.40(F) 84.62(F) 84.39(F) 84.70(F)

fbc11

α = 30

90.62(F) 91.46(D) 80.65(D) 80.41(D) 79.93(D) 79.57(D) 86.06(D) 75.59(D) 71.25(D) 85.58(F) 83.90(F) 85.22(F) 84.50(F) 84.26(F) 84.62(F)

fbc10

Table 5: Critical fluid velocities (m/s) at the middle of Clamped-Clamped conical shells with different vertex angles and their corresponding unstable longitudinal half-wave (m) as well as instability type (F = Flutter or D = Divergence).

Table 6: Critical fluid velocities (m/s) at the middle of Clamped-Free conical shells with different vertex angles and their corresponding unstable longitudinal half-wave (m) as well as instability type (F = Flutter or D = Divergence). α=0

α = 10

α = 30

n

Ref. [38]

fbc11

Ref. [38]

fbc11

Ref. [38]

fbc11

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15

148.50

90.03(F) 87.03(F) 93.03(F) 90.03(F) 82.83(D) 65.25(D) 54.81(D) 50.01(D) 49.56(D) 51.42(D) 53.34(D) 54.69(D) 56.16(D) 58.20(D) 60.78(D)

138.10

119.35(F) 119.30(F) 119.18(F) 107.60(D) 80.33(D) 63.39(D) 53.63(D) 49.58(D) 49.54(D) 52.47(D) 54.84(D) 56.45(D) 58.27(D) 60.79(D) 63.83(D)

97.06

80.82(D) 93.91(F) 96.24(F) 74.43(D) 56.30(D) 46.04(D) 42.08(D) 43.30(D) 47.46(D) 50.66(D) 52.54(D) 55.28(D) 59.25(D) 62.91(D) 69.54(D)

147.67 114.49 82.96 80.47 89.60

91.20

137.70 119.53 92.08 73.86 63.87 61.39 66.11 78.81 80.29 80.52 82.95 86.32 87.93

84.62 59.58 44.79 35.58 32.35 34.00 39.70 48.12 50.46 54.71 58.13 60.40 63.87

Table 7: Critical fluid velocities (m/s) at the middle of Clamped - Simply Supported conical shells with different vertex angles and their corresponding unstable longitudinal half-wave (m) as well as instability type (F = Flutter or D = Divergence). α=0

α = 10

α = 30

n

fbc11

fbc10

fbc11

fbc10

fbc11

fbc10

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15

120.60(F) 120.51(F) 120.33(F) 120.09(F) 120.75(D) 120.30(D) 113.58(D) 100.38(D) 90.21(D) 83.40(D) 79.98(D) 79.59(D) 81.33(D) 82.80(D) 83.55(D)

120.25(F) 119.25(D) 119.25(D) 119.00(D) 118.75(D) 118.50(D) 107.50(D) 94.25(D) 84.75(D) 78.50(D) 75.50(D) 75.25(D) 77.25(D) 79.00(D) 80.50(D)

119.60(F) 119.47(F) 125.43(F) 125.26(F) 124.77(F) 123.51(F) 111.20(D) 97.62(D) 88.05(D) 81.88(D) 79.14(D) 79.51(D) 81.81(D) 82.37(D) 82.74(D)

115.67(F) 127.97(F) 112.38(D) 120.80(F) 119.98(F) 111.36(D) 104.58(D) 92.05(D) 83.00(D) 77.45(D) 75.18(D) 75.80(D) 78.06(D) 79.50(D) 80.74(D)

80.87(D) 92.89(F) 94.89(F) 97.50(F) 87.08(D) 74.44(D) 70.32(D) 65.99(D) 65.36(D) 79.09(F) 85.24(F) 84.44(F) 84.20(F) 84.50(F)

79.21(F) 79.09(D) 78.97(D) 78.49(D) 77.40(D) 71.86(D) 62.19(D) 55.40(D) 52.62(D) 53.71(D) 87.59(D) 59.65(D) 60.13(F) 61.70(D) 63.73(D)

40

Table 8: Critical fluid velocities (m/s) at the middle of Simply Supported - Clamped conical shells with different vertex angles and their corresponding unstable longitudinal half-wave (m) as well as instability type (F = Flutter or D = Divergence). α=0

α = 10

α = 30

n

fbc11

fbc10

fbc11

fbc10

fbc11

fbc10

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15

120.60(F) 120.51(F) 120.33(F) 120.09(F) 120.75(D) 120.30(D) 113.58(D) 100.38(D) 90.21(D) 83.40(D) 79.98(D) 79.59(D) 81.33(D) 82.80(D) 83.55(D)

121.75(D) 120.75(F) 120.50(F) 120.25(F) 120.00(F) 120.5(D) 111.25(D) 98.75(D) 89.25(D) 83.00(D) 79.75(D) 79.75(D) 81.50(D) 82.75(D) 83.50(D)

111.25(D) 111.18(D) 111.00(D) 110.73(D) 110.36(D) 112.06(F) 105.93(D) 94.93(D) 85.58(D) 79.31(D) 76.28(D) 76.10(D) 77.71(D) 79.19(D) 80.28(D)

111.56(D) 111.56(D) 111.36(D) 111.15(D) 110.74(D) 109.92(D) 104.99(D) 93.90(D) 84.85(D) 79.09(D) 76.21(D) 76.21(D) 77.86(D) 79.09(D) 80.33(D)

78.47(D) 78.39(D) 78.22(D) 77.90(D) 77.28(D) 75.82(F) 71.51(D) 64.96(D) 60.00(D) 57.78(D) 57.97(D) 59.31(D) 60.51(D) 61.82(D) 63.56(D)

79.21(D) 79.09(D) 78.97(D) 78.61(D) 77.88(D) 76.20(D) 71.49(D) 64.85(D) 60.01(D) 57.83(D) 58.07(D) 59.40(D) 60.61(D) 61.95(D) 63.76(D)

Table 9: Critical fluid velocities (m/s) at the middle of Simply Supported - Free conical shells with different vertex angles and their corresponding unstable longitudinal half-wave (m) as well as instability type (F = Flutter or D = Divergence). α=0

α = 10

α = 30

n

fbc11

fbc10

fbc11

fbc10

fbc11

fbc10

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15

0.33(D) 1.44(D) 3.78(D) 7.29(D) 12.21(D) 18.63(D) 26.67(D) 36.00(D) 45.45(D) 51.33(D) 52.56(D) 53.16(D) 54.90(D) 57.60(D) 60.66(D)

119.40(D) 119.25(D) 118.95(D) 139.95(D) 119.85(D) 117.30(F) 118.05(F) 114.90(F) 113.70(F) 81.90(D) 120.45(F) -

20.92(D) 6.13(D) 4.75(D) 7.39(D) 12.07(D) 18.42(D) 26.39(D) 35.83(D) 45.79(D) 52.47(D) 53.75(D) 54.59(D) 56.84(D) 60.18(D) 63.73(D)

111.39(D) 111.15(D) 110.91(D) 110.54(D) 109.92(D) 108.81(D) 106.84(D) 106.10(D) 88.59(D) 80.33(D) 77.12(D) 77.36(D) 90.07(F) 80.32(D) 80.70(D)

35.64(D) 11.98(D) 7.06(D) 8.20(D) 12.65(D) 19.21(F) 27.68(D) 37.95(D) 48.51(D) 48.49(D) 49.5(D)1 53.14(D) 60.86(F) 61.20(D) 63.04(D)

79.09(D) 78.66(D) 78.08(D) 77.14(D) 74.82(D) 69.97(D) 64.68(D) 60.83(D) 58.51(D) 57.85(D) 58.51(D) 59.60(D) 60.61(D) 62.00(D) 63.74(D)

Table 10: Critical fluid velocities (m/s) at the middle of Free - Clamped conical shells with different vertex angles and their corresponding unstable longitudinal half-wave (m) as well as instability type (F = Flutter or D = Divergence). α=0

α = 10

α = 30

n

fbc11

fbc10

fbc11

fbc10

fbc11

fbc10

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15

90.03(F) 87.03(F) 93.03(F) 90.03(F) 82.83(D) 65.25(D) 54.81(D) 50.01(D) 49.56(D) 51.42(D) 53.34(D) 54.69(D) 56.16(D) 58.20(D) 60.78(D)

108.00(D) 107.75(D) 107.75(D) 105.25(D) 82.00(D) 65.00(D) 54.75(D) 50.00(D) 49.75(D) 51.50(D) 53.50(D) 54.75(D) 56.25(D) 58.25(D) 61.00(D)

92.22(D) 92.17(D) 91.97(D) 90.96(D) 77.91(D) 61.81(D) 52.12(D) 47.75(D) 47.35(D) 48.78(D) 50.14(D) 51.13(D) 52.42(D) 54.28(D) 56.57(D)

93.48(D) 93.48(D) 93.28(D) 92.05(D) 77.44(D) 61.60(D) 52.13(D) 47.80(D) 47.39(D) 48.83(D) 50.28(D) 51.30(D) 52.54(D) 54.39(D) 56.66(D)

59.50(D) 59.47(D) 59.25(D) 58.09(D) 51.19(D) 42.12(D) 37.45(D) 36.46(D) 36.96(D) 37.38(D) 37.85(D) 38.82(D) 40.33(D) 42.14(D) 44.13(D)

60.74(D) 60.74(D) 60.49(D) 59.04(D) 51.40(D) 42.17(D) 37.55(D) 36.58(D) 37.06(D) 37.43(D) 37.91(D) 39.01(D) 40.47(D) 42.29(D) 44.36(D)

41

(a) S − S, fbc11

(b) S − S, fbc10

(c) C − C, fbc11

(d) C − C, fbc10

(e) C − F, fbc11 Figure 7: Critical velocity variations versus the semi-vertex angle for conical shells at different boundary conditions l cos α = 0.9144m, a = 0.876m , h = 1.5mm , Na = 40, Q9 elements , Mesh = 4961DOF.

42

Table 11: Critical fluid velocities (m/s) at the middle of Free - Simply Supported conical shells with different vertex angles and their corresponding unstable longitudinal half-wave (m) as well as instability type (F = Flutter or D = Divergence). α=0

α = 10

α = 30

n

fbc11

fbc10

fbc11

fbc10

fbc11

fbc10

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15

0.33(D) 1.44(D) 3.78(D) 7.29(D) 12.21(D) 18.63(D) 26.67(D) 36.00(D) 45.45(D) 51.33(D) 52.56(D) 53.16(D) 54.90(D) 57.60(D) 60.66(D)

0.25(D) 1.25(D) 3.50(D) 6.75(D) 11.50(D) 17.75(D) 25.75(D) 35.00(D) 45.00(D) 51.50(D) 52.50(D) 53.00(D) 54.75(D) 57.50(D) 60.75(D)

23.27(D) 6.50(D) 5.00(D) 7.82(D) 12.74(D) 19.34(D) 27.52(D) 36.75(D) 45.03(D) 48.81(D) 49.30(D) 49.87(D) 51.50(D) 53.90(D) 56.49(D)

20.20(D) 5.77(D) 4.53(D) 7.21(D) 11.95(D) 18.55(D) 26.59(D) 36.06(D) 44.92(D) 48.83(D) 49.25(D) 49.66(D) 51.51(D) 53.98(D) 56.66(D)

52.95(D) 14.82(D) 8.48(D) 9.83(D) 15.02(D) 22.32(D) 30.35(D) 35.83(D) 36.90(D) 36.68(D) 37.16(D) 38.49(D) 40.24(D) 42.14(D) 44.13(D)

44.84(D) 13.54(D) 7.93(D) 9.27(D) 14.40(D) 21.59(D) 29.88(D) 35.85(D) 36.94(D) 36.70(D) 37.18(D) 38.64(D) 40.35(D) 42.29(D) 44.36(D)

43

44

221.766 221.760 221.742 221.712 221.670 221.616 221.549 221.470 221.379 221.276 221.160 221.032 220.891

C

816.793 816.791 816.783 816.771 816.754 816.732 816.705 816.673 816.636 816.594 816.547 816.494 816.437

IC

n=1

268.3 268.3 268.3 268.4 268.5 268.5 268.6 268.8 268.9 269.0 269.2 269.4 269.6

% 269.248 269.242 269.223 269.191 269.145 269.087 269.016 268.932 268.835 268.725 268.601 268.465 268.314

C 621.423 621.421 621.415 621.406 621.393 621.376 621.355 621.331 621.302 621.270 621.234 621.195 621.151

IC

n=2

130.8 130.8 130.8 130.8 130.9 130.9 131.0 131.0 131.1 131.2 131.3 131.4 131.5

% 323.375 323.368 323.347 323.312 323.264 323.201 323.125 323.034 322.930 322.811 322.678 322.530 322.367

C 477.027 477.025 477.019 477.008 476.993 476.973 476.949 476.921 476.889 476.852 476.811 476.765 476.715

IC

n=3

47.5 47.5 47.5 47.5 47.6 47.6 47.6 47.6 47.7 47.7 47.8 47.8 47.9

% 361.764 361.757 361.735 361.698 361.646 361.579 361.497 361.399 361.287 361.158 361.014 360.853 360.675

C 376.068 376.066 376.058 376.045 376.027 376.004 375.976 375.942 375.904 375.860 375.811 375.757 375.698

IC

n=4

4.0 4.0 4.0 4.0 4.0 4.0 4.0 4.0 4.0 4.1 4.1 4.1 4.2

% 302.786 302.782 302.770 302.748 302.719 302.681 302.634 302.579 302.515 302.441 302.359 302.267 302.166

C 303.972 303.969 303.960 303.945 303.924 303.897 303.864 303.825 303.780 303.729 303.672 303.609 303.540

IC

n=5

0.4 0.4 0.4 0.4 0.4 0.4 0.4 0.4 0.4 0.4 0.4 0.4 0.5

%

0 25 50 75 100 125 150 175 200 225 250 275 300

Input U ((m/s)

216.937 216.936 216.931 216.923 216.911 216.898 216.881 216.863 216.843 216.822 216.801 216.781 216.763

C

653.119 653.118 653.116 653.113 653.108 653.102 653.095 653.087 653.077 653.066 653.053 653.039 653.024

IC

n=1

201.1 201.1 201.1 201.1 201.1 201.1 201.1 201.2 201.2 201.2 201.2 201.2 201.3

% 266.404 266.402 266.398 266.391 266.381 266.369 266.355 266.339 266.323 266.305 266.289 266.273 266.259

C 495.335 495.334 495.333 495.331 495.327 495.323 495.318 495.312 495.305 495.297 495.288 495.278 495.267

IC

n=2

85.9 85.9 85.9 85.9 85.9 86.0 86.0 86.0 86.0 86.0 86.0 86.0 86.0

% 317.896 317.895 317.890 317.883 317.873 317.862 317.848 317.834 317.819 317.804 317.790 317.777 317.767

C 371.734 371.734 371.732 371.729 371.724 371.719 371.712 371.704 371.695 371.684 371.672 371.659 371.645

IC

n=3

16.9 16.9 16.9 16.9 16.9 16.9 16.9 16.9 17.0 17.0 17.0 17.0 17.0

%

285.113 285.112 285.108 285.103 285.096 285.087 285.078 285.068 285.057 285.048 285.039 285.032 285.026

C

287.094 287.094 287.091 287.087 287.082 287.075 287.066 287.056 287.044 287.031 287.016 286.999 286.981

IC

n=4

0.7 0.7 0.7 0.7 0.7 0.7 0.7 0.7 0.7 0.7 0.7 0.7 0.7

%

227.551 227.550 227.547 227.541 227.534 227.526 227.516 227.506 227.496 227.486 227.478 227.470 227.465

C

227.939 227.938 227.935 227.930 227.923 227.915 227.904 227.892 227.877 227.861 227.843 227.823 227.801

IC

n=5

0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.1

%

Table 13: Compressibility effects on the natural frequencies (Hz) of conical shells (α = 30) at different input velocities with C-C boundary condition, l = 0.9144m, a = 0.876m , h = 1.5mm , M = 40, Q9 elements , Mesh = 4961DOF, fbc11, (C = Compressible model ; IC = Incompressible model).

0 25 50 75 100 125 150 175 200 225 250 275 300

Input U (m/s)

Table 12: Compressibility effects on the natural frequencies (Hz) of cylindrical shells at different input velocities with C-C boundary condition, l = 0.9144m, a = 0.876m , h = 1.5mm , Na = 40, Q9 elements , Mesh = 4961DOF, fbc11, (C = Compressible model ; IC = Incompressible model).

45

221.637 221.631 221.613 221.581 221.538 221.482 221.413 221.332 221.239 221.132 221.013 220.881 220.736

C 520.502 520.501 520.499 520.495 520.490 520.483 520.474 520.464 520.453 520.439 520.425 520.408 520.391

IC

n=1

134.8 134.9 134.9 134.9 134.9 135.0 135.1 135.2 135.2 135.4 135.5 135.6 135.8

% 266.538 266.530 266.504 266.461 266.401 266.323 266.229 266.117 265.987 265.840 265.676 265.494 265.295

C 315.309 315.308 315.304 315.296 315.287 315.274 315.258 315.240 315.218 315.194 315.167 315.137 315.105

IC

n=2

18.3 18.3 18.3 18.3 18.4 18.4 18.4 18.5 18.5 18.6 18.6 18.7 18.8

% 201.775 201.771 201.761 201.744 201.720 201.689 201.651 201.606 201.555 201.496 201.430 201.356 201.276

C 202.303 202.301 202.294 202.283 202.267 202.246 202.221 202.191 202.157 202.118 202.074 202.026 201.973

IC

n=3

0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3

% 137.076 137.072 137.061 137.042 137.016 136.982 136.941 136.892 136.835 136.771 136.699 136.619 136.531

C 137.147 137.144 137.133 137.117 137.093 137.063 137.026 136.983 136.933 136.876 136.812 136.742 136.664

IC

n=4

0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1

% 97.972 97.967 97.953 97.930 97.897 97.855 97.803 97.742 97.671 97.590 97.500 97.400 97.290

C 97.987 97.983 97.969 97.947 97.915 97.875 97.826 97.767 97.700 97.623 97.537 97.442 97.338

IC

n=5

0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0

%

0 25 50 75 100 125 150 175 200 225 250 275 300

Input U ((m/s)

215.826 215.825 215.819 215.810 215.797 215.781 215.763 215.743 215.721 215.699 215.677 215.656 215.637

C 275.641 275.641 275.640 275.639 275.637 275.635 275.632 275.628 275.624 275.619 275.614 275.609 275.602

IC

n=1

27.7 27.7 27.7 27.7 27.7 27.7 27.7 27.8 27.8 27.8 27.8 27.8 27.8

% 152.847 152.846 152.844 152.841 152.837 152.832 152.826 152.821 152.815 152.810 152.805 152.801 152.799

C 153.158 153.158 153.156 153.154 153.150 153.146 153.140 153.134 153.126 153.118 153.109 153.098 153.087

IC

n=2

0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2

% 93.578 93.578 93.575 93.572 93.567 93.561 93.555 93.548 93.542 93.536 93.531 93.527 93.524

C 93.611 93.610 93.608 93.604 93.599 93.593 93.585 93.575 93.565 93.552 93.538 93.523 93.506

IC

n=3

0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0

% 61.841 61.840 61.837 61.832 61.826 61.819 61.811 61.802 61.794 61.787 61.781 61.775 61.772

C

61.846 61.845 61.842 61.837 61.831 61.822 61.811 61.799 61.784 61.768 61.749 61.729 61.706

IC

n=4

0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 -0.1 -0.1 -0.1

%

43.922 43.921 43.917 43.911 43.904 43.895 43.885 43.875 43.865 43.856 43.848 43.842 43.838

C

43.923 43.922 43.918 43.912 43.904 43.893 43.880 43.864 43.846 43.825 43.802 43.777 43.749

IC

n=5

0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.1 0.1 0.1 0.2

%

Table 15: Compressibility effects on the natural frequencies (Hz) of conical shells (α = 30) at different input velocities with C-F boundary condition, l = 0.9144m, a = 0.876m , h = 1.5mm , M = 40, Q9 elements , Mesh = 4961DOF, fbc11, (C = Compressible model ; IC = Incompressible model).

0 25 50 75 100 125 150 175 200 225 250 275 300

Input U (m/s)

Table 14: Compressibility effects on the natural frequencies (Hz) of cylindrical shells at different input velocities with C-F boundary condition, l = 0.9144m, a = 0.876m , h = 1.5mm , M = 40, Q9 elements , Mesh = 4961DOF, fbc11, (C = Compressible model ; IC = Incompressible model).

46

221.762 221.756 221.738 221.708 221.665 221.611 221.544 221.465 221.373 221.269 221.153 221.024 220.882

C

582.797 582.797 582.797 582.797 582.797 582.797 582.797 582.797 582.797 582.797 582.797 582.797 582.797

IC

n=1

162.8 162.8 162.8 162.9 162.9 163.0 163.1 163.2 163.3 163.4 163.5 163.7 163.8

% 269.242 269.236 269.216 269.184 269.139 269.080 269.009 268.924 268.827 268.715 268.591 268.453 268.302

C 616.299 616.297 616.291 616.282 616.269 616.252 616.232 616.208 616.180 616.148 616.113 616.074 616.031

IC

n=2

128.9 128.9 128.9 128.9 129.0 129.0 129.1 129.1 129.2 129.3 129.4 129.5 129.6

% 322.965 322.958 322.937 322.902 322.852 322.789 322.711 322.618 322.512 322.391 322.255 322.104 321.939

C 448.951 448.949 448.942 448.931 448.915 448.895 448.870 448.841 448.807 448.769 448.726 448.679 448.627

IC

n=3

39.0 39.0 39.0 39.0 39.0 39.1 39.1 39.1 39.2 39.2 39.2 39.3 39.4

% 323.309 323.304 323.288 323.262 323.225 323.177 323.118 323.049 322.968 322.876 322.772 322.657 322.529

C 327.389 327.386 327.378 327.363 327.343 327.318 327.286 327.249 327.206 327.158 327.104 327.044 326.978

IC

n=4

1.3 1.3 1.3 1.3 1.3 1.3 1.3 1.3 1.3 1.3 1.3 1.4 1.4

% 243.014 243.010 242.997 242.975 242.945 242.906 242.858 242.801 242.735 242.660 242.576 242.482 242.378

C 243.488 243.484 243.473 243.456 243.431 243.398 243.359 243.313 243.259 243.199 243.131 243.056 242.973

IC

n=5

0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2 0.2

%

0 25 50 75 100 125 150 175 200 225 250 275 300

Input U ((m/s)

216.847 216.845 216.840 216.832 216.821 216.807 216.790 216.771 216.751 216.730 216.709 216.688 216.670

C 472.508 472.508 472.507 472.507 472.505 472.504 472.502 472.499 472.497 472.494 472.490 472.487 472.483

IC

n=1

117.9 117.9 117.9 117.9 117.9 117.9 118.0 118.0 118.0 118.0 118.0 118.0 118.1

% 266.365 266.364 266.359 266.352 266.342 266.330 266.316 266.300 266.283 266.266 266.248 266.232 266.218

C 487.842 487.841 487.840 487.838 487.835 487.830 487.825 487.819 487.812 487.804 487.796 487.786 487.775

IC

n=2

83.1 83.1 83.2 83.2 83.2 83.2 83.2 83.2 83.2 83.2 83.2 83.2 83.2

% 312.277 312.275 312.270 312.263 312.253 312.240 312.226 312.211 312.195 312.179 312.165 312.152 312.142

C 338.713 338.712 338.710 338.706 338.701 338.695 338.687 338.678 338.667 338.655 338.642 338.627 338.611

IC

n=3

8.5 8.5 8.5 8.5 8.5 8.5 8.5 8.5 8.5 8.5 8.5 8.5 8.5

% 234.856 234.855 234.851 234.846 234.838 234.829 234.819 234.808 234.797 234.787 234.778 234.770 234.764

C

235.646 235.645 235.642 235.637 235.631 235.622 235.611 235.599 235.584 235.568 235.549 235.529 235.507

IC

n=4

0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3 0.3

%

169.490 169.489 169.485 169.478 169.470 169.459 169.448 169.436 169.423 169.411 169.401 169.392 169.385

C

169.628 169.626 169.623 169.616 169.608 169.596 169.583 169.567 169.548 169.527 169.503 169.477 169.448

IC

n=5

0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.0

%

Table 17: Compressibility effects on the natural frequencies (Hz) of conical shells (α = 30) at different input velocities with S-S boundary condition, l = 0.9144m, a = 0.876m , h = 1.5mm , M = 40, Q9 elements , Mesh = 4961DOF, fbc11, (C = Compressible model ; IC = Incompressible model).

0 25 50 75 100 125 150 175 200 225 250 275 300

Input U (m/s)

Table 16: Compressibility effects on the natural frequencies (Hz) of cylindrical shells at different input velocities with S-S boundary condition, l = 0.9144m, a = 0.876m , h = 1.5mm , M = 40, Q9 elements , Mesh = 4961DOF, fbc11 , (C = Compressible model ; IC = Incompressible model).