Fluid-structure interaction method for simulating autonomous swimming of flexible organisms based on improved immersed boundary method

Abstract

The application discloses a flow-solid coupling method for simulating autonomous swimming of flexible organisms based on an improved immersed boundary method, and the method comprises the following steps: A, initial parameters are specified at each interface; B, a distribution technique is used to divide the control equation of the immersed boundary-lattice Boltzmann flux solver into a prediction step and a correction step; C, in the prediction step, the intermediate flow field velocity is solved by using the lattice Boltzmann flux solver; D, in the correction step, the intermediate flow field velocity is corrected by using the improved boundary condition forced immersed boundary method to obtain a new flow field velocity; E, if the simulation time is reached, the whole simulation process is ended, and if the simulation time is not reached, the step B is returned to. The flow-solid coupling method can meet the non-slip boundary condition, does not cause the non-physical flow line penetration phenomenon, and eliminates the complicated matrix solving process, so that a large amount of calculation cost can be saved, and the method is suitable for three-dimensional incompressible flow simulation of immersed objects.

Description

TECHNICAL FIELD

[0001] The present application belongs to the technical field of fluid mechanics, and particularly relates to a flow-solid coupling method for simulating autonomous swimming of flexible organisms based on an improved immersed boundary method. BACKGROUND

[0002] The study on the flow coupling motion characteristics of flexible organisms such as swimming or flying organisms is usually applied to bionics and reveals the mechanical and deformation principles of biological motion. Fish autonomous swimming is a complex six-coupling problem, including the active control of the simulation machine motion of fish muscle, the coupling effect of fish body motion and fluid. The traditional flow-solid coupling calculation method has certain limitations in processing large deformation motion of biological tissues, complex structure and calculation efficiency. Therefore, how to establish a theoretical model and a calculation method capable of ideally simulating the interaction between biological machine motion and fluid is the focus of the research in this field.

[0003] The immersed boundary method developed in the 1970s is one of the ideal methods for solving the flow-solid coupling problem of biological muscle. Its main advantages are simple interface processing and mesh generation, high calculation efficiency, and suitability for complex geometric shapes and deformation, especially in the processing of large deformation problems, and it has been widely used in the field of biological fluid mechanics. Fish swimming belongs to a relatively typical fluid-structure interaction problem. Therefore, in order to better simulate fish swimming, a model considering the flow-solid interaction mechanism during fish swimming needs to be established in order to more realistically reflect the thrust generated by the interaction between fish muscle motion and fluid. Due to the complexity of fish swimming problems, previous research work calculated fish motion and fluid separately without considering the flow-solid coupling problem during fish swimming. Or when calculating fish swimming, the fish swimming form and speed are given by using certain control equations, which obviously restricts the freedom and flexibility reflected during fish swimming, and there is still a large difference from the characteristics and dynamic behavior of real fish swimming. Recently, someone has proposed an immersed boundary-lattice Boltzmann flux solver method for simulating three-dimensional incompressible flow with complex geometric shapes and moving boundaries, which has high accuracy for simulating autonomous swimming of flexible organisms. However, in this method, its limitation lies in the implicit solving process and complexity in handling large-scale or moving boundary problems. The no-slip boundary condition is implicitly contained in a well-defined linear equation system, where a square matrix is formed containing Lagrangian-Eulerian field information. When dealing with moving boundaries, it is necessary to assemble and implicitly solve the matrix at each time step, which is very time-consuming. In addition, with the increase of the number of Lagrangian points distributed on the immersed boundary, the calculation cost increases rapidly. SUMMARY

[0004] The application aims at solving the problems in the prior art, and provides a flow-solid coupling method for simulating autonomous swimming of flexible organisms based on an improved immersed boundary method.

[0005] The application aims at solving the problems in the prior art, and provides a flow-solid coupling method for simulating autonomous swimming of flexible organisms based on an improved immersed boundary method.

[0006] The application aims at solving the problems in the prior art, and provides a flow-solid coupling method for simulating autonomous swimming of flexible organisms based on an improved immersed boundary method.

[0007] A, initial parameters are specified at each interface;

[0008] B, a distribution technique is used to divide the control equation of the immersed boundary-lattice Boltzmann flux solver into a prediction step and a correction step;

[0009] C, in the prediction step, the intermediate flow field velocity is solved by using the lattice Boltzmann flux solver;

[0010] D, in the correction step, the intermediate flow field velocity is corrected by using the improved boundary condition forced immersed boundary method to obtain a new flow field velocity;

[0011] E, after the simulation time is reached, the whole simulation process is ended; if the simulation time is not reached, the step B is returned to.

[0012] The simulation time in the step E is the product of the flow direction time step and the iteration step number in the whole iteration process, since the whole flow-solid coupling method for simulating autonomous swimming of flexible organisms is an iteration calculation process around the velocity value, and the determination criterion for ending the iteration calculation process is that the flow field velocity value tends to be stable or the variation fluctuation is stable, for example, the variation curve of the flow field velocity value becomes a straight line or the up and down fluctuation is stable.

[0013] The initial parameters in the step A include the flow direction time step, the fluid density, the fluid velocity, the pressure, the dynamic viscosity and the Reynolds number.

[0014] The control equation of the immersed boundary-lattice Boltzmann flux solver in the step B is:

[0015]

[0016] In the formula (1), ρ is the fluid density, u is the fluid velocity, p is the pressure, μ is the dynamic viscosity, T is the temperature, and f is the restoring force determined by the immersed boundary method and used for explaining the boundary effect of the immersed object.

[0017] According to the multi-scale Chapman-Enskog analysis, the density distribution function f α The relationship between the flux ρu in equation (1) and the density distribution function f

[0018]

[0019] In equation (2): I is the unit tensor; f α is the density distribution function along the α direction; is the equilibrium density distribution function along the α direction; is the non-equilibrium density distribution function along the α direction; τ is the single relaxation parameter; δ t is the flow time step; e α is the particle velocity along the α direction; β and γ are two directions in e α ; p is the pressure.

[0020] In a three-dimensional incompressible flow field, the particle velocity e α along the α direction is given by the following equation (3) using the D3Q15 lattice velocity model:

[0021]

[0022] The equilibrium density distribution function f along the α direction is given by the following equation (4):

[0023]

[0024] In equation (4): w α is the coefficient of the equilibrium density distribution function, c s is the sound speed, and w α and c s values depend on the lattice velocity model;

[0025] The single relaxation parameter τ is related to the dynamic viscosity μ of the fluid and is given by the following equation (5):

[0026]

[0027] The non-equilibrium density distribution function f along the α direction is given by the following equation (6):

[0028]

[0029] The equations (1) and (2) involved in the control equations of the wetting boundary-lattice Boltzmann flux solver are divided into a prediction step and a correction step using the distribution technique, where the prediction step is given by the following equation (7) and the correction step is given by the following equation (8):

[0030]

[0031]

[0032] In formula (7) and formula (8): Δt is a time step; u * is an intermediate flow field velocity; ρ n is a density of a current time step; ρ n+1 is a density of a next time step; u n is a flow field velocity of the current time step; u n+1 is a flow field velocity of the next time step.

[0033] The equation (7) of the prediction step is solved by using a lattice Boltzmann flux solver to obtain the density ρ n+1 and the intermediate flow field velocity u * of the next time step; after the equation (7) of the prediction step is discretized by finite volume in a control volume Ω i , the following is obtained:

[0034]

[0035] In formula (9): W * =(ρ n+1 ,ρ n+1 u * ); ΔV i is a volume of the control volume Ω i ; G k is a flux; ΔS k is an area of the kth interface surrounding the control volume Ω i ; and n is a unit outer normal vector of the kth interface.

[0036] The flux G k in formula (9) is locally reconstructed at each interface, and the specific reconstruction is as follows: the second-order approximation of the non-equilibrium density distribution function is calculated by using Taylor series expansion as follows:

[0037]

[0038] In formula (10): r refers to a position of a particle at a current time t; t refers to a current time; δ t refers to a flow field time step;

[0039] The density and the velocity are locally reconstructed at the interface by using a solution in the lattice Boltzmann method:

[0040]

[0041] Based on the modified equation (8), the intermediate flow field velocity is modified using the submerged boundary method with improved boundary conditions, and a new flow field velocity u is obtained. n+1 The new flow field velocity u n+1 It is determined by the intermediate flow field velocity u * It is obtained by summing and correcting the flow field velocity δu:

[0042] u n+1 =u * +δu (12)

[0043] Wall velocity at Lagrange points The velocity u on the Euler grid n+1 (r j Interpolated from:

[0044]

[0045] In equation (13): h is the mesh size of the Eulerian mesh; N is the number of Lagrange points, M is the number of Eulerian points; r j It is the location of the Euler point. It is the location of the Lagrange point; D ij It is a continuous kernel function, given by the following equation (14):

[0046]

[0047]

[0048] Similarly, the correction velocity δu(r) of the Euler point j Corrected velocity at Lagrange points Interpolation yields:

[0049]

[0050] Based on formulas (12), (13), (14), and (16), the solution equation is obtained:

[0051]

[0052] make If X is the unknown quantity, then the solution equation (17) can be written as a matrix relation (18):

[0053] AX = B (18)

[0054] In equation (18), matrices A and B represent:

[0055]

[0056]

[0057] The matrix relationship (18) is solved by using improved matrix solving techniques, and the specific steps are as follows: the matrix A is converted into a diagonal matrix A diag The diagonal matrix A diag is given by formula (21):

[0058] A diag = diag(d1 d2 … d N ) (21)

[0059]

[0060] According to formula (18)-formula (22), the unknown quantity X is efficiently solved and obtained:

[0061]

[0062] In formula (23), B=(B1 B2 … B N );

[0063] The modified flow field velocity δu is obtained based on the obtained unknown quantity X, and then the new flow field velocity u n+1 is obtained.

[0064] Compared with the prior art, the present application has the following advantages:

[0065] The flow-solid coupling method for simulating the autonomous swimming of flexible organisms provided by the present application improves the original boundary condition forced immersed boundary method by using the improved boundary condition forced immersed boundary method to modify the intermediate flow field velocity on the basis of the original immersed boundary-lattice Boltzmann flux solver method, effectively simulates the incompressible flow with a moving boundary, and saves a large amount of calculation cost.

[0066] The method provided by this invention retains the advantages of the original submerged boundary-lattice Boltzmann flux solver. It not only overcomes the problems of reduced solution efficiency due to the slow convergence of the Poisson equation in the Navier-Stockl solver, and the inconvenience and workload brought to programming by the need for staggered meshes during the solution process, but also overcomes the problem that the lattice Boltzmann method solver can only apply uniform meshes, thus inevitably requiring a large number of mesh points for complex fluid-structure interaction problems, increasing computational costs. By combining the Navier-Stockl solver and the lattice Boltzmann method solver, both of which provide a clear and rigorous physical basis for incompressible flows, and by combining the simplicity, ease of implementation, inherent dynamics, and parallelism of the lattice Boltzmann method solver, it becomes a flexible and efficient solver for simulating incompressible flows on non-uniform meshes. Furthermore, this method incorporates the submerged boundary method, avoiding the tedious process of reconstructing the mesh, saving significant computational costs and coding work. Attached Figure Description

[0067] To more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments recorded in the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0068] Appendix Figure 1 The flowchart shows the fluid-structure coupling method for simulating autonomous swimming of flexible organisms based on the improved immersion boundary method of the present invention.

[0069] Appendix Figure 2 This is a schematic diagram of the lattice velocity model used in this invention;

[0070] Appendix Figure 3 This is a schematic diagram of the partial reconstruction of the LBM solution at the interface between two control units in an embodiment of the present invention;

[0071] Appendix Figure 4 This is a geometric schematic diagram of a stingray in an embodiment of the present invention;

[0072] Appendix Figure 5 This is a comparison chart of drag coefficients at different wave numbers in an embodiment of the present invention. Detailed Implementation

[0073] To enable those skilled in the art to better understand the technical solutions of the present invention, the technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort should fall within the scope of protection of the present invention.

[0074] Appendix Figure 1 This is a flowchart of the fluid-structure interaction method for simulating autonomous swimming of flexible organisms based on the improved immersion boundary method of the present invention. While the present invention provides method operation steps or apparatus structures as shown in the following embodiments or figures, more or fewer operation steps or module units may be included in the method or apparatus based on conventional or non-inventive effort. In steps or structures where there is no logically necessary causal relationship, the execution order of these steps or the module structure of the apparatus is not limited to the execution order or module structure shown in the embodiments or figures of the present invention. When the method or module structure is applied in actual devices or end products, it can be executed sequentially or in parallel according to the method or module structure shown in the embodiments or figures (e.g., in a parallel processor or multi-threaded processing environment, or even a distributed processing implementation environment).

[0075] The idea behind this invention is to improve the matrix solving strategy in the original submerged boundary-lattice Boltzmann flux solver. By using the improved matrix solving technique to solve the submerged boundary method to correct the intermediate speed, the computational efficiency is improved and the computational cost is reduced significantly.

[0076] It should be noted that, unless otherwise specified, the tables below representing the various representations in this embodiment are merely for differentiation and have no special meaning.

[0077] Example

[0078] like Figures 1-3 As shown, this embodiment discloses a fluid-structure interaction method for simulating autonomous swimming of flexible organisms based on an improved immersion boundary method. The steps of this fluid-structure interaction method are as follows:

[0079] A. Specify initial parameters at each interface. These initial parameters include the flow time step, fluid density, fluid velocity, pressure, dynamic viscosity, and Reynolds number. The governing equations for the wetting boundary-lattice Boltzmann flux solver are:

[0080]

[0081] In equation (1): ρ is the fluid density; u is the fluid velocity; p is the pressure; μ is the dynamic viscosity; T is the temperature; f is the restoring force determined by the immersion boundary method, used to explain the boundary action of the immersed object;

[0082] According to the multi-scale Chapman-Enskog analysis, the density distribution function f α The relationship between the flux ρu in formula (1) and the density distribution function f

[0083]

[0084] In formula (2): I is a unit tensor; f α is the density distribution function along the α direction; is the equilibrium density distribution function along the α direction; is the non-equilibrium density distribution function along the α direction; τ is a single relaxation parameter; δ t is the flow direction time step; e α is the particle velocity along the α direction; β and γ are two directions in e α ; p is the pressure.

[0085] In a three-dimensional incompressible flow field, the particle velocity e α along the α direction is given by the following formula (3) by applying the D3Q15 lattice velocity model:

[0086]

[0087] The equilibrium density distribution function f along the α direction is given by the following formula (4):

[0088]

[0089] In formula (4): w α is the coefficient of the equilibrium density distribution function, c s is the sound speed, and the values of w α and c s depend on the lattice velocity model, in this embodiment, w0=2 / 9, w 1-6 =1 / 9, w 7-14 =1 / 72,

[0090] The single relaxation parameter τ is related to the dynamic viscosity μ of the fluid, and is given by the following formula (5):

[0091]

[0092] The non-equilibrium density distribution function f along the α direction is given by the following formula (6):

[0093]

[0094] B. The formulas (1) and (2) involved in the control equation of the immersed boundary-lattice Boltzmann flux solver are divided into a prediction step and a correction step by using distribution technique, wherein the prediction step is given by the following formula (7), and the correction step is given by the following formula (8):

[0095]

[0096]

[0097] In the formulas (7) and (8), △t is a time step; u * is an intermediate flow field velocity; p n is a density of a current time step; p n+1 is a density of a next time step; u n is a flow field velocity of the current time step; u n+1 is a flow field velocity of the next time step.

[0098] C. In the prediction step, the intermediate flow field velocity is solved by the lattice Boltzmann flux solver; specifically, the equation (7) of the prediction step is solved by the lattice Boltzmann flux solver to obtain the density p n+1 and the intermediate flow field velocity u * of the next time step; in a control body Ω i , the equation (7) of the prediction step is discretized by finite volume as follows after discretization:

[0099]

[0100] In the formula (9), W * =(p n+1 , p n+1 u * ); ΔV i is a volume of the control body Ω i ; G k is a flux; ΔS k is an area of the kth interface surrounding the control body Ω i ; and n is a unit outer normal vector of the kth interface.

[0101] The flux G k in the formula (9) is locally reconstructed at each interface, and specifically, the second-order approximation of the non-equilibrium density distribution function is calculated by Taylor series expansion as follows:

[0102]

[0103] In the formula (10), r refers to a position of a particle at a current time t; t refers to a current time; δ t refers to a flow field time step.

[0104] Density and velocity are locally reconstructed at the interface from the solution in the lattice Boltzmann method:

[0105]

[0106] D. In the correction step, the submerged boundary method with improved boundary conditions is used to correct the intermediate flow field velocity, and a new flow field velocity is obtained. Based on equation (8) of the correction step, the submerged boundary method with improved boundary conditions is used to correct the intermediate flow field velocity and a new flow field velocity u is obtained. n+1 The new flow field velocity u n+1 It is determined by the intermediate flow field velocity u * It is obtained by summing and correcting the flow field velocity δu:

[0107] u n+1 =u * +δu (12)

[0108] Wall velocity at Lagrange points The velocity u on the Euler grid n+1 (r j Interpolated from:

[0109]

[0110] In equation (13): h is the mesh size of the Eulerian mesh; N is the number of Lagrange points, M is the number of Eulerian points; r j It is the location of the Euler point. It is the location of the Lagrange point; D ij It is a continuous kernel function, given by the following equation (14):

[0111]

[0112]

[0113] Similarly, the correction velocity δu(r) of the Euler point j Corrected velocity at Lagrange points Interpolation yields:

[0114]

[0115] Based on formulas (12), (13), (14), and (16), the solution equation is obtained:

[0116]

[0117] make If X is the unknown quantity, then the solution equation (17) can be written as a matrix relation (18):

[0118] AX = B (18)

[0119] The matrix A and the matrix B in the formula (18) represent respectively:

[0120]

[0121]

[0122] The matrix relationship formula (18) is solved by using an improved matrix solving technique, and the specific steps are as follows: the matrix A is converted into a diagonal matrix A diag , the diagonal matrix A diag is given by the formula (21):

[0123] A diag = diag(d1 d2 … d N ) (21)

[0124]

[0125] According to the formula (18)-formula (22), the unknown quantity X is efficiently solved and obtained:

[0126]

[0127] In the formula (23), B = (B1 B2 … B N );

[0128] The modified flow field velocity δu is obtained based on the obtained unknown quantity X, and then the new flow field velocity u n+1

[0129] E, after reaching the simulation time, the whole simulation process is ended, and if the simulation time is not reached, it returns to step B.

[0130] Example analysis

[0131] In order to verify the correctness, effectiveness and rapidity of the improved immersed boundary method based flow-solid coupling method for simulating the autonomous swimming of flexible organisms provided in the embodiment, the related numerical simulation operation is carried out in this example to verify the correctness, effectiveness and rapidity of the improved immersed boundary method based flow-solid coupling method for simulating the autonomous swimming of flexible organisms. The steady ball disturbance and the example of the swimming of the yellow eel are selected to simulate and verify the feasibility of the method.

[0132] Wherein, the initial velocity is set to 0.1, the fluid density is 1, the non-uniform grid number of the whole flow field in the calculation of the steady ball disturbance is 171x139x139, the uniform grid number surrounding the ball is 60x60x60, the Lagrange point number is 1333, and δ t 0.015, the drag coefficient under different Reynolds numbers is compared, and the results are shown in Table 1:

[0133] Table 1 Drag coefficient comparison of the circumscribed sphere disturbance at Re = 100, 250, 300

[0134]

[0135] From Table 1 we can see that the results of the two methods are very close, which shows that the method obtains consistent drag distribution with the traditional IB-LBFS method, and proves the effectiveness and accuracy of the new method.

[0136] Further, we select Re = 100, change the uniform grid number surrounding the sphere, which are 40x40x40, 50x50x50, 60x60x60 respectively, to verify the efficiency of IB-LBFS based on the improved matrix solving technique. The CPU time consumed by each iteration step of the two solvers is shown in Table 2 as follows:

[0137] Table 2 Comparison of the calculation efficiency of the original IB-LBFS method and the improved IB-LBFS method under different grid numbers

[0138]

[0139] It can be seen that the current IB-LBFS only spends two percent of the original CPU time in each iteration step, thus greatly saving the calculation cost and improving the calculation efficiency.

[0140] As shown in the geometry of the yellow eel shown in Figure 4 , we give a simple model of the dorsal-ventral motion Z f of the bat-like pectoral fin of the yellow eel as follows: where: r * and θ * represent the positions of each point of the body of the yellow eel, which are read from the existing modeling file.

[0141] wherein the non-uniform grid number of the entire flow field is 254x239x239, and the uniform grid number surrounding the yellow eel is 150x150x150, wherein the reference length L ref = 0.1, the reference frequency F ref = 1, the reference velocity U ref = F ref *L ref = 0.1, the reference amplitude A ref = 0.1*L ref = 0.01, and the Reynolds number force coefficient reference phase reference wave number κ = 2. Then we change the wave number to 1, 2, 3, and 4 respectively, assign the corresponding parameters, and the time step δ tTaking 0.0005, the result graph as shown in Figure 5

[0142] Figure 5 The thrust coefficient comparison chart under different wave numbers is shown, from which it can be seen that the greater the wave number, the smaller the peak of the drag coefficient, and the more stable the drag coefficient, which also shows that the more flexible the wobbler, the more stable the wobbler.

[0143] In summary, through the test and analysis of the above examples and calculation examples, it can be fully illustrated that the new method can correctly simulate the complex flow field of incompressible flow on non-uniform grid, and the new method requires less calculation time than the original IB-LBFS, improves the calculation efficiency, and saves a large amount of calculation cost, which is the advantage of the new method.

[0144] The flow-solid coupling method for simulating the autonomous swimming of flexible organisms provided by the application improves the original boundary condition forced immersed boundary method by using the improved boundary condition forced immersed boundary method to correct the intermediate flow field velocity on the basis of the original immersed boundary-lattice Boltzmann flux solver method, improves the original boundary condition forced immersed boundary method through matrix construction and solving strategy, effectively simulates incompressible flow with a moving boundary, and saves a large amount of calculation cost.

[0145] The above examples only illustrate the technical idea of the application, and cannot limit the protection scope of the application, any modification made according to the technical idea of the application on the basis of the technical scheme falls within the protection scope of the application; the technologies not involved in the application can be realized by the existing technologies.​

Claims

1. A fluid-structure coupling method for simulating autonomous swimming of a flexible organism based on an improved immersed boundary method, characterized in that: The method steps are as follows: A. specifying initial parameters at each interface; B. using distribution technology to divide the control equation of the immersed boundary-lattice Boltzmann flux solver into a prediction step and a correction step; C. in the prediction step, using the lattice Boltzmann flux solver to solve the intermediate flow field velocity; D. in the correction step, using the improved boundary condition forced immersed boundary method to correct the intermediate flow field velocity to obtain a new flow field velocity; E. after reaching the simulation time, the entire simulation process is ended, and if the simulation time is not reached, return to step B; The control equation of the immersed boundary-lattice Boltzmann flux solver in step B is: In equation (1): p is fluid density; u is fluid velocity; p is pressure; m is dynamic viscosity; T is temperature; f is the restoring force determined by the immersed boundary method to account for the boundary effects of the immersed object; According to the multi-scale Chapman-Enskog analysis, the density distribution function f α The relationship between the flux ρu in equation (1) and the density distribution function f In formula (2): I is the unit tensor; f α is the density distribution function along the a direction; is the equilibrium density distribution function along the a direction; is the non-equilibrium density distribution function along the a direction; τ is the single relaxation parameter; δ t is the flow to time step; e α is the particle velocity along the a direction; β and γ are two directions in e α is the pressure; Using distribution technology, the formula (1) and formula (2) involved in the control equation of the immersed boundary-lattice Boltzmann flux solver are divided into a prediction step and a correction step, wherein the prediction step is given by the following formula (7), and the correction step is given by the following formula (8): In formula (7) and formula (8): At is a time step; u * is an intermediate flow field velocity; p n is a density of a current time step; p n+1 is a density of a next time step; u n is a flow field velocity of a current time step; u n+1 is a flow field velocity of a next time step; The equation (7) of the prediction step is solved with the lattice Boltzmann flux solver to obtain the density p of the next time step n+1 and the intermediate flow field velocity u * ; in a control volume Ω i , the equation (7) of the prediction step is discretized by finite volumes as follows: In formula (9): W * = (p n+1 , p n+1 u * ); AV i is the volume of the control volume Ω i ; G k is the flux; ΔS k is the area of the kth interface surrounding the control volume Ω i ; n is the unit outward normal vector of the kth interface. Based on the modified equation (8), the intermediate flow field velocity is modified using the submerged boundary method with improved boundary conditions, and a new flow field velocity u is obtained. n+1 The new flow field velocity u n+1 It is determined by the intermediate flow field velocity u * It is obtained by summing and correcting the flow field velocity δu: u n+1 = u * + δu (12) Wall velocities at Lagrangian points by the velocity u on the Euler grid n+1 (r j ) interpolated In formula (13): h is the grid size of the Euler grid; N is the number of Lagrange points, M is the number of Euler points; r j is the position of the Euler point, is the position of the Lagrange point; D ij is a continuous kernel function, given by the following formula (14): Similarly, the modified velocity δu(r j ) of the Euler point is interpolated from the modified velocities δu(r ) at the Lagrange points. Based on formula (12), formula (13), formula (14) and formula (16), the solving equation is obtained as follows: Let For the unknown quantity X, solving equation (17) can be written as matrix relation (18): AX=B (18) The matrix A and the matrix B in formula (18) are respectively represented as:

2. The flow-solid coupling method for simulating autonomous swimming of a flexible living body based on the improved immersed boundary method according to claim 1, characterized in that: The initial parameters in step A include a flow direction time step, a fluid density, a fluid velocity, a pressure, a dynamic viscosity and a Reynolds number.

3. The flow-solid coupling method for simulating autonomous swimming of a flexible living body based on the improved immersed boundary method according to claim 1, characterized in that: In the three-dimensional incompressible flow field, the particle velocity e α is given by the following equation (3) using the D3Q15 lattice velocity model. Equilibrium density distribution function along the a direction is given by the following equation (4): In formula (4): w α is a coefficient of the equilibrium density distribution function, c s is the sound speed, and w α and c s have values that depend on the lattice velocity model; the single relaxation parameter τ is related to the dynamic viscosity μ of the fluid, given by the following formula (5): Non-equilibrium density distribution function along the a direction is given by the following equation (6):

4. The flow-solid coupling method for simulating autonomous swimming of a flexible living body based on the improved immersed boundary method according to claim 1, characterized in that: At each interface, the flux G in equation (9) k Local reconstruction is performed, specifically as follows: the non-equilibrium density distribution function is expanded using Taylor series. The second-order approximation is calculated as follows: In formula (10): r refers to the position of the particle at the current time t; t refers to the current time; δ t refers to the flow field time step; The density and the velocity are locally reconstructed at the interface by the solution in the lattice Boltzmann method:

5. The flow-solid coupling method for simulating autonomous swimming of a flexible living body based on the improved immersed boundary method according to claim 1, characterized in that: The matrix relational expression (18) is solved by using an improved matrix solving technique. The specific steps are as follows: converting the matrix A into a diagonal matrix A diag The diagonal matrix A diag is given by expression (21). A diag = diag(d1 d2 …d N ) (21) According to formula (18)-formula (22), the unknown quantity X is efficiently solved as follows: In formula (23): B = (B1B2...Bn). N ).

Images