A numerical simulation method for load out-of-cabin process based on dynamic nested grid technology

Abstract

The present application relates to the technical field of air separation, in particular to a load ejection process numerical simulation method based on dynamic nested grid technology, which first establishes a carrier and load geometric model, and then establishes a multi-layer dynamic nested grid model of the separation system, wherein the traction effect of the parachute is simulated based on a multi-particle spring damping model, and the contact and collision process between the load and the adapter during the ejection process is simulated by using a penalty function method; a three-dimensional compressible unsteady air flow model and a multi-rigid body six-degree-of-freedom motion model are established, the load ejection process numerical simulation calculation is carried out, and the motion change of the carrier and the load and the safety of the separation process are monitored based on the Euler coordinate system. The present application simulates the continuous transient process of the load ejection from the tail of the carrier under the action of the traction parachute by using the dynamic nested grid technology, so as to ensure the calculation precision of the ejection process.

Description

Technical Field

[0001] This invention relates to the field of aerial separation technology, and in particular to a numerical simulation method for the payload exit process based on dynamic nested mesh technology. Background Technology

[0002] The payload tail-end ejection system is a type of in-flight multi-body separation system, specifically including air-launched rockets and vehicle airdrops. Ejection methods primarily include parachute-assisted ejection and gravity-based ejection. The separation process is a cross-coupling of aerodynamics and motion, exhibiting unsteady and nonlinear aerodynamic characteristics that pose a significant threat to flight safety. Although the separation process is brief, it determines the flight attitude and velocity of the launch vehicle and payload after separation, significantly impacting the success of subsequent operations. Exploring the aerodynamic and kinematic characteristics during multi-body separation through theoretical research and numerical simulation analysis is crucial for the development of separation systems.

[0003] Current research on numerical simulation of the tail-end exit process of the payload mostly adopts the Dynamic Modeling and Simulation (DMS) method. This method selects some representative positions and attitudes during the separation process to perform steady-state aerodynamic force calculations, fits the calculated discrete aerodynamic data into a continuous aerodynamic data curve through correction coefficients, and then inputs it into the dynamic model for rigid body motion simulation. However, the traditional DMS method still has the following shortcomings:

[0004] I. Aerodynamic data in the DMS method generally originates from: firstly, selecting different relative positions and attitudes of the load in the perturbation flow field of the launch vehicle as calculation points; then, calculating the aerodynamic coefficients at these points using CFD; and finally, fitting these discrete aerodynamic coefficient values ​​into a continuous aerodynamic coefficient curve and database. For example, the aerodynamic coefficients of the separation system when the load exits the cabin at 0m, 3m, and 6m are calculated using CFD, while the aerodynamic coefficients at other states such as 1m and 5m exiting the cabin are obtained through curve fitting. If too few calculation points are selected, the overall fidelity of the aerodynamic database will be very poor; if too many calculation points are selected, the complexity of manual operation and the demand for computing resources will increase dramatically. Furthermore, when there is significant relative motion between multiple bodies, it is often necessary to select a large number of calculation points.

[0005] Second, the aerodynamic data in the DMS method are mostly steady-state aerodynamic data obtained through steady-state calculations, while multibody system separation in mid-air is actually a dynamic and unsteady process. Steady-state calculations cannot take into account the added mass of a three-dimensional object in motion and the unsteady aerodynamic disturbances induced by the relative motion of components during multibody separation, and their values ​​may not necessarily reflect the actual situation.

[0006] In summary, although the DMS method has lower computational resource requirements, it does not take into account the unsteady effects caused by multibody motion during the separation process. Furthermore, the aerodynamic database of the DMS method is fitted from discrete steady-state data, while the actual process of aerial separation is a continuous transient process. Therefore, the computational accuracy of this method is insufficient for quantitative research on the overall separation process. Summary of the Invention

[0007] In view of this, the purpose of this invention is to propose a numerical simulation method for the load exit process based on dynamic nested mesh technology, in order to solve the problems that the DMS method does not take into account the unsteady effects caused by multibody motion during the separation process, and that the aerodynamic database of the DMS method is fitted by discrete steady-state data, while the actual process of separation in the air is a continuous transient process, which makes it difficult to meet the overall quantitative research of the separation process in terms of calculation accuracy.

[0008] To achieve the above objectives, this invention provides a numerical simulation method for the load release process based on dynamic nested mesh technology, comprising the following steps:

[0009] S1. Establish the geometric model of the launch vehicle and its loads;

[0010] S2. The 3D model is meshed using a nested meshing method. Based on the motion of the separation system, the separation system is divided into three mesh layers and several mesh clusters, and the mesh density is increased in the region of concentrated motion.

[0011] S3. Based on the three-dimensional compressible unsteady Navier-Stokes equations and the RNG K-Epsilon turbulence model, a physical model of airflow during the separation process is established.

[0012] S4. Based on the multi-mass spring damping model, simulate the traction effect of the umbrella on the load;

[0013] S5. Based on the penalty function method, simulate the collision between the adapter and the load during the separation process;

[0014] S6. Based on the multi-rigid-body dynamics equations and kinematic equations, establish a six-degree-of-freedom motion model of the vehicle, load, and parachute, and calculate their respective motion states.

[0015] S7. Establish monitoring points on the carrier and the payload to monitor the minimum distance between the carrier and the payload at each time step in order to determine the safety during the separation process;

[0016] S8. Input the initial motion state of the carrier, load, and parachute, as well as the initial and boundary conditions of the separation environment. Solve the calculation and output the pressure cloud map, displacement, velocity, and attitude change curves of the carrier and load, and minimum distance change curve of the carrier and load.

[0017] Preferably, the 3D model is meshed in S2, including the following steps:

[0018] S2.1 Estimate the range of motion of each component during the separation process, and refine the motion trajectory area in the background mesh layer;

[0019] S2.2. The nested mesh is managed in layers. According to the initiative of the movement of each component and the priority of information transmission, the nested mesh is divided into three mesh layers and multiple mesh clusters. Among them, the multi-mass spring damping model representing the traction parachute and the sling has the highest priority. Each mass is a mesh cluster. This layer of mesh is nested in the mesh cluster of the load and provides donor elements to the mesh cluster of the load for information transmission.

[0020] S2.3 The grid clusters of the carrier and the grid clusters of the payload have secondary priority. The two grid clusters are nested together. The boundary between adjacent grids is determined by the distance between the walls. This grid layer provides donor cells to the background grid for data interpolation.

[0021] S2.4 The background mesh layer has the lowest priority and is responsible for providing receptor units to receive information transmitted by the motion mesh layer.

[0022] Preferably, establishing a physical model of airflow during the separation process in S3 includes the following steps:

[0023] S3.1 Establish a three-dimensional compressible unsteady Navier-Stokes equation model, the integral form of which is:

[0024]

[0025] Where W is a conserved variable, F, F v Let Ω be the flow vector term and the dissipation flow vector term, respectively; dS be the area of ​​the infinitesimal element; Re be the Reynolds number; and n be the normal vector outside the control body.

[0026] S3.2. The flow control equations are discretized using the finite volume method. A second-order upwind scheme is used for spatial discretization on the computational grid, and a k-order backward differential discretization is used in the time direction. The discretized equations are:

[0027]

[0028]

[0029] in, j represents the neighboring control bodies of control body i, ij represents the interface between control bodies i and j, nface represents the number of boundary surfaces of the current control body, W is a conserved variable, and φ n-h W n-hIt is the linear interpolation of the conserved variables in the control volume n, where RE is the Reynolds number and F is the linear interpolation of the conserved variables in the control volume n. vij S is the viscous flux at the interface between control volumes i and j. ij This represents the area of ​​the interface.

[0030] S3.3. The RNG K-epsilon turbulence model is used to handle the RANS method. Compared with the standard K-epsilon model, the RNG K-epsilon two-equation turbulence model adds conditions to the epsilon equations and also incorporates turbulent vortices into the calculation, thus improving the accuracy of the calculation. Its equations for turbulent kinetic energy and turbulent dissipation rate ε are as follows:

[0031]

[0032]

[0033] Among them, G k G is the turbulent kinetic energy generated by the average velocity gradient. b Y is the turbulent kinetic energy generated by buoyancy. M The effect of wave expansion on turbulent dissipation rate;

[0034] S3.4 In the standard K-epsilon model, the turbulent viscosity correlation constant is an empirical coefficient obtained from experiments, while the RNG K-epsilon model is derived from theoretical formulas. The formula for calculating turbulent viscosity at high Reynolds numbers is as follows:

[0035]

[0036] Where, μ t It is the turbulent viscosity, ρ is the air density, and C is the turbulent viscosity. μ ε is a constant coefficient, which is 0.09 in the RNG K-epsilon model. k and ε are two fundamental unknowns.

[0037] Preferably, in S4, a multi-mass spring-damped model is used instead of the traction parachute, specifically including the following steps:

[0038] S4.1. A point mass spring-damped model considering mass movement is used to simulate the traction effect of the umbrella on the load.

[0039] S4.2, The point mass spring damping model discretizes the rope into multiple point masses, with adjacent point masses connected by a spring damping force. The spring damping force is:

[0040]

[0041] Where F is the tension of each rope segment, k is the equivalent stiffness of each rope segment, c is the equivalent damping of each rope segment, Δl is the elongation of each rope segment, v is the relative velocity in the direction of the line connecting the endpoints of each rope segment, l0 is the original length of each rope segment, and l is the total length of the rope after elongation.

[0042] Preferably, the stiffness coefficient of each rope segment is:

[0043]

[0044] Where n is the number of ropes connected in parallel, F b δ represents the breaking strength of each rope, and δ represents the breaking elongation of each rope.

[0045] The damping coefficients for each rope segment are:

[0046]

[0047] Where ξ is the damping ratio, and in practice, most cases are underdamped, and generally ξ < 0.2, m p The mass of the traction parachute; in this embodiment, k is 7 × 10⁻⁶. 4 N / m, c is 800 kg / s.

[0048] Preferably, in S5, the collision between the adapter and the load during the separation process is simulated based on the penalty function method, specifically including the following steps:

[0049] Using the inner wall of the carrier aircraft as the primary contact surface and the outer surface of the load as the secondary contact surface, during the calculation, if the distance from the secondary node to the primary contact surface is less than 12mm, a force is applied between the secondary node and the primary contact surface. The formula for the force is:

[0050] F1=k1d

[0051] Where F1 is the contact force, k1 is the adapter stiffness (i.e., the elastic coefficient), and d is the penetration depth.

[0052] Preferably, in S6, the coordinate system is divided into two types: an inertial coordinate system fixed to the earth and a non-inertial coordinate system fixed to the moving object. The non-inertial coordinate system moves with the moving object in the inertial coordinate system with six degrees of freedom. The coordinate transformation relationship from the inertial system to the non-inertial system can be represented by a set of Euler angles. The transformation relationship between the two is as follows:

[0053]

[0054] In this system, the subscript b represents a non-inertial frame of reference, and the subscript I represents an inertial frame of reference. Both are Cartesian coordinate systems that follow the right-hand rule. θ is the pitch angle, φ is the roll angle, and ψ is the yaw angle.

[0055] Preferably, in S6, by solving the kinematic equations of the center of mass of each moving object and the kinematic equations of rotation about the center of mass, the position and motion attitude angle of the center of mass of each moving object within each time step can be obtained. According to Newton's second law, the equation of motion of the center of mass of a rigid body in an inertial coordinate system is:

[0056]

[0057] Where m is the mass of the rigid body, V is the velocity vector of the rigid body's center of mass relative to the inertial frame, d1 is the differential of the velocity vector V with respect to time t, and F2 is the resultant vector of external forces acting on the rigid body.

[0058] Preferably, rotational motion can be described using a non-inertial coordinate system established at the center of mass of each object, and the equation of rotation about the center of mass is:

[0059]

[0060] Where I is the inertia tensor, The rotational angular velocity in body coordinates. Angular acceleration, It is the torque acting on a moving object.

[0061] Preferably, S7 uses the feature point method to monitor the minimum distance between the carrier and the payload in real time to determine the safety during the separation process, specifically including the following steps:

[0062] D1 and D2 on the load within the plumb plane of symmetry, and D3 and D4 on the carrier aircraft, are the most likely points of collision between the load and the carrier aircraft. The distances from D1 to the upper wall of the carrier aircraft's inner cabin, D2 to the lower wall of the carrier aircraft's inner cabin, D3 to the upper wall of the load, and D4 to the lower wall of the load are monitored in real time at each time step. The minimum value among the four distances is the minimum distance between the load and the carrier aircraft's inner cabin at that time step.

[0063] The beneficial effects of this invention are as follows:

[0064] I. The dynamic nested mesh method, as a type of computational fluid dynamics coupled multibody dynamics (CFD-6DOF method), can solve the aerodynamic forces acting on each moving object at every time step in real time, and the accuracy of the aerodynamic data is significantly higher than that of the former. When the calculation conditions such as separation height and velocity need to be changed, the dynamic nested mesh method is also far more flexible in operation than the traditional DMS method.

[0065] II. The dynamic nested mesh method has high accuracy in calculating unsteady transient aerodynamic forces of multibody motion. In recent years, this method has been widely used in numerical simulations of other forms of multibody separation systems, such as aircraft-launched projectile systems and submunition separation systems.

[0066] Third, due to the limitation of the limited number of computational nodes in the dynamic model, the DMS method often only estimates the compression and contact force of the adapter through a subset of feature points when calculating the contact collision between the load and the adapter. In reality, the collision process between the load and the adapter is a process with extremely uneven spatial distribution and violent temporal oscillation. The penalty function method can monitor the degree of contact compression between all mesh nodes on the load surface and the adapter surface at each time step, thus simulating the contact force of the adapter on the load more accurately. 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 embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0068] Figure 1 This is a flowchart illustrating the numerical simulation method for the payload exit process in an embodiment of the present invention.

[0069] Figure 2 This is a schematic diagram of a three-dimensional model of the carrier in an embodiment of the present invention;

[0070] Figure 3 This is a schematic diagram of the three-dimensional load model in an embodiment of the present invention;

[0071] Figure 4 This is a schematic diagram of the separation system and the distribution of lifting points in an embodiment of the present invention;

[0072] Figure 5 This is a flowchart illustrating the process of solving multi-body separation using dynamic nested mesh technology in an embodiment of the present invention.

[0073] Figure 6 This is a schematic diagram of the background grid layer in an embodiment of the present invention;

[0074] Figure 7 This is a schematic diagram of the carrier and load grid layer in an embodiment of the present invention;

[0075] Figure 8 This is a schematic diagram of the mesh layer of the multi-mass spring damping model in an embodiment of the present invention;

[0076] Figure 9 This is a schematic diagram illustrating the effect of the penalty function simulating the adapter buffer in an embodiment of the present invention;

[0077] Figure 10 This is a schematic diagram of the Euler coordinate system and attitude angles in an embodiment of the present invention;

[0078] Figure 11This is a schematic diagram of the carrier and load hazard monitoring points in an embodiment of the present invention;

[0079] Figure 12 These are pressure field cloud diagrams at the initial moment, mid-separation moment, and separation completion moment in this embodiment of the invention;

[0080] Figure 13 The curves showing the variation of adapter friction and minimum distance between the carrier and load in an embodiment of the present invention are shown.

[0081] Figure 14 The figures show the variation curves of the attitude angles and angular velocities of the carrier and the payload in this embodiment of the invention. Detailed Implementation

[0082] To make the objectives, technical solutions, and advantages of the present invention clearer, the present invention will be further described in detail below with reference to specific embodiments and accompanying drawings.

[0083] It should be noted that, unless otherwise defined, the technical or scientific terms used in the embodiments of this invention should have the ordinary meaning understood by those skilled in the art to which this invention pertains. The terms "first," "second," and similar terms used in this invention do not indicate any order, quantity, or importance, but are merely used to distinguish different components. Terms such as "comprising" or "including" mean that the element or object preceding the word encompasses the elements or objects listed following the word and their equivalents, without excluding other elements or objects. Terms such as "connected" or "linked" are not limited to physical or mechanical connections, but can include electrical connections, whether direct or indirect. Terms such as "upper," "lower," "left," and "right" are used only to indicate relative positional relationships; when the absolute position of the described object changes, the relative positional relationship may also change accordingly.

[0084] like Figure 1 , Figure 2 , Figure 3 , Figure 4 , Figure 5 , Figure 6 , Figure 7 , Figure 8 , Figure 9 , Figure 10 , Figure 11 , Figure 12 , Figure 13 , Figure 14 As shown, for reference Figure 1 A numerical simulation method for load excursion process based on dynamic nested mesh technology includes the following steps:

[0085] S1. Establish the geometric model of the launch vehicle and its loads;

[0086] This embodiment uses CATIA software to build the geometric model of the vehicle and load, referring to... Figure 2 and Figure 3 , Figure 2 The image shows a geometric model of a launch vehicle that resembles a missile warhead. Figure 3 The figure shows a load-bearing geometric model with a cylindrical shape. The relevant parameters of the model are shown in Tables 1 and 2.

[0087] Table 1 Geometric parameters of the launch vehicle

[0088]

[0089] Table 2 Load Geometric Parameters

[0090]

[0091]

[0092] Figure 4 (a) is a schematic diagram of the overall separation system in an embodiment of the present invention. The traction parachute is connected to the load via connecting ropes and three slings. A deformable adapter resembling a launch tube is installed on the inner wall of the launch vehicle to prevent the load from colliding with the launch vehicle during separation. The nominal area of ​​the traction parachute is 32m². 2 suspenders Figure 4 (b) As shown, the adapters are uniformly distributed circumferentially within a circle with a radius of 0.18 m; the initial thickness of the adapter is 0.012 m, and the stiffness is 1.0 e. 7 Pa / m, the coefficient of kinetic friction is 0.1.

[0093] S2. The 3D model is meshed using a nested meshing method. Based on the motion of the separation system, the separation system is divided into three mesh layers and several mesh clusters, and the mesh density is increased in the region of concentrated motion.

[0094] Meshing a 3D model in S2 includes the following steps:

[0095] S2.1 Estimate the range of motion of each component during the separation process, and refine the motion trajectory area in the background mesh layer;

[0096] S2.2. The nested mesh is managed in layers. According to the initiative of the movement of each component and the priority of information transmission, the nested mesh is divided into three mesh layers and multiple mesh clusters. Among them, the multi-mass spring damping model representing the traction parachute and the sling has the highest priority. Each mass is a mesh cluster. This layer of mesh is nested in the mesh cluster of the load and provides donor elements to the mesh cluster of the load for information transmission.

[0097] S2.3 The grid clusters of the carrier and the grid clusters of the payload have secondary priority. The two grid clusters are nested together. The boundary between adjacent grids is determined by the distance between the walls. This grid layer provides donor cells to the background grid for data interpolation.

[0098] S2.4 The background mesh layer has the lowest priority and is responsible for providing receptor units to receive information transmitted by the motion mesh layer.

[0099] Figure 5 This is a flowchart of the dynamic nested mesh technology for solving the multibody separation process in this embodiment of the invention. At each time step, the CFD solver calculates the transient aerodynamic forces and torques on the multibody system by solving the unsteady flow field, and transmits the data to the six-degree-of-freedom module. The module calculates the motion of the multibody system based on the forces and torques, and moves it to the corresponding position by the nested mesh, and enters the next time step. The CFD solver then calculates the forces on the multibody system in this state, and so on iteratively to achieve accurate simulation of the dynamic motion process of the load leaving the cabin.

[0100] The hierarchical management of nested meshes has a significant impact on the success of mesh nesting. Cartesian meshes have the advantages of easy mesh generation and high mesh quality, and are suitable for complex structures with large-scale motion. Cartesian meshes are used to perform computational meshing on the separated system, and the mesh model is divided into three mesh layers according to the motion of the separated system.

[0101] Figure 6 This is the first layer of the grid, which is a fixed background grid layer. This layer has only one grid cluster and has the lowest priority in the computational model. It is responsible for providing the receiver unit to receive the data transmitted from the moving grid interpolation. The central area of ​​the grid is the region where the motion trajectory of the separation system is separated. Local densification is performed to ensure the accuracy of the model. The grid gradually becomes sparse from the center to the surrounding areas. The number of grids in this layer is 17.84 million.

[0102] Figure 7 This is the second mesh layer, a moving mesh layer nested within the background mesh. This layer contains two mesh clusters: one for the launcher and one for the payload. These two clusters are nested within each other, with the boundary between adjacent meshes determined by the distance between their walls. The two mesh clusters then act as active unit layers, providing donor units to the background mesh to transmit information. The walls of both mesh clusters are enclosed by prism-structured meshes, and the mesh size smoothly transitions from the wall surface to the far field. The boundary dimensions of the nested meshes are similar to the dimensions of the local refinement zones in the background mesh, ensuring the accuracy of the interpolation. The launcher mesh cluster has 4.09 million meshes, and the payload mesh cluster has 1.88 million meshes.

[0103] Figure 8This is the third mesh layer, embedded within the load's mesh cluster, and has the highest priority in nested mesh information transfer. Since the parachute canopy, lines, and straps are replaced by a multi-mass spring-damped model (see S4 for details), each mass is represented by a small sphere with extremely small mass and volume. Each sphere is a mesh cluster, and all the sphere mesh clusters together constitute this active mesh layer, embedded within the load's mesh cluster, providing information through boundary flow field data exchange and interpolation. This mesh layer has a total of 220,000 meshes.

[0104] S3. Based on the three-dimensional compressible unsteady Navier-Stokes equations and the RNG K-Epsilon turbulence model, a physical model of airflow during the separation process is established.

[0105] S3.1 Establish a three-dimensional compressible unsteady Navier-Stokes equation model, the integral form of which is:

[0106]

[0107] Where W is a conserved variable, F, F v Let Ω be the flow vector term and the dissipation flow vector term, respectively; dS be the area of ​​the infinitesimal element; Re be the Reynolds number; and n be the normal vector outside the control body.

[0108] S3.2. The flow control equations are discretized using the finite volume method. A second-order upwind scheme is used for spatial discretization on the computational grid, and a k-order backward differential discretization is used in the time direction. The discretized equations are:

[0109]

[0110]

[0111] in, j represents the neighboring control bodies of control body i, ij represents the interface between control bodies i and j, nface represents the number of boundary surfaces of the current control body, W is a conserved variable, and φ n-h W n-h It is the linear interpolation of the conserved variables in the control volume n, where RE is the Reynolds number and F is the linear interpolation of the conserved variables in the control volume n. vij S is the viscous flux at the interface between control volumes i and j. ij This represents the area of ​​the interface.

[0112] S3.3. The RNG K-epsilon turbulence model is used to handle the RANS method. The RANS method, a numerical simulation method for fluid mechanics, is a solution method based on the Reynolds-averaged Navier-Stokes equations. Compared with the standard K-epsilon method, the RNG K-epsilon two-equation turbulence model adds conditions to the epsilon equations and also incorporates turbulent vortices into the calculation, thus improving the accuracy of the calculation. Its equations for turbulent kinetic energy and turbulent dissipation rate ε are as follows:

[0113]

[0114]

[0115] Among them, G k G is the turbulent kinetic energy generated by the average velocity gradient. b Y is the turbulent kinetic energy generated by buoyancy. M The effect of wave expansion on turbulent dissipation rate;

[0116] S3.4 In the standard K-epsilon model, the turbulent viscosity correlation constant is an empirical coefficient obtained from experiments, while the RNG K-epsilon model is derived from theoretical formulas. The formula for calculating turbulent viscosity at high Reynolds numbers is as follows:

[0117]

[0118] Where, μ t It is the turbulent viscosity, ρ is the air density, and C is the turbulent viscosity. μ ε is a constant coefficient, which is 0.09 in the RNG K-epsilon model. k and ε are two fundamental unknowns.

[0119] S4. Based on the multi-mass spring damping model, simulate the traction effect of the umbrella on the load;

[0120] Since the canopy, ropes, and straps of a parachute are all deformable flexible bodies, it is difficult to handle such problems in CFD. Therefore, a mass-spring-damped model that considers mass movement is used to simulate the traction effect of the parachute on the load.

[0121] In S4, a multi-mass spring-damped model is used to replace the traction parachute, specifically including the following steps:

[0122] S4.1. A point mass spring-damped model considering mass movement is used to simulate the traction effect of the umbrella on the load.

[0123] S4.2, The point mass spring damping model discretizes the rope into multiple point masses, with adjacent point masses connected by a spring damping force. The spring damping force is:

[0124]

[0125] Where F is the tension of each rope segment, k is the equivalent stiffness of each rope segment, c is the equivalent damping of each rope segment, Δl is the elongation of each rope segment, v is the relative velocity in the direction of the line connecting the endpoints of each rope segment, l0 is the original length of each rope segment, and l is the total length of the rope after elongation.

[0126] The stiffness coefficient of each rope segment is

[0127]

[0128] In the formula, n is the number of ropes connected in parallel, F b δ represents the breaking strength of each rope, and δ represents the breaking elongation of each rope.

[0129] The damping coefficients for each rope segment are:

[0130]

[0131] In the formula, ξ is the damping ratio. In practice, most cases are underdamped, and generally ξ < 0.2, m p The mass of the traction parachute. In this embodiment, k is 7 × 10⁻⁶. 4 N / m, c is 800 kg / s.

[0132] S5. Based on the penalty function method, simulate the collision between the adapter and the load during the separation process;

[0133] Since the collision process between the adapter and the load involves contact between the surfaces, and the adapter is a variable elastic body, such problems are often difficult to solve in CFD. Therefore, in this embodiment, a solid model of the adapter is not established. Instead, the penalty function method is used to simulate the collision process between the load and the adapter, so as to realize the volume coupling effect between the load and the carrier.

[0134] like Figure 9 As shown, the inner wall of the carrier aircraft is the primary contact surface, and the outer surface of the load is the secondary contact surface. During the calculation, at each time step, each node on the secondary contact surface is checked sequentially. The distance from the secondary node to the primary contact surface is calculated to see if it is less than 12mm (the initial thickness of the adapter when it is not compressed). If not, it is assumed that no contact force is generated. If the distance is less than 12mm, the node has penetrated the initial surface of the adapter when it is not compressed, and a force (i.e., contact force) is applied between the secondary node and the primary contact surface. The formula for the force is:

[0135] F1=k1d

[0136] Where F1 is the contact force, k1 is the adapter stiffness (i.e., the elastic coefficient), and d is the penetration depth.

[0137] This approach is equivalent to placing a spring between the slave node and the master surface to buffer the adapter against changes in load attitude. Since the collision process between the load and the carrier is highly nonlinear, the calculated results of the pressure and friction forces acting on the adapter are also oscillating.

[0138] S6. Based on the multi-rigid-body dynamics equations and kinematic equations, establish a six-degree-of-freedom motion model of the vehicle, load, and parachute, and calculate their respective motion states.

[0139] like Figure 10 As shown, the coordinate systems in this embodiment of the invention are divided into two types: an inertial coordinate system fixed to the earth and a non-inertial coordinate system fixed to the moving object (carrier, separation body, parachute). The non-inertial coordinate system moves with the moving object in the inertial coordinate system with six degrees of freedom. The coordinate transformation relationship from the inertial frame to the non-inertial frame can be represented by a set of Euler angles, and the transformation relationship between the two is as follows:

[0140]

[0141] In the formula, the subscript b represents the non-inertial frame of reference, and the subscript I represents the inertial frame of reference; both are Cartesian coordinate systems following the right-hand rule. The inertial frame's X-axis lies in the horizontal plane, with its positive direction pointing from the front to the rear of the vehicle. The Y-axis lies in the vertical plane, with its positive direction pointing vertically upwards. The Z-axis lies in the horizontal plane, with its direction determined by the right-hand rule. The non-inertial frame of reference is fixed to the moving object. Its X-axis coincides with the object's axial direction, with its positive direction pointing from the head to the tail. Its Y-axis lies in the vertical plane of symmetry, coinciding with the object's longitudinal axis, with its positive direction pointing upwards. The Z-axis direction is determined by the right-hand rule. θ is the pitch angle, φ is the roll angle, and ψ is the yaw angle. In this embodiment, which involves planar motion, only the pitch angle changes among the Euler angles.

[0142] By solving the kinematic equations of the center of mass of each moving object and the kinematic equations of rotation about the center of mass, the position and attitude angle of the center of mass of each moving object within each time step can be obtained. According to Newton's second law, the equation of motion of the center of mass of a rigid body in an inertial coordinate system is:

[0143]

[0144] Where m is the rigid body mass, V is the velocity vector of the rigid body's center of mass relative to the inertial frame, d1 is the differential of the velocity vector V with respect to time t, and F2 is the resultant vector of external forces acting on the rigid body. Rotational motion can be described by establishing a body coordinate system (non-inertial frame) at the center of mass of each object. The equation of rotational motion about the center of mass is:

[0145]

[0146] In the formula, I is the inertia tensor. The rotational angular velocity in body coordinates. Angular acceleration, It is the torque acting on a moving object.

[0147] S7. Establish monitoring points on the carrier and the payload to monitor the minimum distance between the carrier and the payload at each time step in order to determine the safety during the separation process;

[0148] like Figure 11 As shown, points D1 and D2 on the load within the plumb plane of symmetry, and points D3 and D4 on the carrier aircraft, are the most likely points of collision between the load and the carrier aircraft. Real-time monitoring is performed on the distances from D1 to the upper wall of the carrier aircraft's inner cabin, D2 to the lower wall of the carrier aircraft's inner cabin, D3 to the upper wall of the load, and D4 to the lower wall of the load at each time step. The minimum of these four distances is the minimum clearance between the load and the carrier aircraft's inner cabin at that time step. When the adapter is not compressed, the initial clearance between the carrier aircraft and the load is 12mm. If the minimum clearance is greater than 4mm throughout the entire exit process, the load exit safety under this condition is considered relatively high.

[0149] S8. Input the initial motion state of the carrier, load, and parachute, as well as the initial and boundary conditions of the separation environment. Perform the calculation and output the pressure cloud map, displacement, velocity, and attitude change curves of the carrier and load, and minimum distance change curve of the carrier and load.

[0150] The relevant parameters calculated in this embodiment are: separation altitude of 25 km, Mach number of 0.8, angle of attack and pitch angle of 0°, and pitch rate of 0 (positive for head-up). The impact of these factors on the safety and engineering applicability of the separation process can be studied individually using the controlled variable method. Using a supercomputing cluster with 256 CPU cores for parallel computing, a single calculation takes approximately 40 hours.

[0151] Figure 12 (ac) are flow field pressure cloud diagrams at 0.1s, 0.3s, and 0.5s, respectively, at the start of separation. In the early to mid-stages of separation, the launch vehicle and the payload are under the same pressure field. As the payload moves backward relative to the launch vehicle, the space originally occupied by the payload in the launch vehicle's interior becomes a near-vacuum local space with extremely low pressure. At 0.5s, the payload and the launch vehicle have completely separated and are separated by a certain distance. Their pressure fields are independent of each other. At this point, it can be considered that the payload is basically unaffected by the aerodynamic interference of the launch vehicle.

[0152] Figure 13 (ab) represent the magnitude of the frictional force on the adapter during the separation process and the minimum distance between the carrier and the load, respectively. Figure 14The figures represent the attitude angles and angular velocities of the launch vehicle and payload during separation. In the initial stage of exoplanetion, due to the sudden pull of the parachute and the maximum contact area between the payload and the adapter, the adapter friction reaches its peak of approximately 200 N. Subsequently, the adapter friction fluctuates with the degree of collision between the payload and the adapter. As the payload's exit distance increases, the overall fluctuation of the adapter friction shows a decreasing trend. At 0.389 s, the payload and launch vehicle completely separate, and the adapter ceases operation. In the early and mid-stages of separation, the adapter has a strong attitude-limiting effect, and the pitch angles of the launch vehicle and payload change in the same direction. As the payload gradually exits the capsule, the equivalent point of application of its force on the launch vehicle shifts backward, increasing the lever arm of the force and increasing the pitch angle of the launch vehicle. The launch vehicle, through the coupling effect of the adapter, also causes the payload's pitch angle to increase. In the later stage of separation, since most of the payload has exited the capsule, the contact area with the adapter is significantly reduced, and the coupling effect between them decreases. Under the continuous action of aerodynamic torque and parachute pull, the payload gradually pitches down, with an accelerating trend. The minimum clearance is positively correlated with the relative pitch angle between the carrier and the payload. The adapter has a good limiting effect in the early and mid-stages of separation. Although the collision between the carrier and the payload is intense, the minimum clearance is always maintained above 10mm. In the later stage of separation, the payload tends to pitch down faster than the carrier, and the relative pitch angle increases, which causes the danger point 3 to be compressed more, and the minimum clearance decreases. However, the minimum clearance is always maintained above 6mm, and the payload can be safely separated.

[0153] Those skilled in the art should understand that the discussion of any of the above embodiments is merely exemplary and is not intended to imply that the scope of the invention is limited to these examples; within the framework of the invention, the technical features of the above embodiments or different embodiments can also be combined, the steps can be implemented in any order, and there are many other variations of the different aspects of the invention as described above, which are not provided in detail for the sake of brevity.

[0154] The embodiments of this invention are intended to cover all such substitutions, modifications, and variations falling within the broad scope of the appended claims. Therefore, any omissions, modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of this invention should be included within the scope of protection of this invention.

Claims

1. A method for numerical simulation of the load out-of-cabin process based on dynamic nested grid technology, characterized in that, Includes the following steps: S1. Establish the geometric model of the launch vehicle and its loads; S2. The 3D model is meshed using a nested meshing method. Based on the motion of the separation system, the separation system is divided into three mesh layers and several mesh clusters, and the mesh density is increased in the region of concentrated motion. Meshing a 3D model in S2 includes the following steps: S2.1 Estimate the range of motion of each component during the separation process, and refine the motion trajectory area in the background mesh layer; S2.

2. The nested mesh is managed in layers. According to the initiative of the movement of each component and the priority of information transmission, the nested mesh is divided into three mesh layers and multiple mesh clusters. Among them, the multi-mass spring damping model representing the traction parachute and the sling has the highest priority. Each mass is a mesh cluster. This layer of mesh is nested in the mesh cluster of the load and provides donor elements to the mesh cluster of the load for information transmission. S2.3 The grid clusters of the carrier and the grid clusters of the payload have secondary priority. The two grid clusters are nested together. The boundary between adjacent grids is determined by the distance between the walls. This grid layer provides donor cells to the background grid for data interpolation. S2.4 The background mesh layer has the lowest priority and is responsible for providing receptor units to receive information transmitted by the motion mesh layer; S3. Based on the three-dimensional compressible unsteady Navier-Stokes equations and the RNG K-Epsilon turbulence model, a physical model of airflow during the separation process is established. S4. Based on the multi-mass spring damping model, simulate the traction effect of the umbrella on the load; In S4, a multi-mass spring-damped model is used to replace the traction parachute, specifically including the following steps: S4.

1. A point mass spring-damped model considering mass movement is used to simulate the traction effect of the umbrella on the load. S4.2, The point mass spring damping model discretizes the rope into multiple point masses, with adjacent point masses connected by a spring damping force. The spring damping force is: in, For the tension of each rope segment, The equivalent stiffness of each rope segment, For the equivalent damping of each rope segment, This refers to the elongation of each rope segment. Let be the relative velocity along the line connecting the endpoints of each rope segment. The original length of each rope segment It is the total length of the rope after it has been stretched. S5. Based on the penalty function method, simulate the collision between the adapter and the load during the separation process; S6. Based on the multi-rigid-body dynamics equations and kinematic equations, establish a six-degree-of-freedom motion model of the vehicle, load, and parachute, and calculate their respective motion states. S7. Establish monitoring points on the carrier and the payload to monitor the minimum distance between the carrier and the payload at each time step in order to determine the safety during the separation process; S8. Input the initial motion state of the carrier, load, and parachute, as well as the initial and boundary conditions of the separation environment. Solve the calculation and output the pressure cloud map, displacement, velocity, and attitude change curves of the carrier and load, and minimum distance change curve of the carrier and load.

2. The numerical simulation method for load exit process based on dynamic nested mesh technology according to claim 1, characterized in that, Establishing a physical model of airflow during the separation process in S3 includes the following steps: S3.1 Establish a three-dimensional compressible unsteady Navier-Stokes equation model, the integral form of which is: in, As a conserved variable, , These are the flow vector term and the dissipation flow vector term, respectively. For any control body, Let the area be the infinitesimal element. The Reynolds number is... It controls the external normal vector; S3.

2. The flow control equations are discretized using the finite volume method. A second-order upwind scheme is used for spatial discretization on the computational grid, and a k-order backward differential discretization is used in the time direction. The discretized equations are: in, , For control body Neighbor control body Indicates control body and control body The interface, This represents the number of boundary surfaces of the current control volume. It is a conserved variable. It is the linear interpolation of the conserved variables in the control volume n, where RE is the Reynolds number. It is the control body and control body Viscous flux at the interface This represents the area of ​​the interface. S3.

3. The RNG K-epsilon turbulence model is used to handle the RANS method. Compared with the standard K-epsilon model, the RNG K-epsilon two-equation turbulence model adds conditions to the epsilon equations and also incorporates turbulent vortices into the calculation, thus improving the calculation accuracy. Its turbulent kinetic energy and turbulent dissipation rate are also improved. The equation is: in, The turbulent kinetic energy generated by the average velocity gradient. The turbulent kinetic energy generated by buoyancy, The effect of wave expansion on turbulent dissipation rate; S3.4 In the standard K-epsilon model, the turbulent viscosity correlation constant is an empirical coefficient obtained from experiments, while the RNG K-epsilon model is derived from theoretical formulas. The formula for calculating turbulent viscosity at high Reynolds numbers is as follows: in, It is turbulent viscosity. air density, The constant coefficient is 0.09 in the RNG K-epsilon model. and These are two fundamental unknowns.

3. The numerical simulation method for load exit process based on dynamic nested mesh technology according to claim 1, characterized in that, The stiffness coefficients of each rope segment are: in, The number of ropes connected in parallel. For the breaking strength of each rope, The breaking elongation of each rope; The damping coefficients for each rope segment are: in, The damping ratio is used, but in reality, most cases are underdamped, and... <0.2, For the mass of the traction parachute; 7×10 4 N / m, It is 800 kg / s.

4. The numerical simulation method for load exit process based on dynamic nested mesh technology according to claim 1, characterized in that, In S5, the collision between the adapter and the load during the separation process is simulated based on the penalty function method, specifically including the following steps: Using the inner wall of the carrier aircraft as the primary contact surface and the outer surface of the load as the secondary contact surface, during the calculation, if the distance from the secondary node to the primary contact surface is less than 12mm, a force is applied between the secondary node and the primary contact surface. The formula for the force is: in, For contact force, This refers to the adapter stiffness, i.e., the elastic coefficient. This represents the penetration depth.

5. The numerical simulation method for load exit process based on dynamic nested mesh technology according to claim 1, characterized in that, In S6, coordinate systems are divided into two types: an inertial coordinate system fixed to the earth and a non-inertial coordinate system fixed to the moving object. The non-inertial coordinate system moves with the moving object in the inertial coordinate system with six degrees of freedom. The coordinate transformation relationship from the inertial frame to the non-inertial frame can be represented by a set of Euler angles. The transformation relationship between the two is as follows: Among them, subscript Indicates a non-inertial frame of reference, subscript Both represent inertial frames of reference and are Cartesian coordinate systems that follow the right-hand rule. For pitch angle, For roll angle, This is the yaw angle.

6. The numerical simulation method for load exit process based on dynamic nested mesh technology according to claim 5, characterized in that, In S6, by solving the kinematic equations of the center of mass of each moving object and the kinematic equations of rotation about the center of mass, the position and motion attitude angle of the center of mass of each moving object in each time step can be obtained. According to Newton's second law, the equation of motion of the rigid body's center of mass in the inertial coordinate system is: in, For rigid body mass, Let be the velocity vector of the rigid body's center of mass relative to the inertial frame. Refers to velocity vector Regarding time The differential, The resultant vector of external forces acting on the rigid body.

7. The numerical simulation method for load exit process based on dynamic nested mesh technology according to claim 6, characterized in that, Rotational motion can be described using a non-inertial coordinate system established at the center of mass of each object. The equation of rotational motion about the center of mass is: in, For inertial tensor, The rotational angular velocity in the body coordinate system. Angular acceleration, It is the torque acting on a moving object.

8. The numerical simulation method for load exit process based on dynamic nested mesh technology according to claim 1, characterized in that, S7 uses a feature point method to monitor the minimum distance between the launch vehicle and the payload in real time to determine the safety during the separation process. This includes the following steps: D1 and D2 on the load within the plumb plane of symmetry, and D3 and D4 on the carrier aircraft, are the most likely points of collision between the load and the carrier aircraft. The distances from D1 to the upper wall of the carrier aircraft's inner cabin, D2 to the lower wall of the carrier aircraft's inner cabin, D3 to the upper wall of the load, and D4 to the lower wall of the load are monitored in real time at each time step. The minimum value among the four distances is the minimum distance between the load and the carrier aircraft's inner cabin at that time step.

Images