A numerical simulation method for predicting and analyzing propeller cavitation

Abstract

This invention discloses a numerical simulation method for predicting and analyzing propeller cavitation, comprising the following steps: S1, establishing a geometric model including the propeller and the surrounding cylindrical flow domain; S2, employing a transient solver based on the finite volume method, selecting a turbulence model and a multiphase flow model; S3, performing non-cavitation steady-state calculations to obtain the initial flow field; S4, after cavitation is triggered, initiating transient calculations and enabling a dynamic adaptive mesh, monitoring and recording the change curve of the total cavitation volume over time; S5, comparing the simulated cavitation morphology with publicly available experimental data from a propeller water tunnel to verify the prediction accuracy of the cavitation morphology. This invention solves the problems of inaccurate cavitation initiation prediction, weak correlation of noise sources, and high computational cost in existing propeller noise simulation methods.

Description

Technical Field

[0001] This invention belongs to the field of fluid dynamics simulation technology, specifically, it relates to a numerical simulation method for predicting and analyzing propeller cavitation. Background Technology

[0002] The CFD (Computational Fluid Dynamics) simulation method for propeller noise mainly combines high-precision flow field solving with acoustic analogy theory or direct acoustic simulation to predict noise generated by mechanisms such as propeller rotation and cavitation.

[0003] Traditional methods rely on a single pressure threshold (such as saturated vapor pressure) to determine cavitation initiation, neglecting the influence of key factors such as the water flow nucleus and turbulent fluctuations. This leads to significant discrepancies between predicted cavitation location and conditions and experimental results. This is because traditional methods use a constant saturated vapor pressure (Pv) as the cavitation initiation criterion, ignoring the impact of turbulent pressure fluctuations. According to cavitation nucleus theory, cavitation initiation requires two conditions to be met simultaneously: a) sufficiently low local pressure; b) the low-pressure state must be maintained for a sufficiently long time for the cavitation nucleus to grow to an observable scale. Turbulence introduces high-frequency pressure fluctuations, which may cause instantaneous pressures below Pv, but the duration is extremely short, insufficient to trigger effective cavitation. Therefore, the effective cavitation initiation pressure threshold should be a "time-averaged" concept; it must be lower than Pv by a certain amount to ensure that, under turbulent fluctuations, the instantaneous pressure has a high probability of remaining below Pv for a sufficiently long time. This results in inaccurate cavitation initiation predictions. Secondly, using fixed meshes or simple static mesh refinement for transient simulations makes it difficult to simultaneously guarantee accuracy and efficiency throughout the entire dynamic evolution of cavitation bubbles (generation, growth, and collapse). This results in weak correlation of noise sources and high computational costs. Summary of the Invention

[0004] To address the problems of inaccurate cavitation initiation prediction, weak correlation of noise sources, and high computational cost in existing propeller noise simulation methods, this invention provides a numerical simulation method for predicting and analyzing propeller cavitation.

[0005] To achieve the above-mentioned technical objectives, the technical solution adopted by the present invention is as follows: A numerical simulation method for predicting and analyzing propeller cavitation includes the following steps: S1. Establish a geometric model that includes the propeller and the outer cylindrical flow domain; S2. A transient solver based on the finite volume method is adopted, and a turbulence model and a multiphase flow model are selected. S3. Perform non-cavitation steady-state calculations to obtain the initial flow field; S4. After cavitation is triggered, start transient calculation and enable dynamic adaptive mesh to monitor and record the change curve of total cavitation volume over time. S5. Compare the simulated cavitation morphology (such as vortex cavitation length) with the publicly available experimental data of the propeller water tunnel to verify the prediction accuracy of the cavitation morphology.

[0006] Furthermore, in step S2, the turbulence model adopts the SST k-ω model, and the multiphase flow model adopts the Schnerr-Sauer cavitation model.

[0007] Furthermore, when starting transient calculation in step S4, the time step is set to the time required for the blade to rotate 0.5 degrees; at the same time, dynamic adaptive mesh technology is enabled to set the encryption criteria.

[0008] Furthermore, the detailed steps of step S3 are as follows: S301. After the flow field stabilizes, not only is the pressure field monitored, but the high vortex core region in the vorticity field (such as the starting position of the blade tip vortex and the hub vortex) is calculated and identified in real time. S302. Modify the triggering condition for cavitation phase transition from local pressure below saturated vapor pressure to local pressure below a dynamic threshold P. th Furthermore, this region also exhibits high vorticity intensity. The dynamic threshold P... th Adjustments are made based on the local turbulence intensity.

[0009] Furthermore, the functional relationship for non-cavitation steady-state calculation in step S3 is as follows: Among them, P v The saturated vapor pressure (constant) at the operating water temperature is T, where ρ is the fluid density and T is the saturated vapor pressure at the operating water temperature. u The local turbulence intensity, T, is provided in real time in each computational cell by a transient solver based on the finite volume method. u =u rms ′ / U local , where u rms ′ represents the root mean square of the pulsation velocity, U ref The characteristic velocity is usually taken as the free flow velocity or the local time-averaged velocity modulus |U local |, α is the core empirical coefficient (dimensionless).

[0010] Furthermore, the core empirical coefficient α is a calibration parameter related to propeller type and operating conditions, and its calibration simulation method is as follows: Select a benchmark case: Select a benchmark model with high-quality cavitation initiation experimental data (such as NACA0015 hydrofoil or PPTC11 propeller at a specific advance coefficient).

[0011] Perform non-cavitation LES simulation: Perform large eddy simulation (LES) on the baseline model to obtain a high-precision turbulent fluctuation field.

[0012] Data Comparison and Inference: Compare the cavitation initiation location observed in the experiment (e.g., at a specific chord length on the blade suction surface). Locate this location in the CFD flow field and read the time-averaged Tu value at that location.

[0013] Inverse calculation of the core empirical coefficient α: Assume that P is satisfied at the initial moment of birth. local =P th The average pressure P at this location is known. local Substitute into formula P local =P v −0.5ρ(α⋅T u ⋅U ref )2, and the value of α under this working condition can be solved in reverse.

[0014] Establish a database: By changing the angle of attack or advance coefficient, multiple core empirical coefficient α values ​​are obtained, which can be fitted into a function related to the operating parameters.

[0015] Furthermore, in step S4, mesh refinement not only depends on the cavitation volume fraction gradient but also correlates with cavitation volume acceleration. During periods and regions of rapid cavitation change (rapid generation or collapse), the mesh is automatically refined to the maximum extent; in regions where cavitation is relatively stable, the mesh is moderately coarsened. This method ensures that computational resources are always focused on the regions where the most significant noise sources are generated.

[0016] Furthermore, the physical details of noise generated during the dynamic process of cavitation are analyzed through mesh encryption, including the cavitation volume acceleration α. v =∂ 2 V / ∂t 2 It is directly proportional to the sound pressure level of the radiated noise; α v Larger time periods and regions are the "burst points" of high-intensity noise, requiring the finest mesh to capture the dramatic deformation and collapse details of the cavitation interface; conversely, coarser meshes can be used to save resources.

[0017] An exponential function is used to describe the target mesh size Δ. target With α v Relationship: Where, Δ min This is the minimum permissible mesh size. It is determined by the minimum cavitation characteristic scale that needs to be analyzed (e.g., the minimum diameter of the cavitation collapse vortex). When cavitation activity is extremely intense (|α) v As |→∞), the mesh size approaches this value. Δ base This represents the background mesh size. It is the basic size of the unrefined region, corresponding to the region in the flow field where no intense cavitation activity occurs (i.e., α). v The mesh size is approximately 0. It is the upper limit of the coarsening of the entire adaptive mesh. α refβ: Scale normalization parameter and shape parameter.

[0018] Calibration method: Calculate the global average α over a complete cavitation shedding cycle. v (t) Time series.

[0019] Analyze its statistical characteristics: determine its maximum value α v ,max, mean, and standard deviation σ.

[0020] Parameter settings: β is a coefficient controlling the decay rate (e.g., β=2). At that time, the mesh size is approximately Δ. base −0.86(Δ base −Δ min This means significant encryption.

[0021] α ref Reference value for cavitation volume acceleration (normalization factor) Definition: A physically meaningful scale used to assess the intensity of current cavitation activity αv.

[0022] Determination method: Typically, a statistical characteristic value of the absolute value of global cavitation volume acceleration is taken from a preliminary simulation or the previous time step. For example:

[0023] Average value: α ref =mean(∣α v |). This is the most common and stable choice, causing the function response to hover around the average activity level.

[0024] Function: To convert α v Divide by α ref Normalization was achieved, making the function exp (−β⋅∣α) v ∣ / α ref It exhibits better robustness and versatility under different operating conditions. When |α v |=α ref When the exponent term is e −β .

[0025] β: Attenuation coefficient (core control parameter) Definition: A dimensionless positive number that controls the mesh size from Δ base Attenuation to Δ min The "speed" or "sensitivity".

[0026] Physical and numerical significance: The larger the β value, the faster the exponential function decays. This means that even if α... v Slightly higher than α ref The grid will also be rapidly encrypted to near Δ minThe calculations are more aggressive and refined, tending to capture all possible high-frequency noise sources, but the computational cost is also higher.

[0027] A smaller β value indicates a smoother function decay. α is needed. v Much greater than α ref This will trigger strong encryption. The calculation is more conservative and economical, but it may miss some fleeting secondary noise events.

[0028] Empirical value range: typically β∈[1,5]. β=1 represents a mild response, β=3 a typical moderate response, and β=5 a radical response. This function ensures that in high-noise source regions where cavitation collapses violently, the mesh is automatically refined to a preset limit of fineness Δ. min In regions where cavitation activity is relatively flat, the background grid size Δ is maintained. base This enables adaptive optimal allocation of computing resources in the spatiotemporal domain.

[0029] Compared with the prior art, the present invention has the following advantages: High prediction accuracy: By using the "pressure-vortex core" collaborative criterion, the physical mechanism of cavitation initiation is captured more accurately, making the simulated initiation location and conditions closer to reality.

[0030] Computational efficiency optimization: The dynamic adaptive mesh strategy precisely allocates computational resources to noise source regions. Compared with a uniformly dense mesh across the entire domain, it is expected to save 30%-50% of computation time while ensuring accuracy.

[0031] The noise correlation is direct and reliable: the noise source model based on cavitation volume acceleration has a clear physical meaning, providing the most direct and clear optimization guidance indicators for propeller noise reduction design (such as modifying the blade profile to suppress violent cavitation collapse). Attached Figure Description

[0032] Figure 1 This is a flowchart of a numerical simulation method for predicting and analyzing propeller cavitation in an embodiment of the present invention; Figure 2 This is a schematic diagram showing the cavitation display of the simulation results in an embodiment of the present invention.

[0033] Explanation of markings in the diagram: Detailed Implementation

[0034] To facilitate understanding by those skilled in the art, the present invention will be further described below with reference to embodiments and accompanying drawings. The content mentioned in the embodiments is not intended to limit the present invention.

[0035] like Figure 1 and 2As shown, this embodiment provides a numerical simulation method for predicting and analyzing propeller cavitation, including the following steps: S1. Establish a geometric model including the propeller and the surrounding cylindrical flow domain; the background mesh is tetrahedral, and multiple prismatic meshes are generated on the blade surface to resolve the boundary layer. Set the initial mesh size to approximately 5 million.

[0036] S2. A transient solver based on the finite volume method is used, and a turbulence model and a multiphase flow model are selected. The SST k-ω model is chosen for the turbulence model because it can better predict boundary layer separation and moderate adverse pressure gradients. The Schnerr-Sauer cavitation model is chosen for the multiphase flow model.

[0037] S3. Perform non-cavitation steady-state calculations to obtain the initial flow field; S4. After cavitation is triggered, start transient calculation and enable dynamic adaptive mesh to monitor and record the change curve of total cavitation volume over time. S5. Compare the simulated cavitation morphology (such as vortex cavitation length) with the publicly available experimental data of the propeller water tunnel to verify the prediction accuracy of the cavitation morphology.

[0038] In step S2, the turbulence model adopts the SST k-ω model, and the multiphase flow model adopts the Schnerr-Sauer cavitation model.

[0039] When starting transient calculation in step S4, the time step is set to the time required for the blade to rotate 0.5 degrees; at the same time, dynamic adaptive mesh technology is enabled to set the encryption criteria.

[0040] Detailed steps of step S3: S301. After the flow field stabilizes, not only is the pressure field monitored, but the high vortex core region in the vorticity field (such as the starting position of the blade tip vortex and the hub vortex) is calculated and identified in real time. S302. Modify the triggering condition for cavitation phase transition from local pressure below saturated vapor pressure to local pressure below a dynamic threshold P. th Furthermore, this region also exhibits high vorticity intensity. The dynamic threshold P... th Adjustments are made based on the local turbulence intensity.

[0041] The functional relationship for non-cavitation steady-state calculation in step S3 is as follows: Among them, P v The saturated vapor pressure (constant) at the operating water temperature is T, where ρ is the fluid density and T is the saturated vapor pressure at the operating water temperature. u The local turbulence intensity, T, is provided in real time in each computational cell by a transient solver based on the finite volume method. u =u rms ′ / U local, where u rms ′ represents the root mean square of the pulsation velocity, U ref The characteristic velocity is usually taken as the free flow velocity or the local time-averaged velocity modulus |U local |, α is the core empirical coefficient (dimensionless).

[0042] The core empirical coefficient α is a calibration parameter related to propeller type and operating conditions, and its calibration simulation method is as follows: Select a benchmark case: Select a benchmark model with high-quality cavitation initiation experimental data (such as NACA0015 hydrofoil or PPTC11 propeller at a specific advance coefficient).

[0043] Perform non-cavitation LES simulation: Perform large eddy simulation (LES) on the baseline model to obtain a high-precision turbulent fluctuation field.

[0044] Data Comparison and Inference: Compare the cavitation initiation location observed in the experiment (e.g., at a specific chord length on the blade suction surface). Locate this location in the CFD flow field and read the time-averaged Tu value at that location.

[0045] Inverse calculation of the core empirical coefficient α: Assume that P is satisfied at the initial moment of birth. local =P th The average pressure P at this location is known. local Substitute into formula P local =P v −0.5ρ(α⋅T u ⋅U ref )2, and the value of α under this working condition can be solved in reverse.

[0046] Establish a database: By changing the angle of attack or advance coefficient, multiple core empirical coefficient α values ​​are obtained, which can be fitted into a function related to the operating parameters.

[0047] In step S4, mesh refinement depends not only on the cavitation volume fraction gradient but also on the cavitation volume acceleration. During periods and regions of rapid cavitation change (rapid generation or collapse), the mesh is automatically refined to the maximum extent; in regions where cavitation is relatively stable, the mesh is moderately coarsened. This method ensures that computational resources are always focused on the regions where the most significant noise sources are generated.

[0048] Physical details of noise generated during the dynamic process of cavitation bubble analysis in mesh encryption, including cavitation bubble volume acceleration α. v =∂ 2 V / ∂t 2 It is directly proportional to the sound pressure level of the radiated noise; α v Larger time periods and regions are the "burst points" of high-intensity noise, requiring the finest mesh to capture the dramatic deformation and collapse details of the cavitation interface; conversely, coarser meshes can be used to save resources.

[0049] An exponential function is used to describe the target mesh size Δ. target With α v Relationship: Where, Δ min This is the minimum permissible mesh size. It is determined by the minimum cavitation characteristic scale that needs to be analyzed (e.g., the minimum diameter of the cavitation collapse vortex). When cavitation activity is extremely intense (|α) v As |→∞), the mesh size approaches this value. Δ base This represents the background mesh size. It is the basic size of the unrefined region, corresponding to the region in the flow field where no intense cavitation activity occurs (i.e., α). v The mesh size is approximately 0. It is the upper limit of the coarsening of the entire adaptive mesh. α ref β: Scale normalization parameter and shape parameter.

[0050] Calibration method: Calculate the global average α over a complete cavitation shedding cycle. v (t) Time series.

[0051] Analyze its statistical characteristics: determine its maximum value α v ,max, mean, and standard deviation σ.

[0052] Parameter settings: β is a coefficient controlling the decay rate (e.g., β=2). At that time, the mesh size is approximately Δ. base −0.86(Δ base −Δ min This means significant encryption.

[0053] α ref Reference value for cavitation volume acceleration (normalization factor) Definition: A physically meaningful scale used to assess the intensity of current cavitation activity αv.

[0054] Determination method: Typically, a statistical characteristic value of the absolute value of global cavitation volume acceleration is taken from a preliminary simulation or the previous time step. For example:

[0055] Average value: α ref =mean(∣α v |). This is the most common and stable choice, causing the function response to hover around the average activity level.

[0056] Function: To convert α v Divide by α ref Normalization was achieved, making the function exp (−β⋅∣α) v ∣ / α ref It exhibits better robustness and versatility under different operating conditions. When |αv |=α ref When the exponent term is e −β .

[0057] β: Attenuation coefficient (core control parameter) Definition: A dimensionless positive number that controls the mesh size from Δ base Attenuation to Δ min The "speed" or "sensitivity".

[0058] Physical and numerical significance: The larger the β value, the faster the exponential function decays. This means that even if α... v Slightly higher than α ref The grid will also be rapidly encrypted to near Δ min The calculations are more aggressive and refined, tending to capture all possible high-frequency noise sources, but the computational cost is also higher.

[0059] A smaller β value indicates a smoother function decay. α is needed. v Much greater than α ref This will trigger strong encryption. The calculation is more conservative and economical, but it may miss some fleeting secondary noise events.

[0060] Empirical value range: typically β∈[1,5]. β=1 represents a mild response, β=3 a typical moderate response, and β=5 a radical response. This function ensures that in high-noise source regions where cavitation collapses violently, the mesh is automatically refined to a preset limit of fineness Δ. min In regions where cavitation activity is relatively flat, the background grid size Δ is maintained. base This enables adaptive optimal allocation of computing resources in the spatiotemporal domain.

[0061] Compared with the prior art, the present invention has the following advantages: High prediction accuracy: By using the "pressure-vortex core" collaborative criterion, the physical mechanism of cavitation initiation is captured more accurately, making the simulated initiation location and conditions closer to reality.

[0062] Computational efficiency optimization: The dynamic adaptive mesh strategy precisely allocates computational resources to noise source regions. Compared with a uniformly dense mesh across the entire domain, it is expected to save 30%-50% of computation time while ensuring accuracy.

[0063] The noise correlation is direct and reliable: the noise source model based on cavitation volume acceleration has a clear physical meaning, providing the most direct and clear optimization guidance indicators for propeller noise reduction design (such as modifying the blade profile to suppress violent cavitation collapse).

[0064] The numerical simulation method for predicting and analyzing propeller cavitation provided in this application has been described in detail above. The specific embodiments are described to help understand the method and its core ideas. It should be noted that those skilled in the art can make various improvements and modifications to this application without departing from its principles, and these improvements and modifications also fall within the protection scope of the claims of this application.

Claims

1. A numerical simulation method for predicting and analyzing propeller cavitation, characterized in that, Including the following steps: S1. Establish a geometric model that includes the propeller and the outer cylindrical flow domain; S2. A transient solver based on the finite volume method is adopted, and a turbulence model and a multiphase flow model are selected. S3. Perform non-cavitation steady-state calculations to obtain the initial flow field; S4. After cavitation is triggered, start transient calculation and enable dynamic adaptive mesh to monitor and record the change curve of total cavitation volume over time. S5. Compare the simulated cavitation morphology with publicly available propeller-driven water tunnel experimental data to verify the accuracy of the cavitation morphology prediction.

2. The numerical simulation method for predicting and analyzing propeller cavitation according to claim 1, characterized in that, In step S2, the turbulence model adopts the SST k-ω model, and the multiphase flow model adopts the Schnerr-Sauer cavitation model.

3. The numerical simulation method for predicting and analyzing propeller cavitation according to claim 2, characterized in that, When starting transient calculation in step S4, the time step is set to the time required for the blade to rotate 0.5 degrees; at the same time, dynamic adaptive mesh technology is enabled to set the encryption criteria.

4. The numerical simulation method for predicting and analyzing propeller cavitation according to claim 3, characterized in that, Detailed steps of step S3: S301. After the flow field stabilizes, not only is the pressure field monitored, but the high vortex core region in the vorticity field is also calculated and identified in real time. S302. Modify the triggering condition for cavitation phase change from local pressure below saturated vapor pressure to local pressure below a dynamic threshold.

5. A numerical simulation method for predicting and analyzing propeller cavitation according to claim 4, characterized in that, The functional relationship for non-cavitation steady-state calculation in step S3 is as follows: Among them, P v The saturated vapor pressure is T at the operating water temperature, ρ is the fluid density, and T is the saturated vapor pressure at the operating water temperature. u U represents the local turbulence intensity. ref The characteristic velocity is taken as the free flow velocity or the local time-averaged velocity modulus, and α is the core empirical coefficient.

6. A numerical simulation method for predicting and analyzing propeller cavitation according to claim 5, characterized in that, The core empirical coefficient α is a calibration parameter related to propeller type and operating conditions, and its calibration simulation method is as follows: Select a baseline model: Select a baseline model with high-quality cavitation initiation experimental data; Perform non-cavitation LES simulation: Large eddy simulation is performed on the benchmark model to obtain a high-precision turbulent fluctuation field; Data comparison and reverse inference: Compare the cavitation initiation location observed in the experiment, find the location in the transient flow field based on the finite volume method, and read the time-averaged local turbulence intensity value at the location; Inverse calculation of the core empirical coefficient α: Assume that P is satisfied at the initial moment of birth. local =P th The average pressure P at a known location local Substitute into formula P local =P v −0.5ρ(α⋅T u ⋅U ref )2, solve for the value of α under this working condition; Establish a database: By changing the angle of attack or advance coefficient, multiple core empirical coefficient α values ​​are obtained, which can be fitted into a function related to the operating parameters.

7. A numerical simulation method for predicting and analyzing propeller cavitation according to claim 6, characterized in that, In step S4, mesh refinement not only depends on the cavitation volume fraction gradient, but also on the cavitation volume acceleration; in the stage and region of drastic cavitation changes, the mesh is automatically refined to the maximum extent; in the region of relatively stable cavitation, the mesh is moderately coarsened.

8. A numerical simulation method for predicting and analyzing propeller cavitation according to claim 7, characterized in that, Physical details of noise generated during the dynamic process of cavitation bubble analysis in mesh encryption, including cavitation bubble volume acceleration α. v =∂ 2 V / ∂t 2 It is directly proportional to the sound pressure level of the radiated noise; an exponential function is used to describe the target mesh size Δ. target With α v Relationship: Where, Δ min For the minimum allowable mesh size, Δ base α is the background grid size. ref β represents the scale normalization parameter and the shape parameter.

Images