A smooth particle simulation method for weakly compressible fluids based on particle position correction

Abstract

The present invention discloses a method for simulating weakly compressible fluid smooth particles based on particle position correction, comprising the following steps: 1. Smooth particle modeling of weakly compressible fluid; 2. Smooth particle position correction; 3. Using a computer to input the corrected position vector of the i-th particle, and using a position-based Verlet algorithm to perform time iteration on each particle in the SPH particle model to obtain an updated particle flow field distribution; and 4. Obtaining fluid state parameters from the updated particle flow field distribution. The method of the present invention has simple steps and a reasonable design. By correcting the position of fluid particles, a uniform distribution of simulated particles is achieved without introducing physical parameters, resulting in better applicability to complex flows and improved stability and accuracy of subsequent simulation calculations.

Description

Technical Field

[0001] The present invention belongs to the technical field of fluid simulation, and in particular relates to a smooth particle simulation method for weakly compressible fluid based on particle position correction. Background Art

[0002] Currently, fluid simulation is a widespread topic in scientific research and engineering. However, due to limitations in hardware and simulation methods, there remains a need for improved accuracy. Fluid simulation can typically be performed using either meshed or meshless methods. The meshless method, known as the SPH method, has attracted increasing attention and has yielded promising results in fluid dynamics and fluid-structure interaction numerical simulations.

[0003] However, due to the Lagrangian characteristics of the SPH method, particles move along trajectories, making it difficult to maintain a uniform distribution during calculations. The calculation accuracy of the SPH method depends on the distribution state of the particles. Uneven distribution will reduce the calculation accuracy of the method and the accuracy of numerical simulations of fluid dynamics problems. To solve this problem, one solution is to modify the calculation method, such as the modified smooth particle method. This method effectively improves the calculation accuracy at the boundary and solves the boundary defect problem. However, these methods cannot guarantee the conservation of the momentum of the system, that is, the interaction between particles is no longer antisymmetric. Another feasible solution is to adjust the particle distribution during the calculation so that the particles remain uniformly distributed during the simulation. For example, the velocity transport method achieves uniform distribution of particles by local average velocity or applying background pressure. Although these methods can improve the smoothness of the flow field, they inevitably introduce unknown physical parameters and need to be adjusted according to different model dynamics, reducing the universality of the method.

[0004] Therefore, there is currently a lack of a rationally designed smooth particle simulation method for weakly compressible fluids based on particle position correction. By correcting the position of fluid particles, a uniform distribution of simulated particles can be achieved without introducing physical parameters, making it more applicable to complex flows and improving the stability and accuracy of subsequent simulation calculations. Summary of the Invention

[0005] The technical problem to be solved by the present invention is to address the deficiencies in the above-mentioned prior art and provide a smooth particle simulation method for weakly compressible fluids based on particle position correction. The method has simple steps and a reasonable design. By correcting the position of fluid particles, a uniform distribution of simulated particles is achieved without introducing physical parameters. The method has better applicability to complex flows and improves the stability and accuracy of subsequent simulation calculations.

[0006] To solve the above technical problems, the present invention adopts a technical solution: a method for simulating smooth particles of weakly compressible fluids based on particle position correction, characterized in that the method comprises the following steps:

[0007] Step 1: Smooth particle modeling of weakly compressible fluid:

[0008] Step 101: Use VC++ software to build a geometric model on a computer;

[0009] Step 102: using a computer and using the SPH framework in VC++ software to set parameters for the geometric model, obtaining a Taylor-Green vortex model based on a weakly compressible fluid; wherein the parameters include initial density, viscosity, Reynolds number, initial velocity, and speed of sound; and applying periodic boundary conditions to the boundaries of the Taylor-Green vortex model based on the weakly compressible fluid;

[0010] Step 103: Generate an SPH particle model based on a Taylor-Green vortex model of a weakly compressible fluid using a computer;

[0011] Step 104: Using a computer, a Verlet neighbor list algorithm is used to traverse each particle in the SPH particle model according to a set search radius to obtain a Verlet linked list.

[0012] Step 105: Using a computer to obtain the neighboring particles of each particle from the Verlet list;

[0013] Step 2: Smooth particle position correction:

[0014] Step A: Use a computer to simulate and analyze the SPH particle model using the SPH method to obtain the SPH basic equation of the field variable space derivative as follows:

[0015] Among them, V j represents the volume of the jth particle among the neighboring particles of the i-th particle, N represents the total number of neighboring particles of the i-th particle, j and N are both positive integers, and 1≤j≤N, represents the gradient operator, f(x i ) represents the field variable function of the i-th particle, f(x j ) represents the field variable function of the jth particle, W ij represents the kernel function between the i-th particle and the j-th particle, and π represents pi, h represents smooth length, x i represents the position of the i-th particle, x j represents the position of the jth particle;

[0016] Step B: Use a computer to set the particle's field variable function to 1, then f(x i )=f(x j )=1, convert the SPH basic equation in step A into the following formula:

[0017]

[0018] Step C: Use a computer to add the correction value of the kernel function derivative to the right side of the equal sign in the formula of step B. So that the right side of the equal sign of the formula in step B in, dx ij =dx i -dx j , dx ij represents the difference in the corrected position vector between the i-th particle and the j-th particle; dx i represents the corrected position vector of the i-th particle; dx j represents the corrected position value vector of the jth particle among the neighboring particles of the i-th particle; e ij represents the unit vector along the direction of the i-th particle and the j-th particle; represents the product operation, 0 represents the zero vector;

[0019] Step D: Use a computer to establish the integral error E of the i-th particle i The expression is as follows

[0020]

[0021] Step E: Using a computer to calculate the integral error E of the i-th particle using the gradient descent method i Calculate and obtain the corrected position vector dx of the i-th particle i ;

[0022] Step 3: Use the computer to input the corrected position vector dx of the i-th particle i , the position-based Verlet algorithm is used to perform time iteration on each particle in the SPH particle model to obtain the updated particle flow field distribution;

[0023] Step 4: Obtain fluid state parameters from the updated particle flow field distribution:

[0024] Step 401: Using a computer to obtain the density ρ, position and velocity of the particle from the updated particle flow field distribution;

[0025] Step 402: Use a computer to calculate the value of P=c 2 (ρ-ρ 0 ), and the pressure P of the particle is obtained; where c represents the speed of sound, ρ 0 Represents the initial density of particles.

[0026] The above-mentioned method for simulating smooth particles of weakly compressible fluids based on particle position correction is characterized in that: Step E, the specific process is as follows:

[0027] Step E01: Use computer to set coefficient matrix Then convert the expression in step D to the following:

[0028]

[0029] Step E02: Use a computer to perform gradient calculation on the integral error of the i-th particle to obtain the gradient of the integral error of the i-th particle. and Where T represents the transpose of the vector;

[0030] Use a computer to set the learning rate to k, then

[0031] Step E03: Use computer to and Substitute into step E01 to obtain the learning rate k;

[0032] Step E04: Use a computer to substitute the learning rate k into Get the corrected position vector dx of the i-th particle i .

[0033] The above-mentioned method for simulating smooth particles of weakly compressible fluids based on particle position correction is characterized in that: Step 3, the specific process is as follows:

[0034] Step 301: Use a computer to calculate the Get the time step Δt of convection control ad ; where h represents the smoothing length, |v max represents the maximum predicted velocity of the flow, λ represents the viscosity, and min(·) represents the minimum value;

[0035] Step 302: Using a computer according to the formula Get the time step Δt of the sound speed control ac ; Where c represents the speed of sound;

[0036] Step 303: Use a computer to control the time step Δt of the convection ad and the time step Δt of the acoustic speed control ac Select a smaller value as the time step Δt;

[0037] Step 304: Use a computer to obtain the initial density of the i-th particle from the SPH particle model. The initial velocity of the i-th particle The initial position of the i-th particle

[0038] Step 305: Use a computer to set the initial time of the pressure relaxation time to zero, and set the time step Δt of the pressure relaxation time and the convection control adMake a judgment, if the pressure relaxation time is less than the convection control time step Δt ad , execute step 306; otherwise, execute step 309;

[0039] Step 306: Using a computer according to the formula Get the density of the i-th particle at the intermediate moment between the n-th time step and the n+1-th time step in, represents the density of the i-th particle at the n-th time step, represents the derivative with respect to time t;

[0040] According to the formula Get the position of the i-th particle at the intermediate moment between the n-th time step and the n+1-th time step in, represents the position of the i-th particle at the n-th time step;

[0041] According to the formula Get the velocity of the i-th particle at the n+1-th time step in, represents the velocity of the i-th particle at the n-th time step, where n is a natural number;

[0042] Step 307: Using a computer Get the density of the i-th particle at the n+1-th time step and the position of the i-th particle at the n+1-th time step Among them, the iteration and The post-pressure relaxation time is updated to (n+1)Δt;

[0043] Step 308: If the iteration yields and The post-pressure relaxation time (n+1)Δt is less than the convection control time step Δt ad , perform the next iteration according to step 306 and step 307 to obtain the updated particle flow field distribution, otherwise, execute step 309;

[0044] Step 309: Using a computer and using the Verlet neighbor list algorithm, the updated particle flow field distribution in the SPH particle model is traversed according to the set search radius to obtain an updated Verlet linked list.

[0045] The above-mentioned method for simulating smooth particles of weakly compressible fluids based on particle position correction is characterized in that: Step 103, the specific process is as follows:

[0046] The grid side length is set by computer, and the Taylor-Green vortex model based on weakly compressible fluid is divided into grids of equal size. A particle is generated in each grid to obtain an SPH particle model; wherein the diameter of each particle is equal to the side length of the grid, the center of each particle is the center of each grid, and the center coordinates of each grid are marked as the position coordinates of the particle.

[0047] Compared with the prior art, the present invention has the following advantages:

[0048] 1. The solution of the particle correction position of the present invention is based on the SPH zero-order consistency equation of the distribution of neighboring particles, and the particle relaxation efficiency is more efficient.

[0049] 2. The gradient descent method is used in the process of solving the particle position correction of the present invention, and the calculation process is more stable.

[0050] 3. The correction of the present invention is based on the local particle distribution, without the introduction of physical parameters, and has better applicability to complex flows.

[0051] In summary, the method of the present invention has simple steps and reasonable design. By correcting the position of fluid particles, a uniform distribution of simulated particles is achieved without introducing any physical parameters. It has better applicability to complex flows and improves the stability and accuracy of subsequent simulation calculations.

[0052] The technical solution of the present invention is further described in detail below through the accompanying drawings and embodiments. BRIEF DESCRIPTION OF THE DRAWINGS

[0053] Figure 1 It is a flowchart of the method of the present invention. DETAILED DESCRIPTION

[0054] like Figure 1 As shown, the present invention provides a method for simulating smooth particles of weakly compressible fluids based on particle position correction, comprising the following steps:

[0055] Step 1: Smooth particle modeling of weakly compressible fluid:

[0056] Step 101: Use VC++ software to build a geometric model on a computer;

[0057] Step 102: using a computer and using the SPH framework in VC++ software to set parameters for the geometric model, obtaining a Taylor-Green vortex model based on a weakly compressible fluid; wherein the parameters include initial density, viscosity, Reynolds number, initial velocity, and speed of sound; and applying periodic boundary conditions to the boundaries of the Taylor-Green vortex model based on the weakly compressible fluid;

[0058] Step 103: Generate an SPH particle model based on a Taylor-Green vortex model of a weakly compressible fluid using a computer;

[0059] Step 104: Using a computer, a Verlet neighbor list algorithm is used to traverse each particle in the SPH particle model according to a set search radius to obtain a Verlet linked list.

[0060] Step 105: Using a computer to obtain the neighboring particles of each particle from the Verlet list;

[0061] Step 2: Smooth particle position correction:

[0062] Step A: Use a computer to simulate and analyze the SPH particle model using the SPH method to obtain the SPH basic equation of the field variable space derivative as follows:

[0063] Among them, V j represents the volume of the jth particle among the neighboring particles of the i-th particle, N represents the total number of neighboring particles of the i-th particle, j and N are both positive integers, and 1≤j≤N, represents the gradient operator, f(x i ) represents the field variable function of the i-th particle, f(x j ) represents the field variable function of the jth particle, W ij represents the kernel function between the i-th particle and the j-th particle, and π represents pi, h represents smooth length, x i represents the position of the i-th particle, x j represents the position of the jth particle;

[0064] Step B: Use a computer to set the particle's field variable function to 1, then f(x i )=f(x j )=1, convert the SPH basic equation in step A into the following formula:

[0065]

[0066] Step C: Use a computer to add the correction value of the kernel function derivative to the right side of the equal sign in the formula of step B. So that the right side of the equal sign of the formula in step B in, dx ij =dx i -dx j , dx ij represents the difference in the corrected position vector between the i-th particle and the j-th particle; dx i represents the corrected position vector of the i-th particle; dx j represents the corrected position value vector of the jth particle among the neighboring particles of the i-th particle; e ijrepresents the unit vector along the direction of the i-th particle and the j-th particle; represents the product operation, 0 represents the zero vector;

[0067] Step D: Use a computer to establish the integral error E of the i-th particle i The expression is as follows

[0068]

[0069] Step E: Using a computer to calculate the integral error E of the i-th particle using the gradient descent method i Calculate and obtain the corrected position vector dx of the i-th particle i ;

[0070] Step 3: Use the computer to input the corrected position vector dx of the i-th particle i , the position-based Verlet algorithm is used to perform time iteration on each particle in the SPH particle model to obtain the updated particle flow field distribution;

[0071] Step 4: Obtain fluid state parameters from the updated particle flow field distribution:

[0072] Step 401: Using a computer to obtain the density ρ, position and velocity of the particle from the updated particle flow field distribution;

[0073] Step 402: Use a computer to calculate the value of P=c 2 (ρ-ρ 0 ), and the pressure P of the particle is obtained; where c represents the speed of sound, ρ 0 Represents the initial density of particles.

[0074] In this embodiment, during specific implementation, step E, the specific process is as follows:

[0075] Step E01: Use computer to set coefficient matrix Then convert the expression in step D to the following:

[0076]

[0077] Step E02: Use a computer to perform gradient calculation on the integral error of the i-th particle to obtain the gradient of the integral error of the i-th particle. and Where T represents the transpose of the vector;

[0078] Use a computer to set the learning rate to k, then

[0079] Step E03: Use computer to and Substitute into step E01 to obtain the learning rate k;

[0080] Step E04: Use a computer to substitute the learning rate k into Get the corrected position vector dx of the i-th particle i .

[0081] In this embodiment, during the specific implementation, step three, the specific process is as follows:

[0082] Step 301: Use a computer to calculate the Get the time step Δt of convection control ad ; where h represents the smoothing length, |v max represents the maximum predicted velocity of the flow, λ represents the viscosity, and min(·) represents the minimum value;

[0083] Step 302: Using a computer according to the formula Get the time step Δt of the sound speed control ac ; Where c represents the speed of sound;

[0084] Step 303: Use a computer to control the time step Δt of the convection ad and the time step Δt of the acoustic speed control ac Select a smaller value as the time step Δt;

[0085] Step 304: Use a computer to obtain the initial density of the i-th particle from the SPH particle model. The initial velocity of the i-th particle The initial position of the i-th particle

[0086] Step 305: Use a computer to set the initial time of the pressure relaxation time to zero, and set the time step Δt of the pressure relaxation time and the convection control ad Make a judgment, if the pressure relaxation time is less than the convection control time step Δt ad , execute step 306; otherwise, execute step 309;

[0087] Step 306: Using a computer according to the formula Get the density of the i-th particle at the intermediate moment between the n-th time step and the n+1-th time step in, represents the density of the i-th particle at the n-th time step, represents the derivative with respect to time t;

[0088] According to the formula Get the position of the i-th particle at the intermediate moment between the n-th time step and the n+1-th time step in, represents the position of the i-th particle at the n-th time step;

[0089] According to the formula Get the velocity of the i-th particle at the n+1-th time step in, represents the velocity of the i-th particle at the n-th time step, where n is a natural number;

[0090] Step 307: Using a computer Get the density of the i-th particle at the n+1-th time step and the position of the i-th particle at the n+1-th time step Among them, the iteration and The post-pressure relaxation time is updated to (n+1)Δt;

[0091] Step 308: If the iteration yields and The post-pressure relaxation time (n+1)Δt is less than the convection control time step Δt ad , perform the next iteration according to step 306 and step 307 to obtain the updated particle flow field distribution, otherwise, execute step 309;

[0092] Step 309: Using a computer and using the Verlet neighbor list algorithm, the updated particle flow field distribution in the SPH particle model is traversed according to the set search radius to obtain an updated Verlet linked list.

[0093] In this embodiment, during specific implementation, step 103, the specific process is as follows:

[0094] The grid side length is set by computer, and the Taylor-Green vortex model based on weakly compressible fluid is divided into grids of equal size. A particle is generated in each grid to obtain an SPH particle model; wherein the diameter of each particle is equal to the side length of the grid, the center of each particle is the center of each grid, and the center coordinates of each grid are marked as the position coordinates of the particle.

[0095] In this embodiment, the SPH method is the Smoothed Particle Hydrodynamics Method.

[0096] In this embodiment, during specific implementation, the geometry in step 101 is a square, and the side length L of the geometric model is 1m;

[0097] In this embodiment, during the specific implementation, the initial density in step 102 is 1 kg / m 3, viscosity is 0.02 Pa·s, Reynolds number is 100 or 10000, initial velocity is 1 m / s, and speed of sound is 10 m / s.

[0098] In this embodiment, during specific implementation, the grid side length in step 103 is L / 50, the number of particles is 2500, and the particle diameter is also L / 50;

[0099] In step 104, a radius of 2 times the smooth length is searched;

[0100] In this embodiment, during specific implementation, the smoothing length is set to 1.5 times the length of the grid side.

[0101] In this embodiment, during specific implementation, either the grid side length is L / 100 and the number of particles is 10,000, or the grid side length is L / 200 and the number of particles is 40,000.

[0102] In this embodiment, during specific implementation, theoretically, the left side of the equation in step B is equal to zero, but in actual practice, due to errors, the right side of the equation is usually not equal to 0, so a correction amount of the kernel function derivative is added to make the right side of the equation in step B equal to 0, so that the equation holds.

[0103] In this embodiment, when it is implemented, the maximum predicted flow velocity |v max The value is 1m / s.

[0104] In summary, the method of the present invention has simple steps and reasonable design. By correcting the position of fluid particles, a uniform distribution of simulated particles is achieved without introducing any physical parameters. It has better applicability to complex flows and improves the stability and accuracy of subsequent simulation calculations.

[0105] The above description is only a preferred embodiment of the present invention and does not limit the present invention in any way. Any simple modification, change and equivalent structural change made to the above embodiment based on the technical essence of the present invention shall still fall within the scope of protection of the technical solution of the present invention.

Claims

1. A smooth particle simulation method for weakly compressible fluid based on particle position correction, characterized in that: The method comprises the following steps: Step 1: Smooth particle modeling of weakly compressible fluid: Step 101: Use VC++ software to build a geometric model on a computer; Step 102: using a computer and using the SPH framework in VC++ software to set parameters for the geometric model, obtaining a Taylor-Green vortex model based on a weakly compressible fluid; wherein the parameters include initial density, viscosity, Reynolds number, initial velocity, and speed of sound; and applying periodic boundary conditions to the boundaries of the Taylor-Green vortex model based on the weakly compressible fluid; Step 103: Generate an SPH particle model based on a Taylor-Green vortex model of a weakly compressible fluid using a computer; Step 104: Using a computer, a Verlet neighbor list algorithm is used to traverse each particle in the SPH particle model according to a set search radius to obtain a Verlet linked list. Step 105: Using a computer to obtain the neighboring particles of each particle from the Verlet list; Step 2: Smooth particle position correction: Step A: Use a computer to simulate and analyze the SPH particle model using the SPH method to obtain the SPH basic equation of the field variable space derivative, as follows: Among them, V j represents the volume of the jth particle among the neighboring particles of the i-th particle, N represents the total number of neighboring particles of the i-th particle, j and N are both positive integers, and 1≤j≤N, represents the gradient operator, f(x i ) represents the field variable function of the i-th particle, f(x j ) represents the field variable function of the jth particle, W ij represents the kernel function between the i-th particle and the j-th particle, and π represents pi, h represents smooth length, x i represents the position of the i-th particle, x j represents the position of the jth particle; Step B: Use a computer to set the particle's field variable function to 1, then f(x i )=f(x j )=1, convert the SPH basic equation in step A into the following formula: Step C: Use a computer to add the correction value of the kernel function derivative to the right side of the equal sign in the formula of step B. So that the right side of the equal sign of the formula in step B in, dx ij =dx i -dx j , dx ij represents the difference in the corrected position vector between the i-th particle and the j-th particle; dx i represents the corrected position vector of the i-th particle; dx j represents the corrected position value vector of the jth particle among the neighboring particles of the i-th particle; e ij represents the unit vector along the direction of the i-th particle and the j-th particle; represents the product operation, 0 represents the zero vector; Step D: Use a computer to establish the integral error E of the i-th particle i The expression is as follows Step E: Using a computer to calculate the integral error E of the i-th particle using the gradient descent method i Calculate and obtain the corrected position vector dx of the i-th particle i ; Step 3: Use the computer to input the corrected position vector dx of the i-th particle i , the position-based Verlet algorithm is used to perform time iteration on each particle in the SPH particle model to obtain the updated particle flow field distribution; Step 4: Obtain fluid state parameters from the updated particle flow field distribution: Step 401: Using a computer to obtain the density ρ, position and velocity of the particle from the updated particle flow field distribution; Step 402: Use a computer to calculate the value of P=c 2 (ρ-ρ 0 ), and the pressure P of the particle is obtained; where c represents the speed of sound, ρ 0 Indicates the initial density of particles.

2. A method for simulating smooth particles in weakly compressible fluids based on particle position correction according to claim 1, characterized in that: Step E: The specific process is as follows: Step E01: Use computer to set coefficient matrix Then convert the expression in step D to the following: Step E02: Use a computer to perform gradient calculation on the integral error of the i-th particle to obtain the gradient of the integral error of the i-th particle. and Where T represents the transpose of the vector; Use a computer to set the learning rate to k, then Step E03: Use computer to and Substitute into step E01 to obtain the learning rate k; Step E04: Use a computer to substitute the learning rate k into Get the corrected position vector dx of the i-th particle i .

3. The method for simulating smooth particles in weakly compressible fluids based on particle position correction according to claim 1, characterized in that: Step 3: The specific process is as follows: Step 301: Use a computer to calculate the Get the time step Δt of convection control ad ; where h represents the smoothing length, |v max represents the maximum predicted flow velocity, λ represents viscosity, and min(·) represents the minimum value; Step 302: Using a computer according to the formula Get the time step Δt of the sound speed control ac ; Where c represents the speed of sound; Step 303: Use a computer to control the time step Δt of the convection ad and the time step Δt of the acoustic speed control ac Select a smaller value as the time step Δt; Step 304: Use a computer to obtain the initial density of the i-th particle from the SPH particle model. The initial velocity of the i-th particle The initial position of the i-th particle Step 305: Use a computer to set the initial time of the pressure relaxation time to zero, and set the time step Δt of the pressure relaxation time and the convection control ad Make a judgment, if the pressure relaxation time is less than the convection control time step Δt ad , execute step 306; otherwise, execute step 309; Step 306: Using a computer according to the formula Get the density of the i-th particle at the intermediate moment between the n-th time step and the n+1-th time step in, represents the density of the i-th particle at the n-th time step, represents the derivative with respect to time t; According to the formula Get the position of the i-th particle at the intermediate moment between the n-th time step and the n+1-th time step in, represents the position of the i-th particle at the n-th time step; According to the formula Get the velocity of the i-th particle at the n+1-th time step in, represents the velocity of the i-th particle at the n-th time step, where n is a natural number; Step 307: Using a computer Get the density of the i-th particle at the n+1-th time step and the position of the i-th particle at the n+1-th time step Among them, the iteration and The post-pressure relaxation time is updated to (n+1)Δt; Step 308: If the iteration yields and The post-pressure relaxation time (n+1)Δt is less than the convection control time step Δt ad , perform the next iteration according to step 306 and step 307 to obtain the updated particle flow field distribution, otherwise, execute step 309; Step 309: Using a computer and using the Verlet neighbor list algorithm, the updated particle flow field distribution in the SPH particle model is traversed according to the set search radius to obtain an updated Verlet linked list.

4. The method for simulating smooth particles in weakly compressible fluids based on particle position correction according to claim 1, characterized in that: Step 103, the specific process is as follows: The grid side length is set by computer, and the Taylor-Green vortex model based on weakly compressible fluid is divided into grids of equal size. A particle is generated in each grid to obtain an SPH particle model; wherein the diameter of each particle is equal to the side length of the grid, the center of each particle is the center of each grid, and the center coordinates of each grid are marked as the position coordinates of the particle.

Images