CN104102783A - Method for forecasting numerical value of cavitation of underwater propeller tip vortex - Google Patents

Abstract

The invention discloses a method for forecasting numerical value of cavitation of an underwater propeller tip vortex and belongs to the field of optimized design of propellers. The method includes steps of (a) determining a basic grid; (b) determining a precision grid; (c) determining an optimal grid; (d) forecasting the cavitation of the underwater propeller tip vortex under the required working condition by means of the optimal grid. The method is used for forecasting the numerical value of the cavitation of the propeller tip vortex, and effectively forecasts tip-vortex cavitation of the E779A type propeller under different working conditions via comparison of results of related experiments, therefore, the method is of great significance in forecasting and evaluating of cavitation performance in the propeller design, design cost and period are effectively reduced, and the method is quite good in application prospect.

Description

A kind of Numerical Prediction Method about underwater propeller tip vortex cavitation

Technical field

The present invention relates to screw propeller optimal design field, specifically, a kind of Numerical Prediction Method about underwater propeller tip vortex cavitation.

Background technology

Propeller cavitation not only can reduce propeller performance, produces cavitation and degrades, and causes ship hull vibration, and can produce the cavitation noise that is enough to stick one's chin out.Along with modern ships is more and more higher to the requirement of load and ship's speed, propeller cavitation phenomenon is difficult to avoid.And tip vortex cavitation is usually the cavitation that screw propeller occurs the earliest, the accurate forecast of its phenomenon can provide important evidence to the judgement of cavitation inception.The numerical forecasting of tip vortex cavitation is also one of cavitation numerical forecasting of the most difficult realization.Numerical value is caught tip vortex cavity, and around pressure, speed and the flow field structure in flow field can further strengthen the understanding to mechanism of noise generation in interior comprehensive cavitating flow information to obtain it.Tip vortex cavitation be always numerous relevant scholars study hotspot and difficult point, be also one of focus of paying close attention to of impeller design manufacture and even shipbuilding industry.

At present, the research of cavitating flow is mainly contained to experiment and two kinds of methods of numerical analysis.The deficiency such as rely on that experiment exists that cost is time-consuming completely, model propeller amendment inconvenience and the installation meeting stream field of sensor in flow field impact.In recent years, along with the fast development of computer technology, the raising greatly of Computing ability, making to utilize modern computational to carry out numerical forecasting to tip vortex cavitation becomes possibility, and numerical forecasting becomes the important means of cavitation forecast.Be comparatively successfully sheet cavitation to propeller cavitation numerical simulation now, and tip vortex cavitation (tip vortex cavitation) numerical forecasting difficulty is larger.Conventional forecasting procedure is at present, by establishing rational computational fields space, set up high-quality grid, select suitable turbulence model and cavitation model, based on mixing the numerical solution that solves eddy stress RANS equation on the basis of multiphase flow model and obtain flow field around screw propeller, can substantially realize the numerical forecasting of propeller cavitation, but the value of forecasting of tip vortex cavitation is not obvious.

At present, abroad to tip vortex cavitation, research mainly lays particular emphasis on hydrofoil.Dynaflow Inc. company of the U.S. combines spherical Bubble dynamics model (Rayleigh – Plesset) and the non-spherical Bubble dynamics model revised, and be built in UnRANS equation solver and solve, the nascent phenomenon of tip vortex cavitation is carried out to numerical forecasting, and they consider size and the impact of space distribution on cavitation inception of flow field cavitation core.South Korea Seoul university research personnel adopt Eulerian – Lagrangian method, and consider that equally the Size Distribution feature of cavity core carried out numerical analysis to tip vortex cavitation and cavitation noise thereof.Norway researchist utilizes business software Fluent software to carry out numerical simulation to the tip vortex cavitation of NACA hydrofoil, and its result and experimental result are more consistent.In addition, Sweden researchist utilizes Large eddy simulation method to carry out numerical forecasting to propeller cavitation, and its screw propeller tip vortex cavitation forecast result and experimental phenomena are more consistent, but sheet cavitation forecast area obviously shows that than experiment result is large.

So, be badly in need of a kind of Numerical Prediction Method that can effectively forecast screw propeller tip vortex cavitation.

Summary of the invention

Screw propeller tip vortex cavitation principal character has:

1. in the various cavitation forms of screw propeller, tip vortex cavitation usually occurs the earliest, and is accompanied by cavitation inception phenomenon;

2. in tail flow field, near tip of propeller blade, produce, present vertically helix shape and distribute, and vestige is longer;

3. there is strong eddying motion in tip whirlpool centre line zone, and produce strong low pressure;

4. tip whirlpool produces at first near blade tip, is departing from after blade, and the helix radius in tip whirlpool can shrink;

5. tip whirlpool can be carried out smooth to wake flow direction along with the vortex sheet structure of screw current;

6. the noise that tip vortex cavitation sends is generally distributed in medium-high frequency section, and its frequency spectrum presents continuous feature.

Principle of the present invention is exactly the principal character according to above-mentioned tip vortex cavitation, Integrated using cavitation model, turbulence model and mixing two-phase flow model carry out numerical solution and receive dimension-Stokes (N-S) system of equations, obtain every physical parameter value such as the vapour phase volume fraction of pressure, speed and reflection cavitation in screw propeller fluid environment, the invention provides a kind of Numerical Prediction Method about underwater propeller tip vortex cavitation, can carry out prediction and evaluation to Design of Propeller cavitation performance.

The present invention adopts following technical scheme:

About a Numerical Prediction Method for underwater propeller tip vortex cavitation, comprise the following steps:

(a) basic grid is determined;

(b) fine grid blocks is determined;

(c) Bestgrid is determined;

(d) utilize Bestgrid to forecast the tip vortex cavitation under required working condition.

The specific implementation of the inventive method comprises the following steps:

Step 1, utilizes modeling software to set up screw propeller 3-D geometric model, and is imported grid division software.

Step 2, divide in software and set up three kinds of alternative grids at grid, and imported calculation procedure:

The computational fields of three kinds of alternative grids is identical, computational fields is cylindrical, its speed frontier distance propeller center that becomes a mandarin is 1.5D, D is airscrew diameter, downstream pressure outlet frontier distance is 5D, propeller center is 2.5D to side periphery distance, and the number of grid of three kinds of alternative grids increases gradually, is about respectively 2,000,000,3,000,000 and 4,000,000.

Step 3, cavitation model and turbulence model are set:

Adopt full cavitation model and Renormalization Group turbulence model (RNG k-ε turbulence model), and its important parameter is revised, in correction to transformation ratio parameter in cavitation model and turbulence model, the correction of turbulent viscosity coefficient adopts C language compilation, and recycling macro call (DEFINE_TURBULENT_VISCOSITY etc.) form embeds calculation procedure.

Step 4, numerical evaluation setting parameter:

Numerical parameter is set the correlation parameter setting that comprises working condition, boundary condition and numerical algorithm; Working condition is mainly set screw propeller rotational speed, environmental pressure and inflow velocity value, thus determine screw propeller dimensionless group, i.e. advance coefficient (J) and cavitation number (σ _n); Set for boundary condition, speed inlet boundary adopts inflow velocity value, and far field boundary condition adopts inflow velocity value to set equally, and the top hole pressure at downstream pressure outlet interface is set to static pressure; For numerical algorithm, it is discrete that in Na Wei-Stokes (N-S) equation, convective term adopts Second-order Up-wind form, diffusion term adopts Using Second-Order Central Difference form discrete, velocity pressure coupling adopts the SIMPLE algorithm that is applicable to non-structured grid, use pointwise Gauss-Seidel iterative discrete equation, utilize the convergence of algebraic multigrid speed-up computation, adopt Sliding mesh computing technique for non-permanent calculating, improve the accuracy of calculating.Because multiphase flow model, cavitation model and Sliding mesh calculating are larger to computer resource usage, adopt parallel computing to shorten computing time;

Mass conservation continuity (continuity) residual error convergence is three rank, in equation, other physical quantity residual error convergence is quadravalence, the initial value that utilizes the convergence solution of single-phase flow to solve as multinomial stream, and the initial value calculating steady state solution as unstable state, in order to ensure calculating convergence, screw propeller rotational speed is progressively increased to predetermined value, and suitably dwindles relaxation factor.

Navier-Stokes equation, also claims to receive dimension-RANS, is called for short N-S equation.Expression is the N-S equation of tensor form, wherein subscript k and j, the same denotation coordination axle of i.

$\begin{array}{c}\frac{\∂}{\∂t}\left(\mathrm{\ρ}{u}_i\right)+\frac{\∂}{\∂x_i}\left(\mathrm{\ρ}{u}_i{u}_j\right)=\mathrm{\ρ}{f}_i-\frac{\∂p}{\∂x_i}+\frac{\∂}{\∂x_j}\mathrm{\μ}[(\frac{\∂u_i}{\∂x_j}+\frac{\∂u_j}{\∂x_i})-\frac{2}{3}{\mathrm{\δ}}_{\mathrm{ij}}\frac{\∂u_k}{\∂x_k}\\ =\mathrm{\ρ}{f}_i-\frac{\∂p}{\∂x_i}+\frac{\∂}{\∂x_i}\left(\mathrm{\μ}\frac{\∂u_i}{\∂x_j}\right)+\frac{1}{3}\frac{\∂}{\∂x_j}\left(\mathrm{\μ}\frac{\∂u_k}{\∂x_k}\right)\end{array}$

In formula, the kinetic viscosity that μ is fluid, δ _ijfor " Kronecker Delta " tensor.

Step 5, carry out numerical evaluation:

Because cavitation model adds after N-S equation, the stability of calculating reduces, and easily occurs unusual appearance.Therefore,, in order to make numerical evaluation steadily carry out, adopt computation process step by step step by step.Specifically, in screw propeller duty parameter, environmental pressure and inflow velocity can directly be set to operating mode value, and revolution speed of propeller adopts classification to increase, until be increased to predetermined operating mode value.In addition, first calculate non-cavitating model Flow Field Distribution, by the time after calculation stability, open again cavitation model.First the parameters such as pressure, density, momentum and vapour phase mark are carried out to single order precision discrete scheme and calculate, after calculation stability, more discrete precision is brought up to second order or QUCIK etc.For the stability that ensures that second order calculates, sub-relaxation factor is suitably reduced.The isoparametric sub-relaxation factor of pressure, momentum, vapour phase mark, Turbulent Kinetic, turbulence dissipation rate and turbulent flow stickiness is set as respectively: 0.25,0.6,0.2,0.7,0.7,0.9.

Step 6, establish basic grid according to hydrodynamic force numerical result:

Utilize these three kinds of alternative grids to carry out numerical evaluation to Propeller, and result of calculation is analyzed; After tending to be steady after the result of calculation of hydrodynamic parameter (thrust coefficient and moment coefficient) increases along with grid number, think that numerical evaluation is along with grid increases and tends towards stability, after selected calculation stability the minimum grid of number of grid as basic grid.

Step 7, imports grid by the file of basic grid and divides software.

Step 8, on the basis of basic grid, according to tip whirlpool shape, utilize grid to divide software and set up tip whirlpool area grid:

First, in professional modeling tool, set up the conduit region geometric model of a helix shape conforming to model parameter according to isometric helix line model, this region is the tip vortex cavitation generation area of supposition, again helix geometry models is imported in basic grid, and carry out grid division, equidistant helix shape area grid is from closing on blade blade tip, equidistant helix radius is preset as 0.82R, the unit grid size in helix region is about 0.0001D, and in former basic grid, the unit size on every limit is constant.

Step 9, the cavitation model of repeating step three and turbulence model are set, the numerical evaluation setting parameter of step 4 and step 5 carry out numerical evaluation.

Step 10, establish fine grid blocks according to tip vortex cavitation numerical forecasting result:

This step is mainly the important parameter r establishing in isometric helix line model.The length of helix radius r is about (0.81-0.83) R, according to tip vortex cavitation numerical result adjust and finally determine this parameter, concrete grammar is, when tip vortex cavitation in the result of calculation of step 9 is arranged in this hypothesis district (being the helix pipeline region that r that step 8 sets sets up is exactly tip vortex cavitation generation area) substantially, fine grid blocks is established, and can carry out next step; Otherwise, when tip vortex cavitation numerical forecasting is not during substantially at this hypothesis district, according to site error, parameter r is adjusted, and return to step 8 and step 9 recalculates, until result conforms to substantially, establish fine grid blocks.

Step 11, imports grid by fine grid blocks file and divides software.

Step 12, set up wake flow apart from the other two kinds of computational fields grids that are respectively 7D and 9D:

On the basis of fine grid blocks, set up other two kinds of computational fields, its wake flow distance is respectively 7D and 9D, and originally in fine grid blocks, the unit size on every limit is constant, divide the wake flow subregion grid increasing, its grid cell is of a size of former fine grid blocks in 5D place grid cell size.

Step 13, the cavitation model of repeating step three and turbulence model are set, the numerical evaluation setting parameter of step 4 and step 5 carry out numerical evaluation.

Step 14, the result that is 5D with wake flow distance compares definite Bestgrid.

These three kinds of computational fields grids are carried out analyzing after numerical evaluation, when along with wake flow is apart from increase, when tip vortex cavitation numerical forecasting precision no longer obviously increases, establish best wake flow distance and Bestgrid.

Step 15, utilizes Bestgrid and numerical evaluation parameter to carry out numerical forecasting to required operating mode.

Further, described full cavitation model and parameter thereof is modified to:

Work as p<p _vtime, steam generation rate is:

$R_e=\mathrm{Ce}\frac{k}{\mathrm{\&gamma;}}{\mathrm{\&rho;}}_l{\mathrm{\&rho;}}_v\sqrt{\frac{2}{3}\frac{p_v-p}{{\mathrm{\&rho;}}_l}}(1-f_v)$

Work as p>p _vtime, vapour phase becomes liquid phase, steam solidification rate R _cfor:

$R_c=\mathrm{Cc}\frac{k}{\mathrm{\&gamma;}}{\mathrm{\&rho;}}_l{\mathrm{\&rho;}}_v\sqrt{\frac{2}{3}\frac{p-p_v}{{\mathrm{\&rho;}}_l}}{f}_v$

Wherein, f _v=α _vρ _v/ ρ _mfor vapour phase massfraction, vaporization coefficient Ce=0.02 and condensation coefficient Cc=0.01 are empirical parameter.

Adopt in transformation ratio (Re and Rc) expression formula according to dimensional analysis k instead of .Under FLUENT software environment, can utilize self-defining function UDF to cavitation model in parameter transformation ratio Re revise, revise after adopting C language compilation and call in calculation procedure.

Further, described Renormalization Group turbulence model and parameter thereof is modified to:

Renormalization Group turbulence model RNG k-ε is the model that the mathematical method of transient state N-S Renormalization Group (Renormalization Group is called for short RNG) for equation is derived.It,, by embody the impact of small scale in Large Scale Motion item and correction viscosity item, is removed from governing equation and make these small scale kinematic systems.Its k equation and ε equation are respectively:

$\frac{\&PartialD;}{\&PartialD;t}\left({\mathrm{\&rho;}}_{m}k\right)+\frac{\&PartialD;}{\&PartialD;x_j}\left({\mathrm{\&rho;}}_{m}k{u}_{\mathrm{mj}}\right)=\frac{\&PartialD;}{\&PartialD;x_j}\left[\left({\mathrm{\&alpha;}}_k\mathrm{\&mu;}\right)\frac{\&PartialD;k}{\&PartialD;x_j}\right]+G-{\mathrm{\&rho;}}_m\mathrm{\&epsiv;}$

$\frac{\&PartialD;}{\&PartialD;t}\left({\mathrm{\&rho;}}_m\mathrm{\&epsiv;}\right)+\frac{\&PartialD;}{\&PartialD;x_j}\left({\mathrm{\&rho;}}_m\mathrm{\&epsiv;}{u}_{\mathrm{mj}}\right)=\frac{\&PartialD;}{\&PartialD;x_j}\left[\left({\mathrm{\&alpha;}}_{\mathrm{\&epsiv;}}\mathrm{\&mu;}\right)\frac{\&PartialD;\mathrm{\&epsiv;}}{\&PartialD;x_j}\right]-R+C_{1\mathrm{\&epsiv;}}\frac{\mathrm{\&epsiv;}}{k}G-C_{2\mathrm{\&epsiv;}}{\mathrm{\&rho;}}_m\frac{{\mathrm{\&epsiv;}}^2}{k}$

In formula, Turbulent Kinetic dissipative shock wave (Turbulent Dissipation Rate) the a reciprocal of effective turbulent prandtl number of Turbulent Kinetic k and dissipative shock wave ε _k=a _ε=1.39; Model parameter C _{1 ε}=1.47, C _{2 ε}=1.68; Viscosity coefficient is μ=μ _t+ μ _m, μ _mfor mixed flow coefficient of viscosity; Revise turbulent viscosity coefficient μ _t=[ρ _v+ α _l ¹⁰(ρ _l-ρ _v)] C _μk ²/ ε, C _μ=0.085, make its non-permanent two-phase simulated flow that is more suitable for high reynolds number, thereby can simulate better propeller cavitation.

Further, described equidistant helix mathematical model and parameter thereof is as follows:

$\left\{\begin{array}{c}x=r\cos \left(\frac{2\mathrm{\&pi;}}{k}z\right)\\ y=r\sin \left(\frac{2\mathrm{\&pi;}}{k}z\right)\end{array}\right.$

Wherein, x, y, z is respectively three coordinate axis of cartesian coordinate system, and constant r is helix radius, and the length of r is about (0.81-0.83) R, and R is propeller radius; Constant k is the distance that often rotates a circle and advance on x axle of helix, k=U _∞/ n, U _∞for inflow velocity, n is revolution speed of propeller.

Further, the method that described grid is divided is: adopt subregion mixed mesh method grid division: screw propeller around flow field regions adopts non-structured grid method to divide, grid is reduced to blade tip gradually by propeller hub, blade tip place surface grids is triangle, size is about 0.001D Jiang Grains place for 0.015D, and D is airscrew diameter; Adopt boundary layer grid at blade surface, boundary layer grid has 4 layers, and its adjacent two layers aspect ratio is 1.1, and ground floor grid cell height is about 0.0008D, makes dimensionless group 20<y ⁺<300, adopts structured grid to divide the computational fields of the peripheral regular shape of screw propeller.

Method of the present invention can realize in general general CFD fluid calculation software (CFX, FLUENT etc.), and grid is divided softwares such as can adopting GAMBIT and realized.First the screw propeller digital model of primary design is imported to grid and divide software, and carry out grid division according to the method in the present invention.Grid model imports computing platform, carries out numerical evaluation, and numerical result is returned to Design of Propeller personnel according to design conditions, carries out cavitation performance assessment, then determines whether designed screw propeller reaches designing requirement.In addition, the present invention also utilizes autoexec to carry out parallel numerical calculating in operating system platform.

Beneficial effect:

(1) main technologies that affects tip vortex cavitation numerical result in tip vortex cavitation forecast has: model (turbulence model and cavitation model), numerical algorithm and grid, because model and numerical algorithm have stronger adaptability and robustness to the numerical evaluation of different objects (oar mould and operating mode), and grid has larger difference to the numerical evaluation of different objects.The present invention passes through model, and the essential condition such as algorithm and grid combine, and find an effective technological approaches to forecast better tip vortex cavitation, and the present invention has significant application value and application prospect to the improvement of Design of Propeller cavitation performance.

(2) the present invention is based on viscous flow theory, Integrated using cavitation model, turbulence model and mix two-phase flow model and carry out numerical solution and receive dimension-Stokes (N-S) system of equations, obtain the pressure in screw propeller fluid environment, every physical parameter value such as the vapour phase volume fraction of speed and reflection cavitation, the technology of the present invention is based on recent development modern computational and Bubble dynamics rapidly, by grid generation technique, numerical computation method and parallel computing are introduced the numerical forecasting of screw propeller tip vortex cavitation, embody the amalgamation of multidisciplinary and multi-field technology.

(3) because hydrofoil does linear translational motion at fluid, screw propeller is to rotate in fluid, and around it, flow field situation is than hydrofoil complexity.The present invention is directed to screw propeller tip vortex cavitation and carry out numerical forecasting; By contrasting with related experiment result, the present invention has more efficiently forecast the tip vortex cavitation of E779A type screw propeller under several different operating modes, therefore, the present invention has vital role to the prediction and evaluation of Design of Propeller cavitation performance, can effectively reduce design cost and design cycle, its application prospect is very good.

Brief description of the drawings

Fig. 1 is the schematic diagram of full runner screw propeller computational fields of the present invention;

Fig. 2 (a) is blade surface schematic diagram;

Fig. 2 (b) is the boundary layer grid schematic diagram at A place in Fig. 2 (a);

Fig. 3 is near the schematic diagram that blade tip, grid distributes;

Fig. 4 is E779A screw propeller geometric model;

Fig. 5 (a) determines process flow diagram for basic grid of the present invention;

Fig. 5 (b) determines process flow diagram for fine grid blocks of the present invention;

Fig. 5 (c) is that Bestgrid of the present invention is established process flow diagram;

Fig. 5 (d) is main flow chart of the present invention;

Fig. 6 is that in blade wake, tip whirlpool seemingly presents helix shape distribution schematic diagram near blade tip

Fig. 7 is the contrast of tip vortex cavitation numerical forecasting of the present invention and experimental result.

Fig. 8 (a) is the tip vortex cavitation numerical simulation result within 1/8 cycle under 0 ° of position E779A oar mould nonlinear inflow condition;

Fig. 8 (b) is the tip vortex cavitation numerical simulation result within 1/8 cycle under 18 ° of position E779A oar mould nonlinear inflow conditions;

Fig. 8 (c) is the tip vortex cavitation numerical simulation result within 1/8 cycle under 36 ° of position E779A oar mould nonlinear inflow conditions;

Fig. 8 (d) is the tip vortex cavitation numerical simulation result within 1/8 cycle under 54 ° of position E779A oar mould nonlinear inflow conditions;

Fig. 8 (e) is the tip vortex cavitation numerical simulation result within 1/8 cycle under 72 ° of position E779A oar mould nonlinear inflow conditions;

Fig. 8 (f) is the tip vortex cavitation numerical simulation result within 1/8 cycle under 81 ° of position E779A oar mould nonlinear inflow conditions.

In figure: 1, the speed border that becomes a mandarin, 2, border, far field, 3, pressure export border.

Embodiment

Below in conjunction with drawings and Examples, the present invention is further detailed explanation.

Embodiment

As shown in Fig. 5 (a), Fig. 5 (b), Fig. 5 (c) and Fig. 5 (d), a kind of Numerical Prediction Method about underwater propeller tip vortex cavitation, specific implementation comprises the following steps:

Step 1, utilizes professional software to set up screw propeller 3-D geometric model, and is imported grid division software.

Step 2, divides in software and sets up three kinds of alternative grids at grid, and imported calculation procedure.

The computational fields of three kinds of alternative grids is identical, according to screw propeller Field Characteristics around, set up cylindrical computational fields, as shown in Figure 1, its speed frontier distance propeller center that becomes a mandarin is 1.5D, and D is airscrew diameter, downstream pressure outlet frontier distance is 5D, propeller center is 2.5D to side periphery distance, and the number of grid of three kinds of alternative grids increases gradually, is about respectively 2,000,000,3,000,000 and 4,000,000; The method that grid is divided is: adopt subregion mixed mesh method grid division: screw propeller around flow field regions adopts non-structured grid method to divide, grid is reduced to blade tip gradually by propeller hub, blade tip place surface grids (as shown in Figure 3) is triangle, size is about 0.001D Jiang Grains place for 0.015D, and D is airscrew diameter; Because cavitation is mainly distributed in blade face, therefore this area grid quality requirements is higher, in order to adapt to better Wall-function, adopt boundary layer grid at blade surface, adopt boundary layer grid to improve the forecast precision to tip vortex cavitation, as shown in Fig. 2 (a), blade surface schematic diagram, Fig. 2 (b) is the boundary layer grid schematic diagram at A place in Fig. 2 (a), blade surface boundary layer grid has 4 layers, its adjacent two layers aspect ratio is 1.1, and ground floor grid cell height is about 0.0008D, makes dimensionless group 20<y ⁺<300, adopts structured grid to divide the computational fields of the peripheral regular shape of screw propeller.

Step 3, cavitation model and turbulence model are set:

Adopt full cavitation model and Renormalization Group turbulence model (RNG k-ε turbulence model), and its important parameter is revised, in correction to transformation ratio parameter in cavitation model and turbulence model, the correction of turbulent viscosity coefficient adopts C language compilation, and recycling macro call (DEFINE_TURBULENT_VISCOSITY etc.) form embeds calculation procedure;

Full cavitation model and parameter thereof are modified to:

Work as p<p _vtime, steam generation rate is:

$R_e=\mathrm{Ce}\frac{k}{\mathrm{\&gamma;}}{\mathrm{\&rho;}}_l{\mathrm{\&rho;}}_v\sqrt{\frac{2}{3}\frac{p_v-p}{{\mathrm{\&rho;}}_l}}(1-f_v)$

In the time of p>pv, vapour phase becomes liquid phase, and steam solidification rate Rc is:

$R_c=\mathrm{Cc}\frac{k}{\mathrm{\&gamma;}}{\mathrm{\&rho;}}_l{\mathrm{\&rho;}}_v\sqrt{\frac{2}{3}\frac{p-p_v}{{\mathrm{\&rho;}}_l}}{f}_v$

Wherein, f _v=α _vρ _v/ ρ _mfor vapour phase massfraction, vaporization coefficient Ce=0.02 and condensation coefficient Cc=0.01 are empirical parameter.

Adopt in transformation ratio (Re and Rc) expression formula according to dimensional analysis k instead of under FLUENT software environment, can utilize self-defining function UDF to cavitation model in parameter transformation ratio Re revise, revise after adopting C language compilation and call in calculation procedure;

Renormalization Group turbulence model RNG k-ε is the model that the mathematical method of transient state N-S Renormalization Group (Renormalization Group is called for short RNG) for equation is derived.It,, by embody the impact of small scale in Large Scale Motion item and correction viscosity item, is removed from governing equation and make these small scale kinematic systems.Its k equation and ε equation are respectively:

$\frac{\&PartialD;}{\&PartialD;t}\left({\mathrm{\&rho;}}_{m}k\right)+\frac{\&PartialD;}{\&PartialD;x_j}\left({\mathrm{\&rho;}}_{m}k{u}_{\mathrm{mj}}\right)=\frac{\&PartialD;}{\&PartialD;x_j}\left[\left({\mathrm{\&alpha;}}_k\mathrm{\&mu;}\right)\frac{\&PartialD;k}{\&PartialD;x_j}\right]+G-{\mathrm{\&rho;}}_m\mathrm{\&epsiv;}$

$\frac{\&PartialD;}{\&PartialD;t}\left({\mathrm{\&rho;}}_m\mathrm{\&epsiv;}\right)+\frac{\&PartialD;}{\&PartialD;x_j}\left({\mathrm{\&rho;}}_m\mathrm{\&epsiv;}{u}_{\mathrm{mj}}\right)=\frac{\&PartialD;}{\&PartialD;x_j}\left[\left({\mathrm{\&alpha;}}_{\mathrm{\&epsiv;}}\mathrm{\&mu;}\right)\frac{\&PartialD;\mathrm{\&epsiv;}}{\&PartialD;x_j}\right]-R+C_{1\mathrm{\&epsiv;}}\frac{\mathrm{\&epsiv;}}{k}G-C_{2\mathrm{\&epsiv;}}{\mathrm{\&rho;}}_m\frac{{\mathrm{\&epsiv;}}^2}{k}$

In formula, Turbulent Kinetic dissipative shock wave (Turbulent Dissipation Rate) the a reciprocal of effective turbulent prandtl number of Turbulent Kinetic k and dissipative shock wave ε _k=a _ε=1.39; Model parameter C _{1 ε}=1.47, C _{2 ε}=1.68; Viscosity coefficient is μ=μ _t+ μ _m, μ _mfor mixed flow coefficient of viscosity; Revise turbulent viscosity coefficient μ _t=[ρ _v+ α _l ¹⁰(ρ _l-ρ _v)] C _μk ²/ ε, C _μ=0.085, make its non-permanent two-phase simulated flow that is more suitable for high reynolds number, thereby can simulate better propeller cavitation.

Step 4, numerical evaluation setting parameter:

Numerical parameter is set the correlation parameter setting that comprises working condition, boundary condition and numerical algorithm; Working condition is mainly set screw propeller rotational speed, environmental pressure and inflow velocity value, thus determine screw propeller dimensionless group, i.e. advance coefficient (J) and cavitation number (σ _n); Set for boundary condition, speed inlet boundary adopts inflow velocity value, and far field boundary condition adopts inflow velocity to set equally, and the top hole pressure at downstream pressure outlet interface is set to static pressure; For numerical algorithm, it is discrete that in Na Wei-Stokes (N-S), convective term adopts Second-order Up-wind form, diffusion term adopts Using Second-Order Central Difference form discrete, velocity pressure coupling adopts the SIMPLE algorithm that is applicable to non-structured grid, use pointwise Gauss-Seidel iterative discrete equation, utilize the convergence of algebraic multigrid speed-up computation, adopt Sliding mesh computing technique for non-permanent calculating, improve the accuracy of calculating.Because multiphase flow model, cavitation model and Sliding mesh calculating are larger to computer resource usage, adopt parallel computing to shorten computing time;

Mass conservation continuity (continuity) residual error convergence is three rank, in equation, other physical quantity residual error convergence is quadravalence, the initial value that utilizes the convergence solution of single-phase flow to solve as multinomial stream, and the initial value calculating steady state solution as unstable state, in order to ensure calculating convergence, screw propeller rotational speed is progressively increased to predetermined value, and suitably dwindles relaxation factor;

Step 5, carry out numerical evaluation:

Because cavitation model adds after N-S equation, the stability of calculating reduces, and easily occurs unusual appearance.Therefore,, in order to make numerical evaluation steadily carry out, adopt computation process step by step step by step.Specifically, in screw propeller duty parameter, environmental pressure and inflow velocity can directly be set to operating mode value, and revolution speed of propeller adopts classification to increase, until be increased to predetermined operating mode value.In addition, first calculate non-cavitating model Flow Field Distribution, by the time after calculation stability, open again cavitation model.First the parameters such as pressure, density, momentum and vapour phase mark are carried out to single order precision discrete scheme and calculate, after calculation stability, more discrete precision is brought up to second order or QUCIK etc.For the stability that ensures that second order calculates, sub-relaxation factor is suitably reduced.The isoparametric sub-relaxation factor of pressure, momentum, vapour phase mark, Turbulent Kinetic, turbulence dissipation rate and turbulent flow stickiness is set as respectively: 0.25,0.6,0.2,0.7,0.7,0.9.

Step 6, establish basic grid according to hydrodynamic force numerical result:

Utilize these three kinds of alternative grids to carry out numerical evaluation to Propeller, and result of calculation is analyzed; After tending to be steady after the result of calculation of hydrodynamic parameter (thrust coefficient and moment coefficient) increases along with grid number, think that numerical evaluation is along with grid increases and tends towards stability, after selected calculation stability the minimum grid of number of grid as basic grid.

Step 7, imports grid by the file of basic grid and divides software.

Step 8, on the basis of basic grid, according to tip whirlpool shape, utilize grid to divide software and set up tip whirlpool area grid:

Here key technical problem is geometric model how reasonably to set up a reflection screw propeller tip vortex cavitation region, and this is also one of important step of the inventive method.Specifically, in wake flow, distribute and present approx helix shape feature according to tip vortex cavitation, set up the computational fields of isometric helix wire shaped in screw current region, centered by helix, divide a long and narrow computational fields space.Equidistant helix shape area grid, from closing on blade blade tip, and extends to wake flow direction.The parameter of equidistant helix is relevant to screw propeller geometric shape parameters and duty parameter.Concrete equidistant helix mathematical model and parameter thereof are as follows:

$\left\{\begin{array}{c}x=r\cos \left(\frac{2\mathrm{\&pi;}}{k}z\right)\\ y=r\sin \left(\frac{2\mathrm{\&pi;}}{k}z\right)\end{array}\right.$

Wherein, x, y, z is respectively three coordinate axis of cartesian coordinate system, and constant r is helix radius, and the length of r is about (0.81-0.83) R, and R is propeller radius; Constant k is the distance that often rotates a circle and advance on x axle of helix, k=U _∞/ n, U _∞for inflow velocity, n is revolution speed of propeller; Fig. 6 seemingly presents helix shape distribution schematic diagram for tip whirlpool in the blade wake providing in external certain list of references near blade tip.

Tip vortex cavitation area grid is made concrete grammar: first, in professional modeling tool, set up the conduit region geometric model of a helix shape conforming to model parameter according to above-mentioned isometric helix line model, this region is the tip vortex cavitation generation area of supposition, again helix geometry models is imported in basic grid, and carry out grid division, equidistant helix shape area grid is from closing on blade blade tip, equidistant helix radius is preset as 0.82R, the unit grid size in helix region is about 0.0001D, in former basic grid, the unit size on every limit is constant.

Step 9, the cavitation model of repeating step three and turbulence model are set, the numerical evaluation setting parameter of step 4 and step 5 carry out numerical evaluation.

Step 10, establish fine grid blocks according to tip vortex cavitation numerical forecasting result:

Fine grid blocks refers to, on the basis of basic grid, tip vortex cavitation generation area be carried out to fine grid blocks division, is mainly the mesh-density that improves this region.

The main important parameter r establishing in isometric helix line model in this step.The length of helix radius r is about (0.81-0.83) R, here do not have this parameter of the disposable establishment of fine method, only have according to tip vortex cavitation numerical result adjust and finally determine this parameter, first adopt specifically the grid of setting up in step 8 (being the tip vortex cavitation area grid of the 0.82R supposition of setting up according to r), and through the numerical evaluation of step 9.When tip vortex cavitation in result of calculation is positioned at this hypothesis district (being that r is region, 0.82R tip whirlpool) substantially, fine grid blocks is established, and can carry out next step.Otherwise, when tip vortex cavitation numerical forecasting is not during substantially at this hypothesis district, according to site error, parameter r is adjusted, and recalculates, until result conforms to substantially.

Step 11, imports grid by fine grid blocks file and divides software.

Step 12, set up wake flow apart from the other two kinds of computational fields grids that are respectively 7D and 9D:

On the basis of fine grid blocks, set up other two kinds of computational fields, its wake flow distance is respectively 7D and 9D, and originally in fine grid blocks, the unit size on every limit is constant, divide the wake flow subregion grid increasing, its grid cell is of a size of former fine grid blocks in 5D place grid cell size;

Step 13, the cavitation model of repeating step three and turbulence model are set, the numerical evaluation setting parameter of step 4 and step 5 carry out numerical evaluation.

Step 14, the result that is 5D with wake flow distance compares definite Bestgrid:

Propeller center is extremely important to the accurate forecast of tip vortex cavitation to the pressure export distance of computational fields.Waste resource that distance is large, distance is little can cause that in numerical procedure, wake flow dissipation increases, thereby reduces the tip vortex strength degree in numerical forecasting, and tip vortex cavitation can not finely be forecast in numerical simulation.And this parameter can change along with the difference of screw propeller shape and duty parameter.Therefore, need in the methods of the invention this parameter to measure.Here set up respectively three kinds of computational fields, its wake flow distance is respectively 5D, 7D and 9D.Wherein the first computational fields (5D) completes in preceding step, after two kinds of computational fields only need increase the grid of wake flow subregion.These three kinds of computational fields grids are carried out analyzing after numerical evaluation.When along with wake flow is apart from increase, when tip vortex cavitation numerical forecasting precision no longer obviously increases, establish best wake flow distance and Bestgrid.

Step 15, utilizes Bestgrid and numerical evaluation parameter to carry out numerical forecasting to other operating mode.

The present invention has more efficiently forecast the tip vortex cavitation of E779A type screw propeller (Fig. 4 is shown in by model) under several different operating modes.Fig. 7 has shown under two kinds of working conditions, the experimental result of reporting in E779A screw propeller tip vortex cavitation numerical simulation result and other document.In figure, adopt the contour surface (black region in figure) of vapour phase volume fraction av=0.2 as cavitation numerical result.In Fig. 7 numerical result, tip vortex cavitation tail is very obvious, and coincide with experimental result.Illustrate that the inventive method is remarkable to tip vortex cavitation numerical forecasting effect.Fig. 8 is the tip vortex cavitation numerical simulation result in nearly 1/4 cycle under nonlinear inflow condition.Fig. 8 (a) (b) (c) (d) (e) (f) in black region be respectively the E779A oar mould being rotated counterclockwise No. 1 blade at 0 °, 18 °, 36 °, 54 °, the cavitation forecast in the position moment such as 72 ° and 81 °.In figure, show that tip vortex cavitation tail is very obvious, its length and thickness, along with blade position of rotation changes and changes, show that in the present invention, using method has remarkable result equally to the numerical simulation of tip vortex cavitation under nonlinear inflow condition.

Claims (10)

1. about a Numerical Prediction Method for underwater propeller tip vortex cavitation, it is characterized in that: comprise the following steps:

(a) basic grid is determined;

(b) fine grid blocks is determined;

(c) Bestgrid is determined;

(d) utilize Bestgrid to forecast the tip vortex cavitation under required working condition.

2. a kind of Numerical Prediction Method about underwater propeller tip vortex cavitation as claimed in claim 1, is characterized in that: described basic grid is determined, comprised following concrete steps:

Step 1, utilizes modeling software to set up screw propeller 3-D geometric model, and is imported grid division software;

Step 2, divide in software and set up three kinds of alternative grids at grid, and imported calculation procedure:

The computational fields of three kinds of alternative grids is identical, computational fields is cylindrical, its speed frontier distance propeller center that becomes a mandarin is 1.5D, D is airscrew diameter, downstream pressure outlet frontier distance is 5D, propeller center is 2.5D to side periphery distance, and the number of grid of three kinds of alternative grids increases gradually, is about respectively 2,000,000,3,000,000 and 4,000,000;

Step 3, cavitation model and turbulence model are set:

Adopt full cavitation model and Renormalization Group turbulence model (RNG k-ε turbulence model), and its important parameter is revised, in correction to transformation ratio parameter in cavitation model and turbulence model, the correction of turbulent viscosity coefficient adopts C language compilation, and recycling macro call (DEFINE_TURBULENT_VISCOSITY etc.) form embeds calculation procedure;

Step 4, numerical evaluation setting parameter:

Numerical parameter is set the correlation parameter setting that comprises working condition, boundary condition and numerical algorithm; Working condition is mainly set screw propeller rotational speed, and environmental pressure and inflow velocity value, determine screw propeller dimensionless group, i.e. advance coefficient (J) and cavitation number (σ _n); Set for boundary condition, speed inlet boundary adopts inflow velocity value, and far field boundary condition adopts inflow velocity value to set, and the top hole pressure at downstream pressure outlet interface is set to static pressure; For numerical algorithm, it is discrete that in Na Wei-Stokes (N-S) equation, convective term adopts Second-order Up-wind form, diffusion term adopts Using Second-Order Central Difference form discrete, velocity pressure coupling adopts the SIMPLE algorithm that is applicable to non-structured grid, use pointwise Gauss-Seidel iterative discrete equation, utilize the convergence of algebraic multigrid speed-up computation, adopt Sliding mesh computing technique for non-permanent calculating, in computation process, adopt parallel computing;

Mass conservation continuity (continuity) residual error convergence in equation is three rank, in equation, other physical quantity residual error convergence is quadravalence, the initial value that utilizes the convergence solution of single-phase flow to solve as multinomial stream, and the initial value calculating steady state solution as unstable state, in order to ensure calculating convergence, screw propeller rotational speed is progressively increased to predetermined value, and suitably dwindles relaxation factor;

Step 5, carry out numerical evaluation:

Adopt computation process step by step step by step, in screw propeller duty parameter, environmental pressure and inflow velocity can directly be set to operating mode value, and revolution speed of propeller adopts classification to increase, until be increased to predetermined operating mode value;

First calculate non-cavitating model Flow Field Distribution, by the time after calculation stability, open again cavitation model, first the parameters such as pressure, density, momentum and vapour phase mark being carried out to single order precision discrete scheme calculates, after calculation stability, again discrete precision is brought up to second order or QUCIK etc., for the stability that ensures that second order calculates, sub-relaxation factor is suitably adjusted to reduction;

Step 6, establish basic grid according to hydrodynamic force numerical result:

Utilize these three kinds of alternative grids to carry out numerical evaluation to Propeller, and result of calculation is analyzed; After tending to be steady after the result of calculation of hydrodynamic parameter (thrust coefficient and moment coefficient) increases along with grid number, think that numerical evaluation is along with grid increases and tends towards stability, after selected calculation stability the minimum grid of number of grid as basic grid.

3. a kind of Numerical Prediction Method about underwater propeller tip vortex cavitation as claimed in claim 1, is characterized in that: described fine grid blocks is determined, comprised following concrete steps:

Step 7, imports grid by the file of basic grid and divides software;

Step 8, on the basis of basic grid, according to tip whirlpool shape, utilize grid to divide software and set up tip whirlpool area grid:

First, in professional modeling tool, set up the conduit region geometric model of a helix shape conforming to model parameter according to equidistant helix mathematical model, this region is the tip vortex cavitation generation area of supposition, again helix geometry models is imported in basic grid, and carry out grid division, equidistant helix shape area grid is from closing on blade blade tip, equidistant helix radius is preset as 0.82R, the unit grid size in helix region is about 0.0001D, and in former basic grid, the unit size on every limit is constant;

Step 9, the cavitation model of repeating step three and turbulence model are set, the numerical evaluation setting parameter of step 4 and step 5 carry out numerical evaluation;

Step 10, establish fine grid blocks according to tip vortex cavitation numerical forecasting result:

This step is mainly the important parameter r establishing in isometric helix line model, the length of helix radius r is about (0.81-0.83) R, according to tip vortex cavitation numerical result adjust and finally determine this parameter, concrete grammar is, when tip vortex cavitation in the result of calculation of step 9 is arranged in this hypothesis district (being the helix pipeline region that r that step 8 sets sets up is exactly tip vortex cavitation generation area) substantially, fine grid blocks is established, and can carry out next step; Otherwise, when tip vortex cavitation numerical forecasting is not during substantially at this hypothesis district, according to site error, parameter r is adjusted, and return to step 8 and step 9 recalculates, until result conforms to substantially, establish fine grid blocks.

4. a kind of Numerical Prediction Method about underwater propeller tip vortex cavitation as claimed in claim 1, is characterized in that: described Bestgrid is determined, comprised following concrete steps:

Step 11, imports grid by fine grid blocks file and divides software;

Step 12, set up wake flow apart from the other two kinds of computational fields grids that are respectively 7D and 9D:

On the basis of fine grid blocks, set up other two kinds of computational fields, its wake flow distance is respectively 7D and 9D, and originally in fine grid blocks, the unit size on every limit is constant, divide the wake flow subregion grid increasing, its grid cell is of a size of former fine grid blocks in 5D place grid cell size;

Step 13, the cavitation model of repeating step three and turbulence model are set, the numerical evaluation setting parameter of step 4 and step 5 carry out numerical evaluation;

Step 14, the result that is 5D with wake flow distance compares definite Bestgrid:

These three kinds of computational fields grids are carried out analyzing after numerical evaluation, when along with wake flow is apart from increase, when tip vortex cavitation numerical forecasting precision no longer obviously increases, establish best wake flow distance and Bestgrid.

5. a kind of Numerical Prediction Method about underwater propeller tip vortex cavitation as described in claim 2,3 or 4, is characterized in that: described full cavitation model and parameter thereof are modified to:

Work as p<p _vtime, steam generation rate is:

$R_e=\mathrm{Ce}\frac{k}{\mathrm{\&gamma;}}{\mathrm{\&rho;}}_l{\mathrm{\&rho;}}_v\sqrt{\frac{2}{3}\frac{p_v-p}{{\mathrm{\&rho;}}_l}}(1-f_v)$

Work as p>p _vtime, vapour phase becomes liquid phase, steam solidification rate R _cfor:

$R_c=\mathrm{Cc}\frac{k}{\mathrm{\&gamma;}}{\mathrm{\&rho;}}_l{\mathrm{\&rho;}}_v\sqrt{\frac{2}{3}\frac{p-p_v}{{\mathrm{\&rho;}}_l}}{f}_v$

Wherein, f _v=α _vρ _v/ ρ _mfor vapour phase massfraction, vaporization coefficient Ce=0.02 and condensation coefficient Cc=0.01 are empirical parameter.

6. a kind of Numerical Prediction Method about underwater propeller tip vortex cavitation as described in claim 2,3 or 4, is characterized in that: described Renormalization Group turbulence model and parameter thereof are modified to:

K equation and the ε equation of Renormalization Group turbulence model are respectively:

$\frac{\&PartialD;}{\&PartialD;t}\left({\mathrm{\&rho;}}_{m}k\right)+\frac{\&PartialD;}{\&PartialD;x_j}\left({\mathrm{\&rho;}}_{m}k{u}_{\mathrm{mj}}\right)=\frac{\&PartialD;}{\&PartialD;x_j}\left[\left({\mathrm{\&alpha;}}_k\mathrm{\&mu;}\right)\frac{\&PartialD;k}{\&PartialD;x_j}\right]+G-{\mathrm{\&rho;}}_m\mathrm{\&epsiv;}$

$\frac{\&PartialD;}{\&PartialD;t}\left({\mathrm{\&rho;}}_m\mathrm{\&epsiv;}\right)+\frac{\&PartialD;}{\&PartialD;x_j}\left({\mathrm{\&rho;}}_m\mathrm{\&epsiv;}{u}_{\mathrm{mj}}\right)=\frac{\&PartialD;}{\&PartialD;x_j}\left[\left({\mathrm{\&alpha;}}_{\mathrm{\&epsiv;}}\mathrm{\&mu;}\right)\frac{\&PartialD;\mathrm{\&epsiv;}}{\&PartialD;x_j}\right]-R+C_{1\mathrm{\&epsiv;}}\frac{\mathrm{\&epsiv;}}{k}G-C_{2\mathrm{\&epsiv;}}{\mathrm{\&rho;}}_m\frac{{\mathrm{\&epsiv;}}^2}{k}$

In formula, Turbulent Kinetic dissipative shock wave (Turbulent Dissipation Rate) the a reciprocal of effective turbulent prandtl number of Turbulent Kinetic k and dissipative shock wave ε _k=a _ε=1.39; Model parameter C _{1 ε}=1.47, C _{2 ε}=1.68; Viscosity coefficient is μ=μ _t+ μ _m, μ _mfor mixed flow coefficient of viscosity; Turbulent viscosity coefficient μ _t=[ρ _v+ α _l ¹⁰(ρ _l-ρ _v)] C _μk ²/ ε, C _μ=0.085.

7. a kind of Numerical Prediction Method about underwater propeller tip vortex cavitation as claimed in claim 3, is characterized in that: described equidistant helix mathematical model and parameter thereof are as follows:

$\left\{\begin{array}{c}x=r\cos \left(\frac{2\mathrm{\&pi;}}{k}z\right)\\ y=r\sin \left(\frac{2\mathrm{\&pi;}}{k}z\right)\end{array}\right.$

Wherein, x, y, z is respectively three coordinate axis of cartesian coordinate system, and constant r is helix radius, and the length of r is about (0.81-0.83) R, and R is propeller radius; Constant k is the distance that often rotates a circle and advance on x axle of helix, k=U _∞/ n, U _∞for inflow velocity, n is revolution speed of propeller.

8. a kind of Numerical Prediction Method about underwater propeller tip vortex cavitation as claimed in claim 2, it is characterized in that: the method that described grid is divided is: adopt subregion mixed mesh method grid division: screw propeller around flow field regions adopts non-structured grid method to divide, grid is reduced to blade tip gradually by propeller hub, blade tip place surface grids is triangle, size is about 0.001D Jiang Grains place for 0.015D, and D is airscrew diameter; Adopt boundary layer grid at blade surface, adopt structured grid to divide the computational fields of the peripheral regular shape of screw propeller.

9. a kind of Numerical Prediction Method about underwater propeller tip vortex cavitation as claimed in claim 8, it is characterized in that: blade surface adopts boundary layer grid, boundary layer grid has 4 layers, its adjacent two layers aspect ratio is 1.1, ground floor grid cell height is about 0.0008D, makes dimensionless group 20<y ⁺<300.

10. a kind of Numerical Prediction Method about underwater propeller tip vortex cavitation as claimed in claim 2, it is characterized in that: in the numerical evaluation of step 5, need to adjust and reduce sub-relaxation factor, the isoparametric sub-relaxation factor of pressure, momentum, vapour phase mark, Turbulent Kinetic, turbulence dissipation rate and turbulent flow stickiness is set as respectively: 0.25,0.6,0.2,0.7,0.7,0.9.