Dynamics of pressure pulsations in thin annular gap

Content from this work may be used under the terms of theCreativeCommonsAttribution 3.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI. Published under licence by IOP Publishing Ltd 30th IAHR Symposium on Hydraulic Machinery and Systems IOP Conf. Series: Earth and Environmental Science 774 (2021) 012145 IOP Publishing doi:10.1088/1755-1315/774/1/012145 1 Dynamics of pressure pulsations in thin annular gap Jan K˚ureˇcka, Daniel Himr and Vladim´ır Hab´an Viktor Kaplan Department of Fluid Engineering, Faculty of Mechanical Engineering, Brno University of Technology, Technick´a 2, Brno, CZ E-mail: [email protected] Abstract. Fluid pulsations can significantly influence behavior of bearings or sealing parts of hydraulic machines and therefore their operation. In this contribution, an analytical approach to solve pressure pulsations dynamics is compared with a computational solution using acoustic simulation in ANSYS Mechanical FEM package. Studied problem is simplified to an annular gap and the pressure pulsations are excited either from one annulus end or excitation is made by a movement of an inner cylindrical surface perpendicularly to its axis. Excitation on an annulus surface can have angle dependent phase. This geometry allows a linearization of governing equations and so the equations are solvable analytically. Compressibility was incorporated to analytical and computational approach so the eigen frequencies for both methods were evaluated and compared. Also the added fluid effects were studied and evaluated. 1. Introduction Thin fluid filled gaps such as seals, balancing drums and other similar devices are widely used so there is need for proper understanding of fluid dynamics in these systems. Behavior of the seals influences a rotor-dynamics [2] and the stability of rotating machines is related to it as well [7]. Effect of fluid on dynamics of solid bodies can be mathematically described by incorporating added fluid effects into an equation of motion [8]. Added effects are also covered in publications aimed at flow and vortex related vibrations [9], [11]. These effects include additional damping that can be reason for energy losses in a system [6]. Depending on the geometry evaluation of added mass, stiffness and damping can be a daunting task, for a simple geometry there are formulas and approaches to compute this values analytically [3], [13]. Computational methods include FSI (fluid structure interaction) analyses consisting of computationally intensive structure (finite element method) and fluid (CFD) connected simulations [5], various effects and assessment of different approaches are thoroughly discussed in [4]. This contribution aims to investigate dynamic properties of pressure waves in a thin annular gap with a two different ways of excitation. Analytic method and finite element method in ANSYS software are both based on an acoustic wave equation used to study of resonating fluid cavities [12]. Simple view of added fluid effects is also presented and this parameter is evaluated across frequency spectrum. 2. Case definition This study is inspired by a previous investigation of a balance drum of feedwater pump [10]. Parameters of fluid and geometry are shown in Figure 1 and Table 1. Dimensions are modified 30th IAHR Symposium on Hydraulic Machinery and Systems IOP Conf. Series: Earth and Environmental Science 774 (2021) 012145 IOP Publishing doi:10.1088/1755-1315/774/1/012145 2 dR L Figure 1. Geometry of balancing drum fluid gap, gap thickness dis exaggerated, domain discretization into 12 blocks for a FEM simulation is shown Table 1. Fluid gap parameters parameter value L 1 m R 150 mm d 1 mm (2 mm) ρ905 kg ·m−3 to a rounded values and density value is for a hot water at 180 degrees. There are three variants of boundary conditions: In section 3. pressure excitation on one of the annular surfaces is used and there is zero pressure boundary condition on the opposite end. This is motivated to check FEM analysis modal results precision compared with the analytical solution for respective modeshapes. Section 4. makes use of excitation on inner cylindrical surface and pressure boundary conditions on annular ends. This harmonic analysis of finite length fluid is compared with analytic formula for infinite length cylinder oscillating in confined space [3]. And third setup consists of the same excitation approach as above but Robin radiation boundary conditions are used on both ends. This non reflecting boundary condition substitutes connection of fluid gap to much larger space so the pressure waves do not propagate back. The default rigid wall boundary condition is assigned to outer cylindrical surface for FEM computation. 3. Angle dependent standing wave excitation To compare pressure field obtained from simulation and computed from analytical formula following setup was used. There is condition of zero pressure on one annular end and there are standing sine pressure waves along the circumference on the opposite end. This means that pressure amplitude value peq.(1) - k is the number of waves along the circumference - on each of twelve surfaces is calculated according to its angular position (discretization of annular surface is in Figure 1). In analytical solution, the gap was considered very small and the solution is then two dimensional variable field with angular and longitudinal position. p=pamp ·cos(k·φ) (1) Analytical and computational solution is presented and compared in following sections. Some assumptions are same for all computational and analytic work in this paper. The fluid is inviscid, there is no bulk (second) viscosity either, the cylindrical annular gap has two concentric surfaces (no eccentricity) and there is not any rotational movement on any surface. 30th IAHR Symposium on Hydraulic Machinery and Systems IOP Conf. Series: Earth and Environmental Science 774 (2021) 012145 IOP Publishing doi:10.1088/1755-1315/774/1/012145 3 3.1. Analytic method Analytical solution was made in thesis [10] supervised by V. Hab´an. Wave equation describing pressure change is obtained by modifying continuity and Navier-Stokes equations. ∂ρ ∂t +∂(viρ) ∂xi = 0 (2) ρ∂vi ∂t +ρ∂vi ∂xj vj−∂Πij ∂xj +∂p ∂xi =ρgi(3) To obtain one equation describing pressure field following steps are taken: Continuity equation (2) is modified by partial derivation with respect to tand an assumption that value of speed of sound is one order of magnitude greater than magnitude flow velocity. Speed of sound is introduced into these equations by formula a2=∂p ∂ρ . Navier-Stokes equation (3) modification omits right hand side since there is assumption of no external forces and on left side leaving out effects of viscosity and convective term. This simplified form is then derived with respect to dimension variable. After these modification, equations are combined and a wave equation for pressure is obtained (4). ∂2p ∂t2−a2∂2p ∂xi2= 0 (4) Considering the geometry of the studied domain, the equation can be written in form with longitudinal and angular dimension. ∂2p ∂t2−a2 ∂2p ∂x2 l +∂2p ∂x2 φ!= 0 (5) Solution of this PDE can be obtained by variable separation method. The solution is assumed in the form: p(t, l, φ) = f(t)·g(l)·h(φ) (6) By incorporating (6) into (5) and simplification with respect to time dependency the resulting equation can be rewritten with functions defined in 6. ∂2f(t) ∂t2 f(t)=a2 ∂2g(l) ∂l2 g(l)+a2 r2 ∂2h(φ) ∂φ2 h(φ)(7) From this form the three terms can be separated with respect to variables t, l, φand auxiliary equations solved. The three terms are for further purposes denoted as variables λ,µand ν. These are eigen values of the system and constants c1to c6can be obtained for boundary conditions. p= (c1·e√λt +c2·e−√λt)(c3·e1 a√µl +c4·e− 1 a√µl)(c5·er a√νφ +c6·e−r a√νφ) (8) Formula for angular frequency of respective eigenmodes (9) can be obtained for boundary conditions of zero pressure on both annular ends and for a λ= (i·ω)2Two more variables are used: nwhich is number of amplitude peaks along the length and knumber of pressure waves along the circumference of annulus (k·2π,k= 0 means constant pressure at any angle). Using boundary conditions of zero pressure amplitude on both ends and continuity of pressure along the circumference (h(0) = h(2·π)) following form can be obtained. 30th IAHR Symposium on Hydraulic Machinery and Systems IOP Conf. Series: Earth and Environmental Science 774 (2021) 012145 IOP Publishing doi:10.1088/1755-1315/774/1/012145 4 ωnk =±asn·π L2+k r2 (9) Geometry of annular studied gap is L= 1 m, R= 150 mm and parameter describing speed of sound in water is a= 1500 m·s−1. Corresponding eigenfrequencies can be calculated for a specific nand k. Evaluated frequencies with these annular gap parameters are listed in Table 2. Table 2. Eigenfrequencies ω 2π[Hz] of pressure modes of annular gap (R = 150 mm, L = 1 m) n↓k→0123 1 750 1759 3270 4833 3 2250 2756 3898 5278 Unfolded modeshape of pressure field is depicted in Figure 2. This is visualisation of pressure on cylindrical surface of gap (in analytic solution the gap thickness is zero) in three dimensional form with angular dimension unrolled along the respective axis. Figure 2. Modeshape, frequency 2756 Hz 3.2. FEM modal analysis The computational analysis with the same geometrical and material properties was performed in the next step. Simulations in this article are made in ANSYS Mechanical which is finite element method package with capabilities of acoustic computations [1]. Acoustic wave equation defined from N-S equation and continuity equation is used to solve fluid mechanics in element nodes. Modal analysis was performed in order to match the analytically obtained eigenfrequencies from previous subsection. Geometry and material parameters used in simulation are the same with one exception there is defined thickness of the gap d= 1mm. It was not needed in analytical solution since no such dimension is included in mathematical model. Computational parameters were compressibility of acoustic medium, material parameters identical to the analytical solution, no viscosity or bulk viscosity, boundary conditions of zero 30th IAHR Symposium on Hydraulic Machinery and Systems IOP Conf. Series: Earth and Environmental Science 774 (2021) 012145 IOP Publishing doi:10.1088/1755-1315/774/1/012145 5 pressure at both ends and boundary condition of rigid wall on cylindrical surfaces. Mesh consisted of 15,000 quadratic hexahedral elements (sufficient mesh density). Resulting natural frequencies acquired by acoustic modal analysis of this fluid cavity are summed up and compared with the analytical results in Table 3. Differences are smaller than 0.5%. Table 3. Eigenfrequencies ω 2π[Hz] of pressure modes of annular gap (R = 150 mm, L = 1 m) FEM modal solution n↓k→0 1 2 3 1 750 1755 3260 4817 % diff to analytic 0 0.25 0.3 0.32 3 2250 2753 3889 5264 % diff to analytic 0 0.11 0.22 0.27 3.3. Harmonic analysis with pressure wave excitation Angle dependent pressure pulsation might create dynamic pressure attenuation if the excitation is close to computed eigenfrequencies. This might be similar to pressure field in a hydraulic machines and i.e. excitation from blade passing frequencies so the presented results could be applicable in a similar case. To create a such standing wave modeshape on a surface, boundary condition is prescribed by equation (10) on one annular end of a the domain. p=pamp ·eiωt ·cos(k·φ) (10) Created pressure wave has an amplitude nodes and anti-nodes at positions defined by angular value at respective circumference. Variable kdefines the character of pressure field that can be excited with this boundary condition. Rest of boundary conditions, material and geometrical properties remain without any changes. FEM Simulation was computed using Harmonic analysis solver across frequency interval from zero up to frequency higher than highest evaluated frequency from modal analysis. Spectral resolution was chosen 53 Hz, this value is rather arbitrary, it is high enough for short computation time and low enough for sensible resolution. Block discretization of domain was used for purposes of defining angle dependent excitation boundary condition and for evaluation method used. Domain was split into twelve parts with 30 degrees angle difference. To evaluate force following analytical formula can be used on inner cylindrical surfaces in calculated domain. Variable nsis vector of surface orientation. Fi=−RZL 0Z2π 0 p·nsdφ dl (11) From this equation force components in both X and Y axis (the force in axial direction of cylinder is zero) can be defined and evaluated. This analytical approach gives result that for specific nand ksome modeshapes don’t have any force effect and some exhibit considerable force as a consequence of pressure attenuation due to excitation on a natural frequency. 30th IAHR Symposium on Hydraulic Machinery and Systems IOP Conf. Series: Earth and Environmental Science 774 (2021) 012145 IOP Publishing doi:10.1088/1755-1315/774/1/012145 6 Fi=−R·C·eiωt ZL 0Z2π 0 (eiπn Ll−e−iπn Ll)·(eikφ −e−ikφ)·nsdφ dl (12) 0 500 1,000 1,500 2,000 2,500 3,000 3,500 4,000 4,500 5,000 5,500 10−3 10−1 101 Frequency [Hz] F [N] Figure 3. Force response in direction of excitation on inner surface to pressure excitation (k=1) In Figure 3, curve shows the influence of eigenmode shapes created by chosen excitation. Modes [n;k], [3;0], [3;1], [1;1] are clearly identifiable. Also from visualized pressure in Figure 2 it is obvious that there might be torque acting on the evaluated surface depending on a pressure field shape. This contribution is limited to first eigen frequencies where there is only force. 4. Inner surface excitation, added fluid effect A simple attempt to asses added fluid effects using FEM harmonic analysis is presented in this section. The system is excited by the normal surface velocity of the inner cylindrical surface. This boundary condition simulates up and down linear movement of a rotor. This boundary condition can be described by equation 13. Value of vamp was chosen 0.2 mm·s−1,k= 1. vs=vamp ·eiωt cos(k·φ) (13) The prescribed surface velocity boundary condition can be explained as a velocity wave with nodes and anti-nodes at static positions (block domain discretization is in Figure 1). Inner cylindrical surface on a block located at 0 degrees has an amplitude of vamp, surface at 180 degrees position same amplitude but negative - it moves in a same direction to the 0 degree surface. Surfaces at 90 and 270 degrees have zero amplitude so there is no movement. Any of the twelve inner surfaces has amplitude with respect to (13). This boundary condition should be an acceptable substitution of a prescribed movement that would require dynamic mesh. Results of this simulation with boundary condition of zero pressure on both ends and inner surface excitation can be used to evaluate added fluid effects. From equation of motion with terms of added fluid effects (14) another set of equation can be derived. −F= (ks+kf)·x+ (bs+bf)·v+ (ms+mf)·a(14) There is no solid domain (lower sindex) in the simulation setup so terms of solid parameters can be omitted. The velocity excitation is purely real and acceleration and displacement are purely imaginary. But there is an assumption of no added fluid stiffness due to simplification of fluid properties. So there is real term of damping, imaginary term with added fluid and complex force. This equation can be separated to imaginary and real components. 30th IAHR Symposium on Hydraulic Machinery and Systems IOP Conf. Series: Earth and Environmental Science 774 (2021) 012145 IOP Publishing doi:10.1088/1755-1315/774/1/012145 7 FIm ·eiωt =1 iω mf·vamp ·eiωt (15) FRe ·eiωt =bf·vamp ·eiωt (16) With this set of equations added mass and added damping effects of fluid can be valuated. FEM results of ANSYS harmonic analysis are in form of a imaginary part and real part of force so there is only one unknown variable of fluid effect in both equations. To further explain the separation of terms to real and imaginary phase: It is obvious that the velocity and so the added damping effect is (in terms of complex numbers) in phase with a real component of force. Whereas stiffness (displacement) and mass (acceleration) effects are in phase with imaginary component. Separation of added stiffness and mass might be problematic, but the investigated model without viscosity makes stiffness zero. Using the velocity excitation and zero pressure boundary condition the resulting response of pressure and velocity fields is two dimensional standing wave. Force on the inner surface has only imaginary part and therefore there is no effect of fluid damping. This might not be good enough model to compare with real machinery so one more method was incorporated to simulation. Radiation boundary condition acts as a non reflecting surface through which the pressure and velocity waves can exit the simulated domain. This is depicted in Figure 4, where there are visible acoustic pressure extremes at both ends of domain. Figure 4. Radiation boundary condition used at both ends Comparison of harmonic analysis results and evaluation of added mass effect can be seen in Figure 5. As expected the zero pressure boundary condition and resulting standing wave in the domain has higher peaks than curve representing radiation boundary condition. The latter boundary condition might be more informative with regards to practical application since it has resemblance to a seal gap connected to bigger fluid volume. 30th IAHR Symposium on Hydraulic Machinery and Systems IOP Conf. Series: Earth and Environmental Science 774 (2021) 012145 IOP Publishing doi:10.1088/1755-1315/774/1/012145 8 0 500 1,000 1,500 2,000 2,500 3,000 3,500 4,000 −4 −2 0 2 4·104 Frequency [Hz] mf[kg] radiation B.C. zero pressure B.C. rad. B.C. 2mm zero press. 2 mm Sinyavskii Sinyavskii 2mm Figure 5. Added mass dependence on the frequency To add context to presented values added mass of fluid can be calculated from analytical formula derived in [13], [3]. This equation contains frequency but it is valid for domain of infinite length, therefore the results at frequencies close to eigenfrequencies differ. Nevertheless there is still certain match to computed values and provides valuable insight on computed dataset. Mf=  1 + γ2 1−γ2+ 2 √2 qω·R2 µf 1 + γ3 (1 −γ2)2 ρ·π·R2(17) It is possible to evaluate damping effect of fluid by incorporating radiation boundary condition from Equation (16) on selected surface. In Figure 6, there is visible effect of eigenfrequencies with regards to maximal values of damping coefficient on inner cylindrical surface. 0 500 1,000 1,500 2,000 2,500 3,000 3,500 4,000 102 104 106 108 Frequency [Hz] bf[N·s·m−1] 1mm 2mm Sin. 1mm Sin. 2mm Figure 6. Added damping depndence on frequency 30th IAHR Symposium on Hydraulic Machinery and Systems IOP Conf. Series: Earth and Environmental Science 774 (2021) 012145 IOP Publishing doi:10.1088/1755-1315/774/1/012145 9 bf= 2 √2 qω·R2 µf 1 + γ3 (1 −γ2)2ρ·π·R2·ω(18) Similar to added mass, an analytical formula (18) for added damping was used from cited papers [3]. Values are limited only to low frequencies, since they differ more than mass because there are differences in both approaches. The FEM added damping does not incorporate viscosity and analytic formula does. There is influence of finite length of FEM domain - difference in modeshapes and eigenfrequencies. Analytic values are still useful as a reference value for comparison 5. Conclusion Analytic method was presented to describe pressure field in a thin annular gap and FEM analysis in ANSYS mechanical was used for comparison. Both methods are in an agreement and our next step is experimental verification. Even so, further examination of dynamics of this system with a different excitation method provides an ability to explore added fluid effects. This way of determination added fluid effects is very quick, since the calculation is not expensive in terms of processing power. This is advantage in comparison to an unsteady CFD simulation which takes significantly more resources and more so for FSI analysis. Choice of investigated geometry was based on a similarity with a previously studied balance drum of hydraulic pump. The results and findings in this paper could be applied to sealing gaps of hydraulic machinery rotors, fluid bearings and other instances of fluid filled thin cavities that could be susceptible to excitation. Further focus to study and to extend this topic would be modify analytical solution to different type of excitation, more complex geometry and travelling waves. Acknowledgments Authors acknowledge the support of their work from following projects: ”Computer Simulations for Effective Low-Emission Energy” funded as project No. CZ.02.1.01/0.0/0.0/16 026/0 008392 by Operational Programme Research, Development and Education, Priority axis 1: Strengthening capacity for high-quality research and from specific research project FSI-S-206235 funded by Faculty of Mechanical Engineering of Brno University of Technology.