Analysis of Free Vibration Characteristics of Cylindrical Shells with Finite Submerged Depth Based on Energy Variational Principle

Abstract

Based on the principle of energy variation, a calculation model for the free vibration characteristics of a cylindrical shell with a finite submerged depth considering the influence of the free liquid surface is established in this paper. First, the Euler beam function is used instead of the shell axial displacement function to obtain the shell kinetic energy and potential energy. Then, by using the mirror image method, the analytical expression of the fluid velocity potential considering the free surface is obtained, and the flow field is added to the system energy functional in the form of fluid work. Then the energy functional is changed to obtain the shell–liquid coupled vibration equation. Solving the equation can obtain the natural frequencies and modes of the structure. The comparison with the finite element calculation results verifies the accuracy of the calculation model in this paper. The research on the influence of the free liquid surface shows that compared to the infinite domain, the free liquid surface destroys the symmetry of the entire system, resulting in a difference in the natural frequency of the positive and negative modes of the shell, and the circumferential mode shapes are no longer mutually uncoupled trigonometric functions. The existence of free liquid surface will also increase the natural frequency of the same order mode, and the closer to the free surface, the natural frequency is greater. As the immersion depth increases, the free vibration characteristics will quickly tend to the result of infinite domain. Additionally, when the immersion depth is equal to or greater than four times the radius of the shell structure, it can be considered that the free liquid surface has no effect. These law and phenomena have also been explained from the mechanism. The method in this paper provides a new analytical solution pattern for solving this type of problem.

1. Introduction

Cylindrical shell structures coupled with heavy fluids are widely used in many engineering fields such as underwater vehicles, pipelines, and submarine cables. Due to the existence of heavy fluid, the vibration and acoustic characteristics of the cylindrical shell structure coupled with it will change significantly. As a classic fluid–solid coupling problem, how to accurately predict the acoustic and vibration characteristics of a shell structure under heavy fluid loading has received extensive attention from predecessors [1,2,3,4,5,6,7]. Although when dealing with heavy fluids, some researchers regard it as a compressible acoustic medium [8], and some scholars regard it as an incompressible fluid medium [9], the heavy fluid can be regarded as an additional mass load on the structure when calculating the low-frequency free vibration characteristics of a cylindrical shell structure. The influence of fluid compressibility is very weak, and the difference between the fluid medium and the acoustic medium is negligible.

If the heavy fluid and the elastic shell are regarded as completely two-way coupled, the calculation of the shell mode under heavy fluid loading (wet mode) and the calculation of the shell mode in vacuum (dry mode) are two completely independent problems [10]. Therefore, it is necessary to find a solution to this type of fluid–solid problem. With the development of numerical calculation methods, numerical methods that combine finite element and boundary elements have also been used to solve the acoustic and vibration characteristics of elastic shells immersed in infinite domains [11,12,13,14]. However, the mixed numerical method also has the disadvantage of high computational cost, so there are still many scholars who are committed to using (semi-)analytic methods to solve this problem. Due to the good axisymmetric properties of the rotating shell immersed in the infinite domain, its circumferential displacement function can be expanded into a Fourier series that is not coupled with each other. Therefore, the integration of the sound (flow) field can be simplified to a one-dimensional integration problem along the shell generating line, which greatly reduces the computational amount [15,16,17]. Based on the wave propagation method (WPA), Li et al. [18] transferred the solution domain to the wavenumber domain and successfully solved the vibration characteristics of the elliptical shell. Caresta et al. [19] approximated the surface of the cone shell as a number of cylindrical shells and used a power series to describe the displacement of the cone shell. The surface sound pressure of the infinite cylindrical shell was used as a fluid load to calculate the low-frequency vibration characteristics of the cone shell. Subsequently, Caresta et al. [20] further extended this method to the calculation of low-frequency vibration characteristics of composite shell structures. Based on the energy principle and potential flow theory, Kwak [9] added the mass of water attached to the fluid in the form of a matrix to the free vibration equation of the cylindrical shell and calculated the natural vibration characteristics of the submerged cylindrical shell. Qu et al. [21] established the structural energy equation of the composite shell based on the modified variational principle and calculated the Helmholtz integral on the surface of the structure by using the collocation method, and the vibro-acoustic characteristics of underwater structures are well predicted. Similarly, Jin et al. [22] used an improved Fourier series to describe the shell displacement, combined with an artificial spring system, and integrated the sound pressure term along the shell surface generation line to calculate the acoustic and vibration characteristics of the underwater composite shell. Wang et al. [23] combined the wave superposition method (WSM) and the fine transfer matrix method (PTMM) to establish a semi-analytical calculation model for the vibration characteristics of an underwater composite shell. Xie et al. [24] regarded the cyclotron shell as a combination of several conical shell elements and assembled them as a whole through structural continuity conditions. The sound pressure terms were treated as a three-node element uniformly and met the condition that is velocity continuity at the joint surface of the shell. Additionally, the sound–structure coupling governing equation of the rotating shell structure is finally obtained.

Although the above publications have presented considerable research results, when the shell structure is close to the free liquid surface or the seabed, it is obviously not appropriate to assume that the flow field is infinite. A common method to deal with the free surface of a completely submerged structure is the mirror image method. The mirror image method was first used to study the radiation characteristics of simple sound sources near the plane boundary and showing Lloyd’s mirror behavior for near-surface monopoles [25,26,27]. The common finite element–boundary element coupling calculation method can also be combined with the mirror image method to modify the basic solution of the three-dimensional acoustic Helmholtz equation, and then realized the fluid–structure coupling vibration performance calculation considering the influence of the free liquid surface [28,29]. Ergin et al. [30] analyzed the vibration characteristics of cylindrical shells under finite immersion depth based on experiments and three-dimensional hydroelastic software and found that the closer the structure is to the free liquid surface, the greater the natural frequency of the same order. The above studies are of considerable value, but the numerical method a will require a greater computational cost with considering the free liquid surface. Therefore, Li et al. [31,32,33,34] have done a lot of work in this field. Their research results show that when the immersion depth of the structure is greater than or equal to four times the radius of the structure, the vibro-acoustic characteristics of the cylindrical shell tend to be in an infinite domain.

Different from the previous method used by our research team [31,32,33,34], this article will obtain the shell–liquid coupled vibration equation based on the energy variation principle, which is a new pattern of solution. First, based on the assumption of no rotation of the shell and the Love shell theory, the Euler beam function is used to replace the axial displacement function of the cylindrical shell, then the kinetic energy and potential energy of the shell are obtained. Subsequently, based on the theory of potential flow, combined with the mirror image method, the potential distribution of the fluid velocity considering the free surface is obtained. Then, according to the continuous condition of the velocity of the shell–liquid junction, the migration matrix of the velocity vector of the flow field and the displacement vector of the shell is obtained, and the flow field is added to the entire system in the form of work. Finally, the energy functional is changed to obtain the shell–liquid coupled vibration equation, solve it, and obtain the free vibration characteristics of the shell. By comparing with the finite element calculation results, the correctness of the method in this paper is verified. Through an in-depth study of the problem, considering the influence of the free liquid surface, the natural frequencies change law and the mode shape coupling phenomenon of the coupled system were discovered and the mechanism was explained.

2. Theoretical Analysis

The cylindrical shell with finite submerged depth considering the influence of the free liquid surface is shown in Figure 1. The length of the cylindrical shell is L, the radius of the middle surface is R, and the thickness of the shell is h. The material parameters of the cylindrical shell are elastic modulus E, Poisson’s ratio μ, density ρ. The density of the fluid is ρ_f. The immersion depth H is defined as the distance between the free liquid surface and the horizontal half-transverse surface of the shell. Take the center of the left end surface of the cylindrical shell as the coordinate origin O, and take the cylindrical coordinate system (r, θ, x) shown in Figure 1 as the coordinate system. Among them, r, θ, x represent radial direction, circumferential direction, and axial direction, respectively. Additionally, u, v, w are, respectively, the displacements of the middle surface of the shell in all directions. At the same time, it is assumed that there are virtual sealing plates at both ends of the cylindrical shell, that is, the internal space of the cylindrical shell is vacuum and fluid-free.

2.1. Shell Kinetic Energy and Potential Energy

The research object of this paper is an isotropic thin-walled cylindrical shell and the midplane displacement of the shell can be assumed as the following form [35]:

$$\{\begin{array}{l}u={\displaystyle \sum _{m=1}^{+\infty }{\displaystyle \sum _{n=-\infty }^{+\infty }{U}_{mn}{X^{\prime }}_m(\eta x)\exp (in\phi )\cos (\omega t)}}\\ v={\displaystyle \sum _{m=1}^{+\infty }{\displaystyle \sum _{n=-\infty }^{+\infty }{V}_{mn}{X}_m(\eta x)\exp (in\phi )\cos (\omega t)}}\\ w={\displaystyle \sum _{m=1}^{+\infty }{\displaystyle \sum _{n=-\infty }^{+\infty }{W}_{mn}{X}_m(\eta x)\exp (in\phi )\cos (\omega t)}}\end{array}$$

(1)

In Equation (1), ω is the circular frequency; m is the order of the axial mode; n is the order of the circumferential mode; U_mn, V_mn, and W_mn are three-way displacement amplitude coefficients; X_m(ηx) is the deflection function of Euler beam. “’” represents the derivative of the deflection function with respect to x. For any boundary condition, the deflection function of the beam can be taken as:

$$X_m\left(\eta x\right)=C_{1}sin(\eta x)+C_{2}cos(\eta x)+C_{1}sinh(\eta x)+C_{1}cosh(\eta x)$$

(2)

In this equation, C₁, C₂, C₃, and C₄ are undetermined coefficients determined by boundary conditions. η is usually called the axial half-wave number k_m and there are corresponding expressions under different boundary conditions. Since trigonometric functions and hyperbolic functions are different in the subsequent treatment of fluid–structure coupling conditions, it is assumed that:

$$\{\begin{array}{l}{X}_m(k_{m}x)=X_{m1}(k_{m}x)+X_{m2}(k_{m}x)\\ X_{m1}(k_{m}x)=C_1\sin (k_{m}x)+C_2\cos (k_{m}x)\\ X_{m2}(k_{m}x)=C_3\sinh (k_{m}x)+C_4\cosh (k_{m}x)\end{array}$$

(3)

As mentioned earlier, this article will study the free vibration characteristics of a finite-length cylindrical shell at a finite submerged depth based on the energy functional variational principle. Therefore, the specific expression form of each part of the energy should be clarified.

The kinetic energy expression of the thin shell ignoring the rotation effect is as follows [36]:

$$T=\frac{\rho h{\omega }^2}{2}{\displaystyle {\iiint }_V({\dot{u}}^2+{\dot{v}}^2+{\dot{w}}^2)\mathrm{d}V}$$

(4)

In Equation (4), “˙” represents the derivative of the displacement component with respect to time. V and dV represent the volume division and volume microelement of the cylindrical shell, respectively.

Substituting Equations (1)–(3) into Equation (4), and due to the orthogonality of the displacement function, the integral shell kinetic energy can be expressed as:

$$T=\frac{1}{2}{\displaystyle \sum _{m=1}^{+\infty }{\displaystyle \sum _{n=-\infty }^{+\infty }\left\{{\xi }_{mn}\right\}[M_{mn}]}}{\left\{{\xi }_{mn}\right\}}^{\mathrm{T}}$$

(5)

In Equation (5), {ξ_mn} is composed of unknown displacement amplitude coefficient, {ξ_mn} = {U_mn, V_mn, W_mn}. The mass matrix [M_mn] is a third-order diagonal matrix.

The expression of the strain energy of thin shells [36] (this article is based on Love shell theory) is as follows:

$$U=\frac{1}{2}{\displaystyle {\iiint }_V{\mathsf{\epsilon }}^{\mathrm{T}}\mathsf{\sigma }\mathrm{d}V}$$

(6)

In Equation (6), ε is the strain vector; σ is the stress vector. The specific forms of strain and stress vector are not listed here, and they can be found in various related documents.

Similarly, the integrated shell strain energy can be expressed as:

$$U=\frac{1}{2}{\displaystyle \sum _{m=1}^{+\infty }{\displaystyle \sum _{n=-\infty }^{+\infty }\left\{{\xi }_{mn}\right\}[K_{mn}]{\left\{{\xi }_{mn}\right\}}^{\mathrm{T}}}}$$

(7)

In this equation, the stiffness matrix [K_mn] is the third-order Hermite matrix.

2.2. Fluid–Structure Coupling Condition

In this paper, the fluid–solid coupling effect will be added to the energy functional of the system in the form of fluid work. First, it is necessary to determine the expression of the velocity potential function at a finite immersion depth (considering the influence of the free surface).

Based on the theory of potential flow, and assuming that the fluid is an incompressible, non-rotational, and inviscid ideal fluid, the velocity potential function ϕ(r,x,θ,t) in the cylindrical coordinate system needs to satisfy the Laplace equation in the cylindrical coordinate system:

$$\frac{1}{r}\frac{\partial }{\partial r}\left(r\frac{\partial \varphi }{\partial r}\right)+\frac{1}{r^2}\frac{{\partial }^2\varphi }{\partial {\theta }^2}+\frac{{\partial }^2\varphi }{\partial x^2}=0$$

(8)

The condition of the velocity potential function at infinity is that the velocity potential is zero.

$${\frac{\partial \varphi }{\partial r}|}_{r=\infty }=0$$

(9)

The free liquid surface can be processed by the principle of image, that is, the velocity potential can be divided into the real source velocity potential directly generated by structural vibration and the virtual source velocity potential generated by the reflection of the free surface. As shown in Figure 2, the virtual source coordinate system (r’, x’, θ’) and the real source coordinate system (r, x, θ) are symmetric about the free liquid surface.

Assuming a random point P in the flow field, its velocity potential can be expressed as:

$${\varphi }_P(r,x,\theta ,t)={\varphi }_P{}^r(r,x,\theta ,t)+{\varphi }_P{}^i(r^{\prime },x^{\prime },{\theta }^{\prime },t)$$

(10)

In Equation (10), ϕ_P^r(r,x,θ,t) is the real source fluid velocity potential; ϕ_Pⁱ(r’,x’,θ’,t) is the virtual source fluid velocity potential.

The velocity potential function that satisfies Equations (8) and (9) is expressed as follows:

$$\{\begin{array}{l}{\varphi }_P{}^r(r,x,\theta ,t)=\\ {\displaystyle \sum _{m=1}^{\infty }{\displaystyle \sum _{n=-\infty }^{+\infty }\left(\begin{array}{l}{\varphi }_{m1,n}{}^r{K}_n(k_{m}r)X_{m1}(k_{m}x)\\ +{\varphi }_{m2,n}{}^r{H}_n{}^{(2)}(k_{m}r)X_{m2}(k_{m}x)\end{array}\right)\exp (in\theta )\sin (\omega t)}}\\ {\varphi }_P{}^i(r^{\prime },x^{\prime },{\theta }^{\prime },t)=\\ {\displaystyle \sum _{m=1}^{\infty }{\displaystyle \sum _{n=-\infty }^{+\infty }\left(\begin{array}{l}{\varphi }_{m1,n}{}^i{K}_n(k_m{r}^{\prime })X_{m1}(k_m{x}^{\prime })\\ +{\varphi }_{m2,n}{}^r{H}_n{}^{(2)}(k_m{r}^{\prime })X_{m2}(k_m{x}^{\prime })\end{array}\right)\exp (in{\theta }^{\prime })\sin (\omega t)}}\end{array}$$

(11)

In Equation (11), K_n() is a modified Bessel function of the second kind of order n; H_n⁽²⁾() is a Hankel function of the second kind of order n.

When point P is above the free liquid surface, it has the following positional relationship under the coordinates of the two sets of coordinate systems:

$$r=r^{\prime },x=x^{\prime },\theta +{\theta }^{\prime }=\pi$$

(12)

If ignores gravity waves and other factors are not considered, the conditions for zero velocity potential at the free liquid surface are as follows:

$${\varphi }_P{}^r(r,x,\theta ,t)+{\varphi }_P{}^i(r^{\prime },x^{\prime },{\theta }^{\prime },t)=0$$

(13)

By substituting Equations (11) and (12) into Equation (13) and orthogonalizing it, the following formulas are obtained:

$$\{\begin{array}{l}{\varphi }_{m1,n}{}^r+{(-1)}^n{\varphi }_{m1,n}{}^i=0\\ {\varphi }_{m2,n}{}^r+{(-1)}^n{\varphi }_{m2,n}{}^i=0\end{array}$$

(14)

The analytical expression of the velocity potential function at any point in the flow field is:

$$\begin{array}{l}{\varphi }_P(r,x,\theta ,t)=\\ {\displaystyle \sum _{m=1}^{\infty }{\displaystyle \sum _{n=-\infty }^{+\infty }\left(\begin{array}{l}{\varphi }_{m1,n}{}^r\left(\begin{array}{l}{K}_n(k_{m}r)\exp (in\theta )\\ -{(-1)}^n{K}_{-n}(k_m{r}^{\prime })\exp (-in{\theta }^{\prime })\end{array}\right)X_{m1}(k_{m}x)\\ +{\varphi }_{m2,n}{}^r\left(\begin{array}{l}{H}_n{}^{(2)}(k_{m}r)\exp (in\theta )\\ -{(-1)}^n{H}_{-n}{}^{(2)}(k_m{r}^{\prime })\exp (-in{\theta }^{\prime })\end{array}\right)X_{m2}(k_{m}x)\end{array}\right)}}\sin (\omega t)\end{array}$$

(15)

According to the Graf plus theorem:

$$\{\begin{array}{l}{K}_{-n}(k_m{r}^{\prime })\exp (-in{\theta }^{\prime })=\\ \begin{array}{l}{\displaystyle \sum _{a=-\infty }^{+\infty }{(-1)}^a{K}_{a+n}(2k_{m}H)I_a(k_{m}r)\exp (ia\theta ) r<2H}\\ {\displaystyle \sum _{a=-\infty }^{+\infty }{(-1)}^a{I}_{a+n}(2k_{m}H)K_a(k_{m}r)\exp (ia\theta ) r\ge 2H}\end{array}\}\\ H_{-n}{}^{(2)}(k_m{r}^{\prime })\exp (-in{\theta }^{\prime })=\\ \begin{array}{l}{\displaystyle \sum _{a=-\infty }^{+\infty }{(-1)}^a{H}_{a+n}{}^{(2)}(2k_{m}H)J_a(k_{m}r)\exp (ia\theta ) r<2H}\\ {\displaystyle \sum _{a=-\infty }^{+\infty }{(-1)}^a{J}_{a+n}(2k_{m}H)H_a{}^{(2)}(k_{m}r)\exp (ia\theta ) r\ge 2H}\end{array}\}\end{array}$$

(16)

In the formula, I_a() is the first-kind modified Bessel function of order a; J_a() is the first-kind Bessel function of order a.

For a cylindrical shell structure with a finite submerged depth, the radius r_c ≈ R < 2H at the outer surface of the structure, so Equation (15) can be further expressed as the following form:

$$\begin{array}{l}{\varphi }_P(r,x,\theta ,t)=\\ {\displaystyle \sum _{m=1}^{\infty }{\displaystyle \sum _{n=-\infty }^{+\infty }\left(\begin{array}{l}{\varphi }_{m1,n}\left(\begin{array}{l}{K}_n(k_{m}r)\exp (in\theta )\\ -{\displaystyle \sum _{a=-\infty }^{+\infty }{(-1)}^{a+n}}{K}_{a+n}(2k_{m}H)I_a(k_{m}r)\exp (ia\theta )\end{array}\right)X_{m1}(k_{m}x)\\ +{\varphi }_{m2,n}\left(\begin{array}{l}{H}_n{}^{(2)}(k_{m}r)\exp (in\theta )\\ -{\displaystyle \sum _{a=-\infty }^{+\infty }{(-1)}^{a+n}}{H}_{a+n}{}^{(2)}(2k_{m}H)J_a(k_{m}r)\exp (ia\theta )\end{array}\right)X_{m2}(k_{m}x)\end{array}\right)}}\sin (\omega t)\end{array}$$

(17)

The coefficient a is equivalent to the coefficient n, so the summation order of the series can be exchanged, and Equation (17) can be further rewritten as:

$$\begin{array}{l}{\varphi }_P(r,x,\theta ,t)=\\ {\displaystyle \sum _{m=1}^{\infty }{\displaystyle \sum _{n=-\infty }^{+\infty }\left(\begin{array}{l}\left(\begin{array}{l}{\varphi }_{m1,n}{K}_n(k_{m}r)\\ -{\displaystyle \sum _{a=-\infty }^{+\infty }{(-1)}^{a+n}}{\varphi }_{m1,a}{K}_{a+n}(2k_{m}H)I_n(k_{m}r)\end{array}\right)X_{m1}(k_{m}x)\\ +\left(\begin{array}{l}{\varphi }_{m2,n}{H}_n{}^{(2)}(k_{m}r)\exp (in\theta )\\ -{\displaystyle \sum _{a=-\infty }^{+\infty }{(-1)}^{a+n}}{\varphi }_{m2,a}{H}_{a+n}{}^{(2)}(2k_{m}H)J_n(k_{m}r)\exp (ia\theta )\end{array}\right)X_{m2}(k_{m}x)\end{array}\right)}}\exp (ia\theta )\sin (\omega t)\end{array}$$

(18)

The junction of the cylindrical shell wall and the flow field needs to meet the velocity continuity condition as follows:

$$-{\frac{\partial {\varphi }_P}{\partial r}|}_{r=R}={\frac{\partial w}{\partial t}|}_{r=R}$$

(19)

Equation (18) is substituted into Equation (19) and it is assumed that the cutoff number of the circumferential wave number is N. After orthogonalization, the relationship of the velocity potential amplitude vector {ξ_MN} and the cylindrical shell radial displacement amplitude vector {W_MN} can be obtained as:

$$\left\{{\xi }_{MN}\right\}=[Q]\left\{W_{MN}\right\}$$

(20)

In Equation (20), [Q] is the migration matrix, {ξ_MN} = {ϕ_m1,_-N, ϕ_m1,_{-N + 1},…, ϕ_m2,_N}^T, {W_MN} = {W_m,_-N, W_m,_{-N + 1},…, W_m,_N}^T.

From the Bernoulli equation, the hydrodynamic pressure P_f at the wall of the cylindrical shell is:

$$P_f={\rho }_f{\frac{\partial {\varphi }_P}{\partial t}|}_{r=R}$$

(21)

Fluid work can be expressed as:

$$W_f=-\frac{1}{2}{\displaystyle {\int }_0^L{\displaystyle {\int }_{-\pi }^{\pi }{P}_{f}wR\mathrm{d}\theta \mathrm{d}}}x$$

(22)

So far, the energy functional expression of the entire fluid–structure coupling system is:

$$\Pi =U-W_f-T$$

(23)

Since the axial functions of the coupled system are orthogonal in the domain, it is only necessary to take different axial half-wave numbers M in sequence. According to the principle of energy variation, the partial derivative of the unknown amplitude coefficient of the displacement is obtained as:

$$\frac{\partial \Pi }{\partial q}=0$$

(24)

In Equation (24), q = [U_M,_-N,…,W_M,_N].

The matrix form of the vibration equations of the coupled system under any axial half-wave number is as follows:

$$\left([K]-{\omega }^2[M]\right)\left\{q\right\}=0$$

(25)

In Equation (25), [K] is the coupling system stiffness matrix, [M] is the coupling system mass matrix.

3. Validation of the Theoretical Method

The calculation model in this paper is the same as that in the literature [32], and the specific parameters are shown in

Table 1. Coupling system parameters.

	Length	L	1.284 m
Shell Geometry	Mean radius	R	0.18 m
	Thickness	h	0.003 m
Shell Material	Density	ρ	7850 kg/m3
	Poisson’s ratio	μ	0.3
	Elastic modulus	E	2.06 × 1011 Pa
Fluid	Density	ρf	1025 kg/m3

.

The beam function can be used to simulate the axial displacement function of a finite-length cylindrical shell under different boundary conditions [35]. This article will take a cylindrical shell with the shear diaphragm boundary condition (SD-SD) at both ends as an example for research. The axial displacement function of the shell under this boundary condition is as follows:

$$\{\begin{array}{l}{X}_m(\eta x)=\sin (k_{m}x)\\ k_m=m\pi /L m=1,2,3\dots \end{array}$$

(26)

3.1. Convergence of the Theoretical Method

When the cylindrical shell is in the infinite domain, the migration matrix [Q] is a diagonal matrix. That is to say, the circumferential modes of the coupled system are also in a decoupled state. At this time, the cut-off number N of the circumferential wave number only affects whether certain modes appear and does not affect the modal results that have already appeared. However, when considering the influence of the free surface, [Q] is no longer a diagonal matrix, that is, there may be a coupling relationship between the circumferential modes of the coupled system, and the cutoff number N will have an impact on the natural frequency calculation results. In order to discuss the convergence of the method in this paper, taking the working condition of latent depth H = 0.2 m as an example, the first 10 natural frequencies of the coupled system are calculated, and the results are shown in

Table 2. The first 10 natural frequencies of the coupled system under different cutoff numbers N when H = 0.2 m.

Order	N = 4	N = 6	N = 8	N = 10	N = 12	N = 14
1	103.41	103.40	103.40	103.40	103.40	103.40
2	103.47	103.47	103.47	103.47	103.47	103.47
3	117.55	117.29	117.26	117.26	117.26	117.26
4	117.68	117.41	117.39	117.38	117.38	117.38
5	210.36	209.69	209.60	209.59	209.59	209.59
6	210.43	209.76	209.67	209.65	209.65	209.65
7	223.95	223.77	223.76	223.76	223.76	223.76
8	224.04	223.86	223.85	223.85	223.85	223.85
9	251.61	250.63	250.51	250.50	250.50	250.50
10	251.62	250.64	250.53	250.51	250.51	250.51

.

It can be seen from Table 2 that the method in this paper has good convergence. When the cutoff number N reaches about 10, the first 10 natural frequency values of the coupled system have converged. Unless otherwise specified, the cutoff number is 10 for all calculation examples in this paper.

3.2. Accuracy of the Theoretical Method

In order to verify the accuracy of the method in this paper, the calculation results in this paper and the FEM calculation results under different diving depth conditions are compared, as shown in

Table 3. Comparison of the method in this paper and the natural frequency of FEM at different depths.

Mode	H = 0.21 m			H = 0.45 m			Infinite Domain
	Present	FEM	Δ	Present	FEM	Δ	Present	FEM	Δ
1	103.05	103.12	0.07	99.10	99.09	0.01	98.86	98.83	0.03
2	103.11	103.14	0.03	99.13	99.11	0.02	98.86	98.84	0.02
3	115.26	115.71	0.39	109.30	109.87	0.52	109.26	109.83	0.52
4	115.36	115.77	0.35	109.30	109.87	0.52	109.26	109.83	0.52
5	207.40	209.30	0.91	202.46	204.45	0.97	202.45	204.44	0.97
6	207.43	209.33	0.91	202.46	204.45	0.97	202.45	204.44	0.97
7	222.65	223.83	0.53	217.02	218.16	0.52	216.98	218.10	0.51
8	222.74	223.90	0.52	217.03	218.16	0.52	216.98	218.11	0.52
9	247.66	249.60	0.78	241.57	243.63	0.85	241.56	243.62	0.85
10	247.68	249.61	0.77	241.57	243.64	0.85	241.56	243.63	0.85

. The finite element software Nastran is used for FEM numerical simulation, and the fluid part of the coupling system will be simulated by the virtual mass module. When the virtual mass module is used for calculation, the effect of fluid on the structure is equivalent to additional water mass, so there is no need to build fluid element additionally, and the diving depth can be set according to the draught. The number of elements is 2400. This number of elements can meet the requirements for the number of nodes at any wavelength in the vibration calculation of the structure below 1000 Hz to be greater than or equal to 6, so as to ensure the calculation accuracy. The finite element mesh model is shown in Figure 3.

The error percentage of the calculation results of these two methods is defined as:

$$\Delta =\frac{\left|f_F-f_P\right|}{f_F}\times 100$$

(27)

In Equation (27), f_F and f_P represent the finite element and the natural frequency calculation results of the method in this paper, respectively.

It can be seen from Table 3 that under different immersion depths, the calculation results of the method in this paper and the result of FEM are in good agreement, and the maximum percentage error does not exceed 1, which shows that the method in this paper is accurate and reliable. In addition, as the immersion depth increases, the natural frequency of the same order of the coupled system gradually decreases, and eventually it will tend to the calculation result under infinite domain conditions, that is to say, the influence of the free liquid surface will disappear. In addition to the influence of the frequency value, the free liquid surface will also make the natural frequency of the positive and negative modes different, which is also different from a cylindrical shell under a vacuum or an infinite flow field.

4. Numerical Examples and Discussion

4.1. Influence of Free Liquid Surface on Free Vibration Characteristics

The foregoing qualitatively gives the influence of the free liquid surface on the natural frequency based on the content of Table 3. To further clarify the quantitative change law of natural frequency and immersion depth, the following takes the previous four-order natural frequencies as an example, and plots the natural frequency change curves under different dimensionless immersion depths (the ratio of immersion depth to structure radius H/R), and the results are shown in Figure 4.

It can be seen from Figure 4 that the first four-order natural frequency of the shell is significantly affected by the immersion depth when the dimensionless immersion depth H/R is less than 2, and the closer the free liquid surface is, the natural frequency of the shell will increase sharply. After that, the curve changes smoothly, indicating that the influence of the free liquid surface on the coupling system is weakened, and the vibration characteristics of the shell tend to be stable.

In order to determine the specific value of the dimensionless immersion depth H/R in which the flow field can be regarded as an infinite domain, the percentage deviation of the natural frequency of the shell at any immersion depth from the same order natural frequency of the shell in the infinite domain is defined as κ = (f_H − f_∞)/f_∞ × 100, f_H is the natural frequency of the shell at any immersion depth, and f_∞ is the same order natural frequency of the shell in the infinite domain. This paper takes the first 10 natural frequencies of the shell and analyzes the magnitude of the deviation κ of the natural frequency of each order when the dimensionless immersion depth H/R changes. The results are shown in

Table 4. Percent deviations κ of natural frequencies of each order of the shell and infinite domain natural frequencies under different immersion depths.

Order	H/R = 1	H/R = 2	H/R = 3	H/R = 4	H/R = 5
1	5.27	0.67	0.10	0.02	0.01
2	5.25	0.70	0.12	0.03	0.01
3	12.46	0.15	0.01	0.00	0.00
4	12.63	0.15	0.01	0.00	0.00
5	7.09	0.03	0.00	0.00	0.00
6	7.21	0.03	0.00	0.00	0.00
7	4.11	0.10	0.00	0.00	0.00
8	4.15	0.11	0.00	0.00	0.00
9	7.32	0.03	0.00	0.00	0.00
10	7.35	0.03	0.00	0.00	0.00

.

It can be seen from Table 4 that when the dimensionless immersion depth H/R is greater than or equal to 2, the relative percentage deviation κ is all within 1; when H/R is greater than or equal to 4, the relative percentage deviation κ is further reduced to within 0.1. Therefore, it can be considered that when the immersion depth of the cylindrical shell is four times or more than its structural radius, the influence of the free liquid surface on the vibration characteristics of the system can be ignored. That is to say, it can be regarded as a fluid–solid coupling system of a cylindrical shell under an infinite flow field.

In fact, this phenomenon has a clear physical meaning. As the immersion depth increases, the distance between the virtual source and the actual structure will increase by twice the distance away. Therefore, the virtual source fluid load acting on the outer surface of the structure is also rapidly reduced, and the vibration characteristics of the coupled system will naturally tend to the characteristics in infinite domain.

The above research is aimed at the influence of the natural frequency of the free liquid surface. However, the free liquid surface will increase the corresponding natural frequency and will also cause the difference between the positive and negative modal frequencies of the shell that originally appeared in pairs in the infinite domain. In order to explain the reason for this phenomenon, it is also necessary to analyze the mode vibration of the shell. By taking the first four modal circumferential modes in infinite domain or in infinite domain obtained by the method in this paper and FEM as an example, the result is as follows:

It can be seen from Figure 5 that the coupling system still maintains strictly circumferential symmetry in the infinite domain, and the shell circumferential mode is a regular circumferential mode that appears in pairs. The first four circumferential modes are cos(2φ), sin(2φ), cos(3φ), and sin(3φ), the circumferential modes are orthogonal to each other and are not coupled to each other, and with the same circumferential wavenumber n, the natural frequencies of different modes are the same. When the immersion depth is small and the influence of the free liquid surface is not negligible, there is a clear difference of shell mode shape between small immersion depth and the infinite domain, and the circumferential symmetry of the coupled system is broken, so the shell shape is no longer a regular circumferential wave shape. The circumferential mode will become a linear superposition of the mode of different circumferential wave numbers n. However, the non-coupling of the positive and negative modes still exists, that is, the cos system function and the sin system function will not be coupled, so the natural frequencies of the positive and negative modes are different. In addition, the calculation method in this paper is in good agreement with the mode shape of FEM, which further verifies the correctness of the method in this paper.

4.2. The Influence of Shell Parameters on Free Vibration Characteristics

In the analysis of the vibration characteristics of a cylindrical shell structure, the dimensionless axial wave number k_mR is a sensitive physical quantity, and the variable of k_mR can be further defined as the length–diameter ratios L/R. The above conclusion that the watershed can be regarded as infinite when the dimensionless immersion depth H/R is greater than or equal to 4 needs to be verified under different length–diameter ratios L/R. Then, keeping the cylindrical shell radius R = 0.18m and the dimensionless immersion depth H/R = 4 unchanged, and when the aspect ratio L/R is 5, 10, 15, 20, we calculated and obtained the percentage deviation κ of first four natural frequencies of the cylindrical shell and the corresponding natural frequencies in infinite domain; the results are shown in

Table 5. The relative percentage deviation κ of the first four natural frequencies of the shell under different length–diameter ratios L/R.

Order	L/R = 5	L/R = 10	L/R = 15	L/R = 20
1	0.00	0.04	0.05	0.03
2	0.01	0.04	0.05	0.03
3	0.01	0.00	0.05	0.00
4	0.01	0.00	0.05	0.00

.

It can be seen from Table 5 that for different aspect ratios, the percentage relative deviations of the first four natural frequencies of the shell and the corresponding natural frequencies in infinite domain are all less than 0.1. Therefore, it can be considered that for any finite submerged-depth cylindrical shell structure, when the immersion depth is four times or more the radius of the structure, the influence of the free liquid surface can be ignored.

5. Conclusions

Based on the energy variation principle and the image method, this paper proposed an analytical method for solving the natural vibration characteristics of a cylindrical shell with finite submerged depth. First, assuming that the shell is an isotropic thin-walled shell, the Euler beam function under arbitrary boundary conditions was used as the axial displacement function of the shell under the same boundary, and the expression form of each displacement component of the shell was determined. Subsequently, based on the Irrotational hypothesis and the Love shell theory, the kinetic energy and potential energy of the shell structure were obtained, respectively. For the flow field with free liquid surface, assuming that the flow field is non-rotational and non-viscous, the image method was used to derive the fluid velocity potential function considering the influence of the liquid surface. Then, according to the continuous condition of the velocity at the fluid–solid interface, the transformation matrix of the velocity potential function and the shell displacement function was established, and the flow field was added to the entire system in the form of fluid work. According to the comparison with the finite element simulation calculation results, the accuracy of the method in this paper was verified. Further research can draw the following conclusions:

• The existence of the free liquid surface will increase the natural frequency of the shell in the same order mode, and the smaller the immersion depth, that is, the closer the shell is to the free liquid surface, the more obvious the increase in natural frequency.
• When considering the influence of the free liquid surface, due to the destruction of symmetry, there is a clear difference between the shell mode shape and the result in infinite domain. On the one hand, the circumferential waves are coupled, and the shell vibration shape is no longer regular; on the other hand, the natural frequency of the shell’s positive and negative modes is different.
• The influence of free liquid on the shell will quickly diminish as the immersion depth increases. Additionally, when the immersion depth is four times or more than the radius of the shell structure, the influence of the free liquid surface can be ignored.

The structural material studied in this paper is the positive isotropic homogeneous material. In the future, structural materials can be replaced with other materials such as composite materials, functionally graded materials, or carbon nanotubes materials. The method in this paper has a broad application prospect.