A Parachute Fluid-Structure Interaction Calculation Method Based on Immersion Boundary and Finite Element Method

Abstract

This invention provides a parachute fluid-structure interaction (FSI) calculation method based on submerged boundary and finite element methods, comprising: Step 1, establishing a parachute FSI model; Step 2, describing the parachute system using the Neo-Hookean model and the Saint Venant–Kirchhoff model for hyperelastic nonlinear materials, respectively; Step 3, achieving information transfer between Eulerian-Lagrange grids through Dirac trigonometric functions; Step 4, establishing geometric models of the canopy and lines and generating the corresponding grids required for finite element calculations; Step 5, setting the computational domain; and Step 6, analyzing the deformation of various parachutes and the changes in the surrounding flow field. This invention, based on the submerged boundary method, eliminates the need for dynamically generated body-fitted grids, offering significant advantages over traditional body-fitted grid solutions in simulating the large deformations involved in parachute inflation.

Description

Technical Field

[0001] This invention belongs to the fields of computational fluid dynamics and fluid-structure interaction, and particularly relates to a parachute fluid-structure interaction calculation method based on immersion boundary and finite element method. Background Technology

[0002] The parachute inflation process involves complex fluid-structure interaction (FSI) phenomena, which are highly nonlinear and unsteady, making numerical simulation of this process very challenging.

[0003] Over the past few decades, numerical methods for solving parachute inflation FSI problems have seen tremendous development. One of the most popular methods is the Deformation Domain / Steady Space-Time (DSD / SST) method, which can obtain the flow field and structural response during parachute inflation. LS-DYNA, a commercial explicit dynamics solver, is also widely used in engineering practice for FSI simulation of parachute inflation. The Immersed Boundary (IB) method, because it does not require dynamic generation of body-fitted meshes, is suitable for solving problems involving large structural deformations or displacements, and has gradually become a popular method for parachute inflation FSI simulation in recent years. Kim and Peskin used the IB method to study the inflation process of two-dimensional and three-dimensional semi-open parachutes. Liu et al. proposed an IB-lattice Boltzmann (LB) based FSI solver to study the inflation process of different types of parachute systems.

[0004] IBAMR is a distributed, memory-parallel implementation of the IB method, supporting adaptive mesh refinement on Cartesian meshes. IBAMR's core functionality relies on several high-quality open-source libraries, including: SAMRAI, an application infrastructure for adaptive structured mesh refinement; PETSc, a portable, scalable toolkit for scientific computing; libMesh, a finite element library written in C++; and hypre, a high-performance preprocessing library based on parallel multimeshing methods. Currently, IBAMR is widely used to simulate various FSI problems involving large deformations. Summary of the Invention

[0005] Objective of the Invention: The technical problem to be solved by this invention is to address the shortcomings of existing technologies by providing a parachute fluid-structure interaction (FSI) calculation method based on immersion boundary and finite element methods to accurately simulate the inflation process of various parachutes. Numerical experiments verify that the computational framework of this invention is an accurate and reliable FSI simulation platform for the parachute inflation process, and it is expected to be further extended to the simulation of the parachute system inflation process in various aerospace missions.

[0006] This invention utilizes the open-source IBAMR library to establish a partitioned fluid-structure interaction (FSI) model of a parachute based on the submerged boundary / finite element method. Fluid motion is solved using the diffusion interface IB method, with adaptive mesh refinement employed to reduce computational costs while maintaining the required mesh resolution near the object surface. The dynamic response of the parachute system is solved using the finite element method. The parachute canopy and lines are modeled as hyperelastic materials. Rigid constraints are applied to the fixed ends of the parachute system by introducing tethering forces. The established FSI model is used to simulate various parachute inflation processes. By analyzing and comparing the deformation of the parachute canopy and the changes in the surrounding flow field, the reliability of the established FSI computational framework for the parachute is verified. The method of this invention specifically includes the following steps:

[0007] Step 1: Establish a parachute fluid-structure interaction model;

[0008] Step 2: The parachute system is described using the Neo-Hookean model and the Saint Venant–Kirchhoff model for hyperelastic nonlinear materials, respectively.

[0009] Step 3: Information transfer between Euler-Lagrange grids is achieved through Dirac trigonometric functions;

[0010] Step 4: Establish the geometric model of the canopy and parachute lines and generate the corresponding mesh required for finite element calculation;

[0011] Step 5, set the computation domain;

[0012] Step 6: Analyze the deformation of various parachutes and the changes in the surrounding flow field.

[0013] Step 1 includes: Describing the entire computational domain as Ω = Ω f ∪Ω s , where Ω f and Ω s Let X represent the fluid domain and the solid domain respectively, x represent the Eulerian coordinates, and χ(X,t) represent the physical position of the Lagrange point X at time t. The fluid-structure interaction governing equation of the parachute system is expressed as:

[0014]

[0015]

[0016]

[0017]

[0018] in, It is the gradient operator. dX is the differential operator at the Lagrange point X, dX is the Lagrange point differential, and ρ is the density. It is the partial derivative of the physical location of the Lagrange point with respect to time. Let p(x,t) be the mass derivative, p(x,t) be the pressure, u(x,t) be the velocity field described by Euler, f(x,t) be the elastic force described by Euler, μ be the viscosity coefficient of the fluid, δ(x) be the Dirac trigonometric function, N(x) be the unit outward normal vector of the structural surface, dA be the area differential, and P be the mass derivative. s (X,t) is the first Piola-Kirchhoff stress used to describe the structural response.

[0019] Step 2 includes:

[0020] Step 2-1, for the Neo-Hookean model, the strain energy function Ψ after volume stability correction is:

[0021]

[0022] Where J is the Jacobian determinant of the deformation gradient, and I1 is the first invariant of the right Cauchy-Green deformation tensor. ν is the shear modulus, E is Young's modulus, ν is Poisson's ratio, and κ is... stab For numerical bulk modulus, the expression is: ν stab The numerical Poisson's ratio;

[0023] The corresponding first Piola-Kirchhoff stress P s Represented as:

[0024]

[0025] Where T represents the matrix transpose and F is the transformed gradient;

[0026] Step 2-2, for the Saint Venant–Kirchhoff model, the strain energy function Ψ after volume stability correction is:

[0027]

[0028] in, Represents the Green-Lagrange strain tensor, C = F T F is the right Cauchy-Green deformable tensor, and I is the unit tensor. and Represents Lamé's constant;

[0029] The corresponding first Piola-Kirchhoff stress representation P s for:

[0030] P s=F(λtr(E)I+2μE) (8)

[0031] Where tr is the trace of the matrix;

[0032] The fixed-end constraint at the paracord junction is satisfied by applying a tethering force, which is approximately represented by the Lagrange multiplier F(X,t):

[0033]

[0034] Where κ and η are the stiffness and damping penalty parameters, respectively.

[0035] Step 3 includes: For an Eulerian mesh, the spacing is defined as (Δx1, Δx2), where Δx1 is the mesh size in the x-direction and Δx2 is the mesh size in the y-direction; for a given mesh element (i, j), where i is the mesh number in the x-direction and j is the mesh number in the y-direction, the pressure is approximately at the center of the element and denoted as p. i,j The velocity is approximately at the center of the element edge and is denoted as . exist Place and exist At this point, u1 and u2 are the velocity components at the grid point in the x-direction and the y-direction, respectively; Indicates the grid number The coordinates of the grid point at that location, This represents the velocity component in the x-direction;

[0036] The incompressible Navier-Stokes equations are discretized using a second-order staggered grid finite difference scheme. The structure is superimposed on an Eulerian grid and discretized by a series of Lagrange points. The governing equations of the structure are discretized using the nodal finite element method. Information transfer between the Eulerian and Lagrange grids is achieved through Dirac trigonometric functions.

[0037] Step 4 includes: using the open-source mesh generation software Gmsh to establish the geometric models of the canopy and parachute lines and generate the corresponding meshes required for finite element calculations; using shell elements to mesh the canopy and cable elements to mesh the parachute lines.

[0038] Step 5 includes: setting the computational domain to an 8m×4m×4m cuboid, setting the left side of the computational domain as a Dirichlet velocity inlet boundary condition, and setting the incoming flow velocity U. ∞ =1 m / s, the right side is set as the Neumann velocity exit boundary condition, and the other four sides are set as free slip boundary conditions; the Reynolds number Re is set to Where the air density ρ = 1 kg / m³ 3 μ fLet L be the viscosity coefficient of air and L be the characteristic length. The computational domain is discretized by an adaptively refined Cartesian mesh with 3 layers. The ratio of the mesh scale between adjacent layers is 4. The number of elements in the coarsest layer is 16, and the mesh scale Δx in the densest layer is Δx = 1.5 × 10⁻⁶. -2 m and the Lagrange grid scale ratio with the structure is 1:1, and the computation time Δt is Δt = 1 × 10 -4 s.

[0039] Step 6 includes: using the established fluid-structure interaction model to simulate the inflation process of square, cross, and double cross parachutes, and using the open-source post-processing software VisIt to analyze the deformation of various parachutes and the changes in the surrounding flow field, as well as the drag coefficient C of the parachutes. D Represented as:

[0040]

[0041] Where F D A and A' represent the drag force on the canopy and the nominal area of ​​the canopy, respectively.

[0042] Beneficial effects: This invention can accurately and reliably simulate the inflation process of various parachutes. The framework has the following advantages: (1) Based on the submerged boundary method, it does not require dynamic generation of body-fitted meshes. Compared with the traditional body-fitted mesh solution method, it has a significant advantage in simulating the large deformations involved in the parachute inflation process; (2) For parachute fluid-structure interaction simulation, unlike the current reliance on commercial software or laboratory-developed software, this invention establishes a corresponding calculation framework based on the open-source IBAMR library for the first time. Moreover, the entire pre-processing and post-processing process also relies on open-source software, providing a new approach for parachute fluid-structure interaction simulation. Attached Figure Description

[0043] The present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments, and the advantages of the present invention in the above and / or other aspects will become clearer.

[0044] Figure 1 This is a two-dimensional schematic diagram of an Eulerian-Lagrange grid layout.

[0045] Figure 2a This is a geometric model and finite element mesh diagram of a square umbrella.

[0046] Figure 2b This is a geometric model and finite element mesh diagram of the cross umbrella.

[0047] Figure 2c This is a geometric model and finite element mesh diagram of a double cross umbrella.

[0048] Figure 3 This is a schematic diagram of the computational domain and boundary condition settings.

[0049] Figure 4a This is a schematic diagram showing the time history results of the drag coefficient of a square umbrella.

[0050] Figure 4b This is a schematic diagram showing the time history of the displacement of the center point of the canopy of a square umbrella.

[0051] Figure 5a A schematic diagram of the deformation of a square umbrella at t = 0.25s.

[0052] Figure 5b A schematic diagram of the deformation of a square umbrella at t = 0.5s.

[0053] Figure 5c A schematic diagram of the deformation of a square umbrella at t = 0.75s.

[0054] Figure 5d A schematic diagram of the deformation of a square umbrella at t=1s.

[0055] Figure 5e A schematic diagram of the deformation of a square umbrella at t = 1.25s.

[0056] Figure 5f A schematic diagram of the deformation of a square umbrella at t = 1.5s.

[0057] Figure 6a This is a diagram showing the transient vorticity isosurface (Q=40) of a square umbrella at t=2s and the effect of dynamic adaptive mesh refinement.

[0058] Figure 6b This is a diagram showing the transient vorticity isosurface (Q=40) of a square umbrella at t=3s and the effect of dynamic adaptive mesh refinement.

[0059] Figure 6c This is a diagram showing the transient vorticity isosurface (Q=40) of a square umbrella at t=4s and the effect of dynamic adaptive mesh refinement.

[0060] Figure 6d This is a diagram showing the transient vorticity isosurface (Q=40) of a square umbrella at t=5s and the effect of dynamic adaptive mesh refinement.

[0061] Figure 7a This is a schematic diagram showing the time history of the drag coefficient of the cross-shaped umbrella.

[0062] Figure 7b This is a schematic diagram showing the time history of the displacement of the center point of the canopy of the cross-shaped umbrella.

[0063] Figure 8a This is a schematic diagram of the deformation of the cross-shaped umbrella at t = 0.2s.

[0064] Figure 8b This is a schematic diagram of the deformation of the cross-shaped umbrella at t = 0.4s.

[0065] Figure 8c This is a schematic diagram of the deformation of the cross-shaped umbrella at t = 0.6s.

[0066] Figure 8d This is a schematic diagram of the deformation of the cross-shaped umbrella at t = 0.8s.

[0067] Figure 8e This is a schematic diagram of the deformation of the cross-shaped umbrella at t=1s.

[0068] Figure 8f This is a schematic diagram of the deformation of the cross-shaped umbrella at t = 1.2s.

[0069] Figure 9a This is a schematic diagram of the deformation of the double-cross umbrella at t=0.1s.

[0070] Figure 9b This is a schematic diagram of the deformation of the double-cross umbrella at t=0.2s.

[0071] Figure 9c This is a schematic diagram of the deformation of the double-cross umbrella at t=0.3s.

[0072] Figure 9d This is a schematic diagram of the deformation of the double-cross umbrella at t=0.4s.

[0073] Figure 9e This is a schematic diagram of the deformation of the double-cross umbrella at t=0.5s.

[0074] Figure 9f This is a schematic diagram of the deformation of the double-cross umbrella at t=0.6s.

[0075] Figure 10a This is a slice view of the velocity field around the canopy of a double-cross umbrella at t=0.4s.

[0076] Figure 10b This is a slice view of the velocity field around the canopy of the double-cross umbrella at t=0.8s.

[0077] Figure 10c This is a slice view of the velocity field around the canopy of a double-cross umbrella at t=1.3s.

[0078] Figure 10d This is a slice view of the velocity field around the canopy of a double-cross umbrella at t=2s.

[0079] Figure 10e This is a slice view of the velocity field around the canopy of a double-cross umbrella at t=3s.

[0080] Figure 10f This is a slice view of the velocity field around the canopy of a double-cross umbrella at t=4s. Detailed Implementation

[0081] This invention provides a parachute fluid-structure interaction calculation method based on immersion boundary and finite element method, including the following steps:

[0082] Step 1: Denote the entire computational domain as Ω = Ω f ∪Ω s , where Ω f and Ω s Let represent the fluid domain and the solid domain respectively, x represent the Eulerian coordinates, and χ(X,t) represent the physical position of the Lagrange point X at time t. The fluid-structure interaction governing equation of the parachute system is expressed as:

[0083]

[0084]

[0085]

[0086]

[0087] Where ρ is density, Let p(x,t) be the mass derivative, p(x,t) be the pressure, u(x,t) be the velocity field described by Euler, f(x,t) be the elastic force described by Euler, μ be the viscosity coefficient of the fluid, δ(x) be the Dirac trigonometric function, N(x) be the unit outward normal vector of the structural surface, dA be the area differential, and P be the mass derivative. s (X,t) is the first Piola-Kirchhoff stress used to describe the structural response.

[0088] Step 2: The parachute system is described using the Neo-Hookean model and the Saint Venant–Kirchhoff model for hyperelastic nonlinear materials, respectively.

[0089] For the Neo-Hookean model, the strain energy function after volume stability correction is:

[0090]

[0091] Where J is the Jacobian determinant of the deformation gradient, and I1 is the first invariant of the right Cauchy-Green deformation tensor. ν is the shear modulus, E is Young's modulus, ν is Poisson's ratio, and κ is... stab For numerical bulk modulus, the expression is: ν stab This is the numerical Poisson's ratio.

[0092] The corresponding first Piola-Kirchhoff stress is expressed as:

[0093]

[0094] For the Saint Venant–Kirchhoff model, the strain energy function after volume stability correction is:

[0095]

[0096] in, Represents the Green-Lagrange strain tensor, C = F T F is the right Cauchy-Green deformable tensor, and I is the unit tensor. and This represents the Lamé constant.

[0097] The corresponding first Piola-Kirchhoff stress is expressed as:

[0098] P s =F(λtr(E)I+2μE) (8)

[0099] The fixed-end constraint at the parachute line junction is satisfied by applying a tethering force. This tethering force is approximately represented by the Lagrange multiplier F(X,t):

[0100]

[0101] Where κ and η are the stiffness and damping penalty parameters, respectively.

[0102] Step 3: A two-dimensional schematic diagram of the Eulerian-Lagrange grid layout in spatial discretization is shown below. Figure 1 As shown. For an Eulerian mesh, the spacing is defined as (Δx1, Δx2); for a given mesh element (i, j), the pressure is approximately at the element center and denoted as p. i,j The velocity is approximately at the center of the element edge and is denoted as . exist Place and exist The incompressible Navier-Stokes equations are discretized using a second-order staggered-grid finite-difference scheme. The structure is superimposed on an Eulerian grid and discretized by a series of Lagrange points. The governing equations of the structure are discretized using the nodal finite element method. Information transfer between the Eulerian and Lagrange grids is achieved through Dirac trigonometric functions.

[0103] Step 4: Using the open-source mesh generation software Gmsh, the geometric models of the umbrella canopy and lines are created, and the corresponding meshes required for finite element calculations are generated. The geometric models and finite element meshes for square umbrellas, cross umbrellas, and double cross umbrellas are shown below. Figure 2a , Figure 2b , Figure 2cAs shown, the canopy length of the square umbrella is L = 0.8128m, and the length of the parachute ropes is l = 1.27m. The canopy length of the cross umbrella is the same as that of the square umbrella, and the distance from the point where the parachute ropes intersect to the surface of the canopy is the same as that of the square umbrella. Shell elements are used to divide the canopy, and cable elements are used to divide the parachute ropes.

[0104] Step 5: As Figure 3 As shown, the computational domain is set as an 8m × 4m × 4m cuboid. The left side of the computational domain is set as a Dirichlet velocity inlet boundary condition, with an incoming flow velocity of U. ∞ =1 m / s, the right side is set as the Neumann velocity exit boundary condition, and the other four sides are set as free slip boundary conditions, where x1, x2, and x3 represent the x, y, and z directions respectively, and u1, u2, and u3 are the velocity components in the x, y, and z directions respectively. The Reynolds number is set to... Where the air density ρ = 1 kg / m³ 3 μ f Let be the viscosity coefficient of air. The computational domain is discretized by an adaptively refined Cartesian mesh with 3 layers. The ratio of the mesh scale between adjacent layers is 4. The number of elements in the coarsest layer is 16, and the mesh scale in the densest layer is Δx = 1.5 × 10⁻⁶. -2 m and the Lagrange grid scale ratio with the structure is 1:1, and the computation time is Δt = 1 × 10 -4 s.

[0105] Step 6: The established fluid-structure interaction model is used to simulate the inflation process of square, cross, and double cross parachutes, and the deformation of various parachutes and the changes in the surrounding flow field are analyzed using the open-source post-processing software VisIt. The drag coefficient of the parachute is expressed as:

[0106]

[0107] Where F D A and A' represent the drag force on the canopy and the nominal area of ​​the canopy, respectively.

[0108] For square umbrellas, Figure 4a , Figure 4b The drag coefficient C corresponding to different material models is shown. D The result of the change in displacement (dX) of the center point of the parachute canopy with time. Figure 5a , Figure 5b , Figure 5c , Figure 5d , Figure 5e , Figure 5f The umbrella canopy deformation at a typical moment is shown. Figure 6a , Figure 6b , Figure 6c , Figure 6dThe transient vorticity isosurface at typical moments and the effect of dynamic adaptive mesh refinement are shown. The calculated steady-state drag coefficient of the square umbrella is 1.1, which is close to the reference result.

[0109] For the cross umbrella, Figure 7a , Figure 7b The drag coefficient C corresponding to different material models is shown. D The result of the change in displacement (dX) of the center point of the parachute canopy with time. Figure 8a , Figure 8b , Figure 8c , Figure 8d , Figure 8e , Figure 8f The deformation of the canopy at a typical moment is shown. The calculated steady-state drag coefficient of the cross-shaped parachute is 0.8, which is close to the reference result.

[0110] For double cross umbrellas, Figure 9a , Figure 9b , Figure 9c , Figure 9d , Figure 9e , Figure 9f The umbrella canopy deformation at a typical moment is shown. Figure 10a , Figure 10b , Figure 10c , Figure 10d , Figure 10e , Figure 10f The image shows a sliced ​​view of the velocity field (U_magnitude) around the canopy at a typical moment. The calculated steady-state drag coefficient of the double-cross umbrella is 1.0, which is 1.25 times that of the single-cross umbrella. This value is also close to the reference result.

[0111] These results demonstrate that the parachute fluid-structure interaction calculation framework based on the immersion boundary / finite element method proposed in this invention can accurately and reliably simulate the inflation process of various parachutes, and therefore is expected to be further extended to the simulation of the inflation process of parachute systems in various aerospace missions.

[0112] In its specific implementation, this application provides a computer storage medium and a corresponding data processing unit. The computer storage medium is capable of storing a computer program, which, when executed by the data processing unit, can run the invention's content regarding the parachute fluid-structure interaction calculation method based on immersion boundary and finite element methods, as well as some or all of the steps in various embodiments. The storage medium can be a magnetic disk, optical disk, read-only memory (ROM), or random access memory (RAM), etc.

[0113] Those skilled in the art will clearly understand that the technical solutions in the embodiments of the present invention can be implemented using computer programs and their corresponding general-purpose hardware platforms. Based on this understanding, the technical solutions in the embodiments of the present invention, or the parts that contribute to the prior art, can be embodied in the form of computer programs, i.e., software products. These computer program software products can be stored in a storage medium and include several instructions to cause a device containing a data processing unit (which may be a personal computer, server, microcontroller, MUU, or network device, etc.) to execute the methods described in various embodiments or certain parts of the embodiments of the present invention.

[0114] This invention provides a parachute fluid-structure interaction calculation method based on immersion boundary and finite element method. Many methods and approaches exist for implementing this technical solution; the above description is merely a preferred embodiment of this invention. It should be noted that those skilled in the art can make various improvements and modifications without departing from the principles of this invention, and these improvements and modifications should also be considered within the scope of protection of this invention. All components not explicitly stated in this embodiment can be implemented using existing technologies.

Claims

1. A parachute fluid-structure interaction calculation method based on immersion boundary and finite element method, characterized in that, Includes the following steps: Step 1: Establish a parachute fluid-structure interaction model; Step 2: The parachute system is described using the Neo-Hookean model and the Saint Venant–Kirchhoff model for hyperelastic nonlinear materials, respectively. Step 3: Information transfer between Euler-Lagrange grids is achieved through Dirac trigonometric functions; Step 4: Establish the geometric model of the canopy and parachute lines and generate the corresponding mesh required for finite element calculation; Step 5, set the computation domain; Step 6: Analyze the deformation of various parachutes and the changes in the surrounding flow field; Step 1 includes: denoting the entire computational domain as ,in and Representing the fluid domain and the solid domain respectively. Represents Euler coordinates, express Lagrange point of time The physical location of the parachute system, and the fluid-structure interaction control equations of the parachute system, are expressed as: （1）， （2）， （3） ， （4）， in, It is the gradient operator. It is the differential operator of the Lagrange point X. It is the Lagrange point differential. For density, It is the partial derivative of the physical location of the Lagrange point with respect to time. For the matter derivative, For pressure, The velocity field described by Euler. The elastic force described by Euler, Let be the viscosity coefficient of the fluid. For Dirac trigonometric functions, Let be the unit outward normal vector of the structural surface. For the area differential, The first Piola-Kirchhoff stress is used to describe the structural response; Step 2 includes: Step 2-1, for the Neo-Hookean model, the strain energy function after volume stability correction. for: （5）， in, Let Jacobian determinant be the deformation gradient. Let be the first invariant of the right Cauchy-Green deformation tensor. Shear modulus For Young's modulus, Poisson's ratio, For numerical bulk modulus, the expression is: , The numerical Poisson's ratio; The corresponding first Piola-Kirchhoff stress Represented as: （6）， Where T represents the matrix transpose and F is the transformed gradient; Step 2-2, for the Saint Venant–Kirchhoff model, the strain energy function after volume stability correction. for: （7）， in, Represents the Green-Lagrange strain tensor. For the right Cauchy-Green deformation tensor, For unit tensors, and Represents Lamé's constant; The corresponding first Piola-Kirchhoff stress representation for: （8）， in It is the trace of the matrix; The fixed-end constraint at the paracord junction is satisfied by applying a tethering force, which is determined by the Lagrange multiplier. Approximate representation: （9）， in, and These are the stiffness and damping penalty parameters, respectively.

2. The method according to claim 1, characterized in that, Step 3 includes: For Euler meshes, the spacing is defined as ( ), It is the grid size in the x-direction. It is the mesh size in the y-direction; for a given mesh cell , It is the grid number in the x-direction. It is the grid number in the y-direction, and the pressure is approximately at the cell center and denoted as . The velocity is approximately at the center of the element edge and is denoted as . exist Place and exist place, These represent the velocity components at grid points in the x-direction and the velocity components at grid points in the y-direction, respectively. Indicates the grid number The coordinates of the grid point at that location, This represents the velocity component in the x-direction; The incompressible Navier-Stokes equations are discretized using a second-order staggered grid finite difference scheme. The structure is superimposed on an Eulerian grid and discretized by a series of Lagrange points. The governing equations of the structure are discretized using the nodal finite element method. Information transfer between the Eulerian and Lagrange grids is achieved through Dirac trigonometric functions.

3. The method according to claim 2, characterized in that, Step 4 includes: using the open-source mesh generation software Gmsh to establish the geometric models of the canopy and parachute lines and generate the corresponding meshes required for finite element calculations; using shell elements to mesh the canopy and cable elements to mesh the parachute lines.

4. The method according to claim 3, characterized in that, Step 5 includes: setting the computation domain to... The cuboid is used, and the left side of the computational domain is set as the Dirichlet velocity inlet boundary condition, with the incoming flow velocity... The right side is set as the Neumann velocity exit boundary condition, and the other four sides are set as free slip boundary conditions; Reynolds number Set as air density , The viscosity coefficient of air. The feature length is specified; the computational domain is discretized by an adaptively refined Cartesian grid with 3 layers, a grid scale ratio of 4 between adjacent layers, 16 cells in the coarsest layer, and a grid scale of [missing information - likely a specific value]. for Furthermore, the Lagrangian mesh scale ratio with the structure is 1:1, and the computation time is... for .

5. The method according to claim 4, characterized in that, Step 6 includes: using the established fluid-structure interaction model to simulate the inflation process of square, cross, and double cross parachutes, and using the open-source post-processing software VisIt to analyze the deformation of various parachutes and the changes in the surrounding flow field, as well as the drag coefficient of the parachutes. Represented as: （10）， in and These represent the drag force experienced by the canopy and the nominal area of ​​the canopy, respectively.

Images