SPH particle distribution optimization method based on gradient particle migration
Through the gradient particle migration optimization method, the problems of particle stratification and low computational efficiency in the WCSPH method are solved, and the simulation accuracy and efficiency of engineering problems such as ship water impact are improved.
Patent Information
- Application Number
- CN202510967851.1
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2025-07-14
- Publication Date
- 2025-09-26
AI Technical Summary
In the weakly compressible smooth particle hydrodynamics (WCSPH) method, particle stratification, aggregation, or fracture causes pressure oscillations and distortions in the calculation of motion resistance. Existing particle migration technology has low computational efficiency and cannot effectively handle boundary particles, making it difficult to meet the requirements of refined simulation of large-scale fluid-structure interaction problems.
A localized migration strategy based on gradient particle migration and a migration gradient correction mechanism are adopted. By calculating the Euclidean distance from the particle to the migration area boundary and the correction coefficient, the particle distribution is optimized. Combined with the density field update, the calculation accuracy and efficiency are improved.
The pressure calculation accuracy in the structural mutation area is significantly improved, the calculation time is reduced, and the applicability of the WCSPH method in engineering problems such as ship water impact is enhanced.
Smart Images
Figure CN120706327A_ABST
Abstract
Description
Technical Field
[0001] The present invention belongs to the field of particle distribution optimization, and in particular relates to an SPH particle distribution optimization method based on gradient particle migration. Background Art
[0002] In the weakly compressible smooth particle hydrodynamics (WCSPH) method, the fluid computational domain is discretized into particles with physical properties for simulation. However, its purely Lagrangian characteristics lead to non-physical particle stratification, aggregation, or fracture near the free liquid surface and in areas of structural mutations, resulting in pressure oscillations and distortion in the calculation of motion resistance, seriously affecting the simulation accuracy of engineering problems such as ship entry impact and real ship surge effects.
[0003] Although the existing particle transfer technology (PST) based on the improved scheme of Fick's law can optimize particle distribution, it still has significant limitations: on the one hand, traditional PST cannot effectively handle boundary particles whose support domain is truncated, and it needs to globally solve the eigenvalues and accurately identify the free surface particles, resulting in low computational efficiency; on the other hand, the migration process only relies on the particle density distribution within the region, ignoring the influence of particles outside the region, which is prone to interface fault problems. Although subsequent studies have proposed the δ+-SPH method to introduce corrective aggregation forces and the VEM / VCS scheme to improve continuity, the global migration strategy still causes the computational burden to increase sharply with the particle size. When the number of particles reaches the 100,000 level, the computational resource consumption of the existing PST technology far exceeds the feasible range of engineering practice, making it difficult to meet the requirements of refined simulation such as the pressure peak prediction at the corner of the ship, which seriously restricts the application of the WCSPH method in large-scale fluid-structure interaction problems.
[0004] In response to the above technical problems, the present invention proposes an SPH particle distribution optimization method based on gradient particle migration. Summary of the Invention
[0005] To address the aforementioned technical issues of pressure oscillation and computational efficiency caused by uneven particle distribution in the WCSPH method, a new SPH particle distribution optimization method based on gradient particle migration is proposed. This method aims to significantly improve the pressure calculation accuracy in structural abrupt regions (such as hull corners and sharp angles) through a localized migration strategy and a migration gradient correction mechanism, reducing computational time and enhancing the applicability of the WCSPH method to engineering problems such as ship impact.
[0006] The technical means adopted in the present invention are as follows: A SPH particle distribution optimization method based on gradient particle migration includes the following steps: A weakly compressible smooth particle fluid dynamics calculation step is performed, the fluid domain and rigid body motion boundaries are initialized in the dynamics solver, the center coordinates are defined according to the structural characteristics of the solution object, the migration area is delineated based on the center coordinates, and the particles in the marked area are traversed after the calculation step is completed, and the center coordinates are updated in real time when the solution object moves; the Euclidean distance from the marked particle to the migration area boundary is calculated, and the migration amount correction coefficient is calculated, and a synthetic migration vector is calculated in combination with the position vector of the marked particle; the position of the marked particle is updated according to the synthetic migration vector to obtain the updated particle position, and the density field is updated based on the updated particle position to complete the particle migration.
[0007] Furthermore, the calculation formula of the Euclidean distance from the marked particle to the migration area boundary is:
[0008] Among them, 𝑟𝑟 is the Euclidean distance from the marked particle to the boundary of the migration area, is the coordinate of the particle in the x direction, is the x-direction coordinate of the center coordinate, is the migration area range parameter, is the coordinate of the particle in the y direction, is the y-direction coordinate of the center coordinate, To mark the particle number, the calculation formula of the migration correction coefficient is:
[0009] in, is the migration correction factor, is the absolute value function.
[0010] Furthermore, the calculation formula of the synthetic migration vector is:
[0011] in, is the synthetic migration vector, is the migration correction factor, is the position vector of the marker particle, and the calculation formula of the position vector of the marker particle is:
[0012] in, is the scalar coefficient, is the identity matrix, is the normal vector of the free surface of the i-th particle, are the first correction displacement vector, the second correction displacement vector and the third correction displacement vector respectively, is the fluid domain, is the particle number, For the transition domain, is the interface domain, is the regional characteristic value, and the calculation formula of the scalar coefficient is:
[0013] in, is the neighboring particle of the i-th marked particle, is the normal vector of the free surface of the i-th particle, is the normal vector of the free surface of the jth particle, is the virtual particle domain.
[0014] Furthermore, the calculation formula for the position of the updated particle is:
[0015] in, is the updated particle position, is the particle's position vector, is the migration correction factor.
[0016] Furthermore, the migration area boundary meets the following conditions:
[0017] in, is the x-direction coordinate of the center coordinate, is the migration area range parameter, is the y-direction coordinate of the center coordinate, is the x-coordinate of any point in the migration area, is the y-coordinate of any point in the migration area.
[0018] Furthermore, the updating of the density field based on the updated particle positions to complete particle migration includes: Based on the updated particle positions, the density field is updated, and the particles are migrated to the updated particle positions to optimize the calculation. During the migration execution, the virtual particles are coupled to the fixed boundary, the sponge layer absorbs the shock wave reflection, and the characteristic line interpolation non-reflection boundary.
[0019] Compared with the prior art, the present invention has the following advantages: This method utilizes a localized migration strategy to avoid global eigenvalue calculations. Combined with a boundary distance gradient correction mechanism, it improves computational efficiency while maintaining peak pressure accuracy at sharp corners. Its dynamic migration region delineation mechanism adapts to complex scenarios such as ship corners and sudden changes in cabin structure, significantly enhancing the feasibility of simulating engineering problems such as the inrush effect of water inflow on real ships, demonstrating its outstanding industrial application value.
[0020] Based on the above reasons, the present invention can be widely promoted in fields such as particle distribution optimization. BRIEF DESCRIPTION OF THE DRAWINGS
[0021] In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the following briefly introduces the drawings required for use in the embodiments or the description of the prior art. Obviously, the drawings described below are some embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying any creative labor.
[0022] Figure 1 Schematic diagram of the SPH particle distribution optimization method based on gradient particle migration of the present invention.
[0023] Figure 2 Schematic diagram of the process of the SPH particle distribution optimization method based on gradient particle migration of the present invention.
[0024] Figure 3 This is the wedge-shaped body water entry model in the embodiment.
[0025] Figure 4 Schematic diagram of the pressure measurement points of the wedge-shaped body entering water in the embodiment.
[0026] Figure 5 1 is the pressure curve of each measuring point of the wedge-shaped body water entry pressure in the embodiment. DETAILED DESCRIPTION
[0027] In order to enable those skilled in the art to better understand the solutions of the present invention, the technical solutions in the embodiments of the present invention will be clearly and completely described below in conjunction with the drawings in the embodiments of the present invention. Obviously, the embodiments described are only part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without making creative efforts should fall within the scope of protection of the present invention.
[0028] It should be noted that the terms "first", "second", etc. in the description and claims of the present invention and the above-mentioned drawings are used to distinguish similar objects and are not necessarily used to describe a specific order or sequence. It should be understood that the numbers used in this way can be interchanged where appropriate, so that the embodiments of the present invention described herein can be implemented in an order other than those illustrated or described herein. In addition, the terms "including" and "having" and any variations thereof are intended to cover non-exclusive inclusions. For example, a process, method, system, product or device that includes a series of steps or units is not necessarily limited to those steps or units clearly listed, but may include other steps or units that are not clearly listed or inherent to these processes, methods, products or devices.
[0029] like Figure 1-5 As shown, the present invention provides an SPH particle distribution optimization method based on gradient particle migration, and the specific steps are as follows: S1. Perform a weakly compressible smooth particle hydrodynamics (WCSPH) calculation step. Initialize the fluid domain and rigid body motion boundaries in the dynamics solver. Define the center coordinates based on the structural characteristics of the solution object. Delineate the migration region based on the center coordinates. After the calculation step, traverse the particles in the marked region and update the center coordinates in real time as the solution object moves.
[0030] Among them, WCSPH is a common improvement method of the SPH method. The full name of SPH is Smoothed Particle Hydrodynamics, that is, smooth particle hydrodynamics, and the full name of WCSPH is Weakly Compressible Smoothed Particle Hydrodynamics, that is, weakly compressible smooth particle hydrodynamics, which is a relatively common method.
[0031] At the end of each computational step, a localized migration region is dynamically defined and its range parameters are determined based on particle distribution density and structural characteristics (including free surface locations and sharp corners or bends). Only particles within the region are marked, excluding global particles from participating in the migration. This region covers areas with significant particle delamination or fracture, adaptively adjusting the region's range based on particle spacing, while also marking particles within the region to preserve boundary topology.
[0032] S2. Calculate the Euclidean distance from the marker particle to the boundary of the migration area, and calculate the migration amount correction coefficient. Combined with the position vector of the marker particle, a synthetic migration vector is calculated.
[0033] Specifically, the migration amount correction is performed on the marked particles, and the distance from the target particle to the boundary of the migration area is first calculated.
[0034] The calculation formula of the Euclidean distance from the marked particle to the boundary of the migration area is:
[0035] Among them, 𝑟𝑟 is the Euclidean distance from the marked particle to the boundary of the migration area, is the coordinate of the particle in the x direction, is the x-direction coordinate of the center coordinate, is the migration area range parameter, is the coordinate of the particle in the y direction, is the y-direction coordinate of the center coordinate, is the particle number.
[0036] Then determine the correction factor through the formula. The calculation formula for the migration correction factor is:
[0037] in, is the migration correction factor, It is an absolute value function equivalent to the absolute value of rr-drpst, not a coefficient.
[0038] According to the particle position vector, the synthetic migration vector is calculated to make the boundary particle migration amount present a linear attenuation distribution to avoid interface fault. The calculation formula of the synthetic migration vector is:
[0039] in, is the synthetic migration vector, is the migration correction factor, is the position vector of the marker particle.
[0040] The position vector of the marker particle is calculated by Fick's law to partition the migration into fluid domain, transition domain and interface domain, and a distance weighted correction mechanism is introduced to eliminate interface discontinuities.
[0041] The calculation formula of the position vector of the marker particle is:
[0042] in, is the scalar coefficient, is the identity matrix, is the normal vector of the free surface of the i-th particle, are the first correction displacement vector, the second correction displacement vector and the third correction displacement vector respectively, is the fluid domain, is the particle number, For the transition domain, is the interface domain, is the regional characteristic value.
[0043] The scalar coefficient is calculated as:
[0044] in, is the neighboring particle of the i-th marked particle, is the normal vector of the free surface of the i-th particle, is the normal vector of the free surface of the jth particle, is the virtual particle domain.
[0045] S3. Update the position of the marked particle according to the synthesized migration vector, obtain the updated position of the particle, and update the density field to complete the particle migration.
[0046] Specifically, the density field is updated by combining the WCSPH density field data, and particle migration is performed according to the corrected synthetic migration vector to optimize the calculation, improve the calculation accuracy, and ensure the accuracy of the simulation.
[0047] During the migration process, the virtual particles are coupled to the fixed boundary, the sponge layer absorbs the shock wave reflection, and the characteristic line interpolation non-reflecting boundary ensures the free surface condition and pressure stability.
[0048] The calculation formula for the updated particle position is:
[0049] in, is the updated particle position, is the particle's position vector, is the migration correction factor.
[0050] Example S1. First, the migration area is dynamically delineated and the particle is marked. The migration area is delineated, such as Figure 3 As shown, according to the wedge structure characteristics, the wedge bottom rise angle θ = 25°, the top side length w = 1.2m, and the particle spacing d x =0.004m, mass M is 94kg, initial velocity v0 is 5.05m / s, artificial sound speed c0 is 50m / s, and the fluid domain and rigid body motion boundary are initialized in the WCSPH solver. c ,y c ) as the center point, define a rectangular migration area with a range parameter drpst = 0.02m. The migration area boundary satisfies: At the end of each WCSPH calculation step, the fluid particles are traversed and the particles in the region are marked. Then a dynamic adaptive mechanism is established. When the wedge movement causes the position of the tip to change, the center coordinates of the migration region (x c ,y c ), ensuring that the region covers the structural mutation area.
[0051] S2. Based on S1, calculate the migration gradient correction. First, calculate the boundary distance. For each marked particle i, its position vector r i =(x i ,y i ), calculate its Euclidean distance to the nearest migration area boundary, such as the particle position is (x c +0.015,y c), rr = |0.015 - 0.02| = 0.005 m. The correction factor is then determined to obtain the migration correction factor dr. For example, substituting rr = 0.005 m, the Euclidean distance to the nearest migration zone boundary in the previous example, yields dr = 0.75. The actual migration amount is then calculated using the particle position vector, which is calculated differently based on particle type to ensure particle continuity.
[0052] S3. Apply S1 and S2 to the wedge-shaped body entering the water, perform migration and verify its effect. First, update the particle position and synchronously correct the WCSPH density field. Select four pressure measurement points on the impact surface of the wedge, the specific locations are as follows: Figure 4 As shown, Figure 5 Pressure curves at the measuring points under different conditions are given. By comparison, the impact pressure values at the four measuring points near the bottom lift angle reach over 100,000 Pa, indicating that the local pressure in this area is very high. The pressures of the traditional weakly compressible SPH method and particle migration technology deviate from the experimental results, and the accuracy of the peak pressure calculation results is not reliable enough. Including measuring points 5 times the particle spacing away from the sharp corner, the GPST method can control the relative error within 15%, and the calculation accuracy is significantly improved compared to the traditional SPH.
[0053] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention, rather than to limit it. Although the present invention has been described in detail with reference to the above embodiments, those skilled in the art should understand that they can still modify the technical solutions described in the above embodiments, or replace some or all of the technical features therein with equivalents. However, these modifications or replacements do not cause the essence of the corresponding technical solutions to deviate from the scope of the technical solutions of the embodiments of the present invention.
Claims
1. A SPH particle distribution optimization method based on gradient particle migration, characterized in that: The following steps are involved: Performing a weakly compressible smooth particle fluid dynamics calculation step, initializing a fluid domain and a rigid body motion boundary in the dynamics solver, defining a center coordinate based on structural characteristics of the solution object, demarcating a migration region based on the center coordinate, traversing particles within the marked region after the calculation step, and updating the center coordinate in real time as the solution object moves; Calculate the Euclidean distance from the marker particle to the boundary of the migration area, and calculate the migration correction coefficient. Combined with the position vector of the marker particle, the synthetic migration vector is calculated. According to the synthetic migration vector, the position of the marked particle is updated to obtain the updated position of the particle. Based on the updated position of the particle, the density field is updated to complete the particle migration.
2. The SPH particle distribution optimization method based on gradient particle migration according to claim 1, characterized in that: The calculation formula of the Euclidean distance from the marked particle to the migration area boundary is: Among them, 𝑟𝑟 is the Euclidean distance from the marked particle to the boundary of the migration area, is the coordinate of the particle in the x direction, is the x-direction coordinate of the center coordinate, is the migration area range parameter, is the coordinate of the particle in the y direction, is the y-direction coordinate of the center coordinate, is the particle number, The calculation formula of the migration correction coefficient is: in, is the migration correction factor, is the absolute value function.
3. The SPH particle distribution optimization method based on gradient particle migration according to claim 1, characterized in that: The calculation formula of the synthetic migration vector is: in, is the synthetic migration vector, is the migration correction factor, is the position vector of the marker particle, The calculation formula of the position vector of the marker particle is: in, is the scalar coefficient, is the identity matrix, is the normal vector of the free surface of the i-th particle, are the first correction displacement vector, the second correction displacement vector and the third correction displacement vector respectively, is the fluid domain, is the particle number, For the transition domain, is the interface domain, is the regional characteristic value, The calculation formula of the scalar coefficient is: in, is the neighboring particle of the i-th marked particle, is the normal vector of the free surface of the i-th particle, is the normal vector of the free surface of the jth particle, is the virtual particle domain.
4. The SPH particle distribution optimization method based on gradient particle migration according to claim 1, characterized in that: The calculation formula for the position of the updated particle is: in, is the updated particle position, is the particle's position vector, is the migration correction factor.
5. The SPH particle distribution optimization method based on gradient particle migration according to claim 1, characterized in that: The migration area boundary meets the following conditions: in, is the x-direction coordinate of the center coordinate, is the migration area range parameter, is the y-direction coordinate of the center coordinate, is the x-coordinate of any point in the migration area, is the y-coordinate of any point in the migration area.
6. The SPH particle distribution optimization method based on gradient particle migration according to claim 1, characterized in that: The updating of the density field based on the updated particle position to complete particle migration includes: Based on the updated particle positions, the density field is updated, and the particles are migrated to the updated particle positions to optimize the calculation. During the migration execution, the virtual particles are coupled to the fixed boundary, the sponge layer absorbs the shock wave reflection, and the characteristic line interpolation non-reflection boundary.