A shock capturing method and system for aircraft aerodynamic characteristic analysis

Abstract

The application discloses a shock capturing method and system for analyzing aerodynamic characteristics of an aircraft, and relates to the field of aircraft design, and comprises the following steps: S1, problem modeling and discretization; S2, characteristic analysis and transformation preparation; S3, characteristic space dimension reduction; S4, nonlinear reconstruction in the dimension-reduced subspace; S5, full-dimension flux synthesis and physical space return; and S6, time advancement and flow field solving. By reducing the characteristic projection from five dimensions to two dimensions, the number of times of nonlinear reconstruction and the amount of matrix operation are significantly reduced, the calculation efficiency is improved, meanwhile, the shock capturing ability, high precision and non-oscillation characteristics of the original format are maintained, the method has strong universality and is easy to integrate. The method is suitable for fine simulation of complex flow structures in the flow around a supersonic aircraft.

Description

Technical Field

[0001] This invention relates to the field of aircraft design, and more specifically, to a shock wave capture method and system for analyzing the aerodynamic characteristics of aircraft. Background Technology

[0002] In the field of aircraft design, extremely high accuracy is required for the prediction of aerodynamic forces (such as lift and drag) and aerothermal forces (such as surface heat flux). Traditional characteristic projection methods, due to their high computational cost, force engineers to compromise between mesh resolution (affecting accuracy) and computation time (affecting design cycle time), limiting the refined design capabilities based on high-fidelity CFD. For supersonic aircraft, involving complex flows such as intense shock waves, even slight distortions in numerical simulations can lead to misjudgments in the design of thermal protection systems.

[0003] High-order, high-resolution shock capture schemes (such as WENO and TENO) are standard tools for accurately simulating supersonic flows involving complex phenomena such as shock waves and contact discontinuities. To ensure numerical stability and the absence of physical oscillations (ENO properties), these schemes require nonlinear reconstruction (i.e., eigenvector projection) in the characteristic space. This process involves frequent global matrix transformations of the flow field variables (left-multiplying by the left eigenvector matrix L, and then right-multiplying by the right eigenvector matrix R after reconstruction). For a three-dimensional problem, the eigenvector matrix is ​​5x5, and each projection and back-projection involves a large number of floating-point multiplication and addition operations, resulting in extremely high computational costs, which can account for more than 50% of the total computation time of the entire scheme, becoming a bottleneck restricting high-performance simulations.

[0004] Existing techniques such as Co-WENO (Common-weights WENO) partially reduce the number of reconstructions by finding a common set of weights for all variables, but do not reduce the computational complexity of the feature projection itself. Other methods attempt to disable feature projection in smooth regions, but determining smooth regions in complex flow fields introduces additional overhead and uncertainty. Therefore, there is an urgent need in this field for a new method that can fundamentally reduce the computational complexity of feature projection without sacrificing the format shock capture capability and oscillation-free characteristics.

[0005] Therefore, developing a shock wave capture method that can both ensure accuracy and significantly improve efficiency is of great significance for aircraft design and development. Summary of the Invention

[0006] The purpose of this invention is to significantly reduce the computational cost of aircraft design while maintaining computational accuracy and robustness, and ultimately achieve efficient and high-precision prediction of aerodynamic and aerothermal loads in aircraft design.

[0007] To achieve the above-mentioned objective, this invention provides a shock wave capture method for aerodynamic characteristic analysis of aircraft, the method comprising:

[0008] Step S1: Based on the incoming flow conditions and shape parameters of the aircraft, determine the local regions that may generate shock waves or flow discontinuities; obtain the geometric parameters of the local regions and construct the corresponding local computational grid; on the local computational grid, spatially discretize the gas dynamics control equations describing the flow to obtain a discretized set of equations, wherein the convective flux term in the discretized set of equations represents the interface numerical flux to be determined.

[0009] Step S2: According to the requirements of the shock wave capture format, on the interface of the grid cell where flux reconstruction is required, perform eigenvalue decomposition on the local flow Jacobian matrix to obtain the eigenvalues ​​of the local flow Jacobian matrix, the full-dimensional left eigenvector matrix L, and the full-dimensional right eigenvector matrix R; wherein, the full-dimensional left eigenvector matrix L is used to project the variables in the physical space to the eigenspace, and the full-dimensional right eigenvector matrix R is used to inversely transform the variables in the eigenspace back to the physical space;

[0010] Step S3: Construct the first dimension-reduced eigenvector l using the eigenvectors in the full-dimensional left eigenvector matrix L. c1 Second dimensionality reduction eigenvector l c2 , where l c1 It is a characteristic vector representing the convection mode, composed of the characteristic vectors corresponding to entropy waves and shear waves. c2 It is a weighted combination of the corresponding eigenvectors in the full-dimensional left eigenvector matrix L and the full-dimensional right eigenvector matrix R, representing the eigenvectors of the acoustic mode. c1 and l c2 This forms a two-dimensional reduced feature subspace, resulting in a reduced projection matrix L. red =[l c1 ;l c2 ];

[0011] Step S4: The physical flux F at the grid cell interface is transformed by the dimension reduction projection matrix L. red Projecting onto a two-dimensional reduced feature subspace yields two reduced feature variables F. a and F c For the characteristic variable F a and F c Interpolation is performed using a nonlinear reconstruction scheme, and two sets of nonlinear weighting factors ω, corresponding to the convection mode and the acoustic mode respectively, are extracted during the reconstruction process. k,a and ω k,c And obtain the reconstructed feature variable values ​​on both sides of the interface in the two-dimensional reduced feature subspace;

[0012] Step S5: Adjust the nonlinear weighting factor ω k,a and ω k,cExtending to the full-dimensional feature space, a full-dimensional diagonal weight matrix is ​​constructed; in the full-dimensional diagonal weight matrix, the values ​​corresponding to the first and fifth diagonal positions of the convection mode are assigned as ω. k,a The values ​​ω are assigned to the second to fourth diagonal positions corresponding to the acoustic mode. k,c If the full-dimensional diagonal weight matrix has other diagonal positions, then it is assigned a value of 1; based on the physical flux reconstruction value on each sub-template, the synthesis formula is derived by decomposing the full-dimensional diagonal weight matrix to calculate the numerical flux in the physical space.

[0013] Step S6: Repeat steps S2 to S5 for all grid cell interfaces that need to be processed until the final physical interface numerical flux of all grid cell interfaces in the entire field is calculated; substitute the final physical interface numerical flux of all grid cell interfaces in the entire field into the discretized equation set, perform time integration, and advance the flow field evolution until a convergent or transient flow field solution is obtained, thereby obtaining the shock wave structure or aerodynamic / thermal load distribution of the aircraft.

[0014] This method, based on physical mechanisms, aggregates five-dimensional feature patterns into two dominant modes: convection modes (entropy waves and shear waves) and acoustic modes. This allows for accurate capture of flow discontinuities even after dimensionality reduction, while significantly reducing computational complexity. By classifying feature patterns into two main categories, the left-hand feature variable is reduced from 5-dimensional to 2-dimensional, and from a 5×5 matrix to a 2×5 matrix. The number of floating-point multiplication and addition operations required for projection is drastically reduced from 25 to 10. Simultaneously, the number of feature variables requiring independent nonlinear reconstruction is reduced from five to two. The method is designed in two dimensions because the main discontinuities in the flow can be characterized by acoustic and convection modes, and the linear combination of the two feature vectors maximizes the preservation of the original shock wave capture performance. This significantly reduces the computational cost of aircraft design while maintaining computational accuracy and robustness, ultimately achieving efficient and high-precision prediction of aerodynamic and aerothermal loads in aircraft design.

[0015] Wherein, the two-dimensional reduced eigenspace is composed of the first reduced eigenvector l c1 Second dimensionality reduction eigenvector l c2 Composition, the first dimensionality-reduced eigenvector l c1 The second dimensionality-reduced eigenvector l is formed by a linear combination of eigenvectors related to entropy waves and shear waves. c2 It is composed of a linear combination of eigenvectors related to the left and right acoustic waves. The purpose of this design is to aggregate the five-dimensional feature patterns into two dominant modes based on physical mechanisms: the convection mode (entropy wave and shear wave) and the acoustic wave mode. This allows for accurate capture of flow discontinuities even after dimensionality reduction, while significantly reducing computational cost. The design is two-dimensional because the main discontinuities in the flow can be characterized by the acoustic and convection modes, and the linear combination of the two eigenvectors maximizes the preservation of the shock wave capture performance of the original scheme.

[0016] Among them, the second dimensionality reduction eigenvector l c2 The construction method is as follows:

[0017] ;

[0018] Where θ is the local covariant velocity, c is the speed of sound, and l1 and l5 are the row vectors of the left and right sound waves in the full-dimensional left eigenvector matrix L, respectively. The construction method of the second dimensionality reduction eigenvector aggregates the left and right sound waves through adaptive weighting, which not only greatly reduces the amount of computation, but also ensures the robustness and high accuracy of the method when capturing complex discontinuities such as shock waves, thus achieving a balance between efficiency and accuracy in the aerodynamic analysis of aircraft.

[0019] The nonlinear reconstruction format is either WENO, TENO, or WENN.

[0020] The expression for the full-dimensional diagonal weight matrix is ​​as follows:

[0021] ;

[0022] in, It is a diagonal weight matrix. It is an identity matrix.

[0023] The specific synthesis formula is as follows:

[0024] ;

[0025] in, Let q be the numerical flux in physical space. k Let be the reconstructed value of the physical flux F on the k-th sub-template, r1 and r5 be the column vectors corresponding to the left and right sound waves in the full-dimensional right eigenvector matrix R, l1 and l5 be the row vectors corresponding to the left and right sound waves in the full-dimensional left eigenvector matrix L, i represent the i-th grid point, and represent the grid point index.

[0026] The present invention also provides a shock wave capture system for aerodynamic characteristic analysis of aircraft, the system comprising:

[0027] The modeling and discretization unit is used to determine local regions that may generate shock waves or flow discontinuities based on the incoming flow conditions and shape parameters of the aircraft; to obtain the geometric parameters of the local regions and construct the corresponding local computational grid; and to spatially discretize the gas dynamics control equations describing the flow on the local computational grid to obtain a set of discretized equations, wherein the convective flux term in the set of discretized equations represents the interface numerical flux to be determined.

[0028] The feature analysis and transformation preparation unit is used to perform eigenvalue decomposition on the local flow Jacobian matrix on the grid cell interface where flux reconstruction is required, according to the requirements of the shock wave capture format, to obtain the eigenvalues ​​of the local flow Jacobian matrix, the full-dimensional left eigenvector matrix L, and the full-dimensional right eigenvector matrix R; wherein, the full-dimensional left eigenvector matrix L is used to project the variables in the physical space to the feature space, and the full-dimensional right eigenvector matrix R is used to inversely transform the variables in the feature space back to the physical space;

[0029] The dimensionality reduction eigenspace construction unit is used to construct the first dimensionality reduction eigenvector l using the eigenvectors in the full-dimensional left eigenvector matrix L. c1 Second dimensionality reduction eigenvector l c2 , where l c1 It is a characteristic vector representing the convection mode, composed of the characteristic vectors corresponding to entropy waves and shear waves. c2 It is a weighted combination of the corresponding eigenvectors in the full-dimensional left eigenvector matrix L and the full-dimensional right eigenvector matrix R, representing the eigenvectors of the acoustic mode. c1 and l c2 This forms a two-dimensional reduced feature subspace, resulting in a reduced projection matrix L. red =[l c1 ;l c2 ];

[0030] The dimensionality reduction space reconstruction and weight extraction unit is used to convert the physical flux F at the grid cell interface into the dimensionality reduction projection matrix L. red Projecting onto a two-dimensional reduced feature subspace yields two reduced feature variables F. a and F c For the characteristic variable F a and F c Interpolation is performed using a nonlinear reconstruction scheme, and two sets of nonlinear weighting factors ω, corresponding to the convection mode and the acoustic mode respectively, are extracted during the reconstruction process. k,a and ω k,c And obtain the reconstructed feature variable values ​​on both sides of the interface in the two-dimensional reduced feature subspace;

[0031] The interface numerical flux synthesis unit is used to synthesize the nonlinear weighting factor ω. k,a and ω k,c Extending to the full-dimensional feature space, a full-dimensional diagonal weight matrix is ​​constructed; in the full-dimensional diagonal weight matrix, the values ​​corresponding to the first and fifth diagonal positions of the convection mode are assigned as ω. k,a The values ​​ω are assigned to the second to fourth diagonal positions corresponding to the acoustic mode. k,cIf the full-dimensional diagonal weight matrix has other diagonal positions, then it is assigned a value of 1; based on the physical flux reconstruction value on each sub-template, the synthesis formula is derived by decomposing the full-dimensional diagonal weight matrix to calculate the numerical flux in the physical space.

[0032] The time-progression and flow field solution unit is used to process all the mesh cell interfaces that need to be processed by the feature analysis and transformation preparation unit, the dimension reduction feature subspace construction unit, the dimension reduction space reconstruction and weight extraction unit, and the interface numerical flux synthesis unit, and calculate the final physical interface numerical flux of all mesh cell interfaces in the entire field. The final physical interface numerical flux of all mesh cell interfaces in the entire field is substituted into the discretized equation system, and time integration is performed to advance the flow field evolution until a convergent or transient flow field solution is obtained, thereby obtaining the shock wave structure or aerodynamic / thermal load distribution of the aircraft.

[0033] Wherein, the two-dimensional reduced eigenspace is composed of the first reduced eigenvector l c1 Second dimensionality reduction eigenvector l c2 Composition, the first dimensionality-reduced eigenvector l c1 The second dimensionality-reduced eigenvector l is formed by a linear combination of eigenvectors related to entropy waves and shear waves. c2 It is composed of a linear combination of eigenvectors related to the left and right sound waves.

[0034] Among them, the second dimensionality reduction eigenvector l c2 The construction method is as follows:

[0035] ;

[0036] Where θ is the local covariant velocity, c is the speed of sound, and l1 and l5 are the row vectors corresponding to the left and right sound waves in the full-dimensional left eigenvector matrix L, respectively.

[0037] The nonlinear reconstruction format is either WENO, TENO, or WENN.

[0038] One or more technical solutions provided by this invention have at least the following technical effects or advantages:

[0039] The computational efficiency of aircraft design is significantly improved: by reducing 5D reconstruction to 2D, the number of nonlinear reconstructions is reduced by 60%. Simultaneously, the matrix operation dimension of feature projection is reduced from 5×5 to 2×5, significantly reducing floating-point operations. Real-world measurements show an efficiency improvement of over 20% for the WENO5-JS format and 40%~50% for the TENO5-A and WENN4-O formats. This allows for more refined mesh simulations of aircraft flows or the completion of numerous simulations within shorter design cycles, all within the same computational resources.

[0040] Maintaining the excellent properties of the format: The dimensionality reduction strategy based on physical mechanisms retains the shock wave capture capability, high accuracy and non-oscillation characteristics of the original format to the greatest extent. In particular, it ensures the non-oscillation design characteristics of strong discontinuities such as shock waves at the aircraft nose and leading edge of the wing, thereby obtaining accurate aerodynamic load distribution.

[0041] High versatility: This method is a framework improvement that can be used in conjunction with a variety of mainstream nonlinear formats (WENO-JS, WENO-Z, TENO, WENN, etc.) and is easy to integrate into existing aircraft CFD code.

[0042] Accuracy loss is controllable: Dimensionality reduction is based on the physical classification of flow modes, mainly affecting the modes that are reasonably aggregated, and has minimal impact on the accuracy of the overall solution.

[0043] Significant engineering application value: This invention directly addresses a core pain point in aircraft design—the contradiction between high-precision aerodynamic / thermal prediction and computational efficiency. Through this method, designers can obtain high-resolution flow field data sufficient to support refined design at an acceptable computational cost. In particular, the accurate capture of key parameters such as heat flux peaks, shock wave locations, and separation zone sizes significantly improves the reliability and economy of aircraft design. Attached Figure Description

[0044] The accompanying drawings, which are provided to further illustrate embodiments of the invention and constitute a part of this invention, are not intended to limit the scope of the invention.

[0045] Figure 1 This is a flowchart illustrating a shock wave capture method for analyzing the aerodynamic characteristics of an aircraft.

[0046] Figure 2 This is a comparative diagram showing the density distribution obtained by simulating the Lax shock tube problem using the method of this invention (WENO5-JS-RO), standard feature reconstruction (WENO5-JS), and featureless reconstruction (WENO5-JS-NCR);

[0047] Figure 3 This is a schematic diagram comparing the density cloud maps obtained by simulating the two-dimensional Riemann problem using the traditional method (left half) and the method of the present invention (right half). Detailed Implementation

[0048] To better understand the above-mentioned objectives, features, and advantages of the present invention, the present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments. It should be noted that, where there is no conflict, the embodiments of the present invention and the features thereof can be combined with each other.

[0049] Many specific details are set forth in the following description in order to provide a full understanding of the invention. However, the invention may also be practiced in other ways different from those described herein, and therefore the scope of protection of the invention is not limited to the specific embodiments disclosed below.

[0050] Example 1;

[0051] Please refer to Figure 1 , Figure 1 This is a flowchart illustrating a shock wave capture method for aerodynamic characteristic analysis of an aircraft. Embodiment 1 of the present invention provides a shock wave capture method for aerodynamic characteristic analysis of an aircraft, the method comprising:

[0052] Step S1, Problem Modeling and Discretization: Based on the incoming flow conditions and shape parameters of the aircraft, determine the local regions that may generate shock waves or flow discontinuities; obtain the geometric parameters of the local regions and construct the corresponding local computational grid; on the local computational grid, spatially discretize the gas dynamics governing equations (Navier-Stokes equations) describing the flow to obtain a discretized set of equations, wherein the convective flux term in the discretized set of equations represents the interface numerical flux to be determined;

[0053] Step S2, Feature Analysis and Transformation Preparation: According to the requirements of the shock wave capture format, on the interface of the grid cell where flux reconstruction is required, the local flow Jacobian matrix is ​​eigenvalued to obtain the eigenvalues, the full-dimensional left eigenvector matrix L, and the full-dimensional right eigenvector matrix R of the local flow Jacobian matrix; wherein, the full-dimensional left eigenvector matrix L is used to project the variables in the physical space to the feature space, and the full-dimensional right eigenvector matrix R is used to inversely transform the variables in the feature space back to the physical space;

[0054] Step S3, Feature Space Dimensionality Reduction: Using the eigenvectors of the full-dimensional left eigenvector matrix L, construct two new dimensionality-reduced eigenvectors l with clear physical meaning. c1 and l c2 Among them, l c1 It is a characteristic vector representing the convection mode, composed of the characteristic vectors corresponding to entropy waves and shear waves. c2 The eigenvectors representing the acoustic mode are formed by a weighted combination of the corresponding eigenvectors in the full-dimensional left eigenvector matrix L and the full-dimensional right eigenvector matrix R. These two vectors form a two-dimensional reduced-dimensional eigenspace and a reduced-dimensional projection matrix L. red =[l c1 ;l c2 ];

[0055] Step S4, Nonlinear Reconstruction in the Dimensionally Reduced Subspace: The physical flux F at the grid cell interface is transformed through the dimensionality-reduced projection matrix L. redDirectly projecting onto the two-dimensional feature subspace yields two dimensionality-reduced feature variables F. a and F c Subsequently, these two variables were independently interpolated using a nonlinear reconstruction scheme. During the reconstruction process, two sets of nonlinear weighting factors ω, corresponding to the convection mode and the acoustic mode respectively, were extracted. k,a and ω k,c And obtain the reconstructed feature variable values ​​on both sides of the interface in the two-dimensional reduced feature subspace;

[0056] Step S5, flux synthesis: The nonlinear weighting factor ω k,a and ω k,c Extending to the full-dimensional feature space, a full-dimensional diagonal weight matrix is ​​constructed; in the full-dimensional diagonal weight matrix, the values ​​corresponding to the first and fifth diagonal positions of the convection mode are assigned as ω. k,a The values ​​ω are assigned to the second to fourth diagonal positions corresponding to the acoustic mode. k,c If the full-dimensional diagonal weight matrix has other diagonal positions, it is assigned a value of 1; based on the physical flux reconstruction value on each sub-template, the synthesis formula is derived by decomposing the full-dimensional diagonal weight matrix to calculate the numerical flux in the physical space; the efficient synthesis formula derived by the above weight matrix decomposition directly calculates the numerical flux in the physical space, which is different from the traditional method of calculating the numerical flux by weighted fusion in the full-dimensional feature space and then back-projection, thus avoiding the expensive full-dimensional matrix operation;

[0057] Step S6, Time Advancement and Flow Field Solution: Repeat steps S2 to S5 for all grid cell interfaces that need to be processed until the final physical interface numerical flux of all grid cell interfaces in the entire field is calculated; substitute the final physical interface numerical flux of all grid cell interfaces in the entire field into the discretized equation set, perform time integration, advance the flow field evolution until a convergent or transient flow field solution is obtained, and obtain the shock wave structure or aerodynamic / thermal load distribution of the aircraft.

[0058] In one embodiment of the present invention, the determination of the local region where shock waves or flow discontinuities may occur in step S1 can be achieved in the following way:

[0059] First, theoretical predictions are made. Based on the aircraft's design Mach number and angle of attack, and using its nose radius of curvature, airfoil leading edge wedge angle, and other shape parameters, the oblique shock wave theory, normal shock wave relationship, or empirical formulas are used to make preliminary predictions of the nose bow shock wave separation distance, wing leading edge shock wave angle, and possible flow separation zone, thereby initially defining the main spatial range of shock wave generation and influence.

[0060] Secondly, rapid flow field diagnostic calculations are performed. Within the theoretically predicted range, or by establishing a relatively sparse initial grid for the entire computational domain, a first-order upwind scheme or other low-order scheme is used to perform rapid calculations of the steady flow field. After the calculations are completed, the gradient fields of key physical quantities in the flow field, such as the pressure gradient, are extracted. Field and density gradient field Set a gradient threshold empirically determined based on the flow field scale. , will satisfy or The physical region containing the grid cells is marked as a candidate local region where shock waves, contact sparses, or shear layers may exist. This region is the target area for subsequent high-resolution local mesh construction and fine shock wave capture.

[0061] Finally, mesh generation is guided. Using the candidate local regions identified in the diagnostic calculations above as the core, the mesh is locally refined to form a final computational mesh with high resolution in critical regions and moderate mesh density in smooth regions, for high-precision calculations.

[0062] In a preferred embodiment of the present invention, the spatial discretization of the gas dynamics equations describing the flow (i.e., the Navier-Stokes equations) is achieved using the finite difference method. The core idea of ​​this method is to directly discretize and approximate the spatial derivative terms in the governing equations at the grid nodes of the structured mesh. Further details regarding the finite difference method are not provided in this embodiment of the invention.

[0063] In step S2, for the three-dimensional gas dynamics problem, the Euler equation in its conserved form can be written as:

[0064] ;

[0065] Among them, the conserved variables It is a five-dimensional vector (mass, momentum components in three directions, and energy), where, For density, , and Three directions , ,and The speed on, Let t represent energy and T represent temperature. In the finite difference method, it is necessary to calculate the numerical flux of a physical flux (such as F) at the grid interface, where t represents time, x, y, and z represent the three spatial coordinate directions, and F, G, and H are the physical vectors in the x, y, and z directions, respectively. For the flux Jacobian matrix... Perform feature decomposition into form ,in:

[0066] It is an eigenvalue matrix. , representing the propagation speed of the five wave patterns, and These are left-moving sound waves and right-moving sound waves, respectively. Entropy wave, and This represents two shear waves;

[0067] R is a full-dimensional right eigenvector matrix, and the column vectors represent the shape of each wave pattern in physical space.

[0068] L is the full-dimensional left eigenvector matrix, which is the inverse matrix of R. Its row vectors are used to project physical quantities onto the feature space.

[0069] In step S3, for the three-dimensional Euler equation, the flux Jacobian matrix... Given a 5×5 matrix, after eigenvalue decomposition of the Jacobian matrix, the five eigenvalues ​​of the eigenvalue matrix are:

[0070] ;

[0071] Here, u is velocity, and c is the local speed of sound. The standard feature projection steps include:

[0072] Projection: Projecting the physical flux at the i-th grid point Projecting onto the feature space yields the feature variables. ;

[0073] ;

[0074] Reconstruction: Components of the five feature variables Perform nonlinear reconstruction;

[0075] Taking a 5-point template as an example, the template is: Where, the subscript i represents the index of the current grid cell, and i-2, i-1, i+1, and i+2 represent the indices of the second and first adjacent grid cells to the left and the first and second adjacent grid cells to the right of the current grid cell, respectively. F i , F i-2 , F i-1 , F i+1 , F i+2 These represent the physical flux at the corresponding grid points. These represent the feature variables projected onto the feature space at the corresponding grid points.

[0076] Smoothness index ;

[0077] in, Let k be the smoothness index of the k-th sub-template, where k is the index of the sub-template. This indicates the coordinates of the right boundary of the current mesh cell interface. It is the interpolation polynomial of the k-th sub-template, and r is the order of the reconstructed polynomial, r=3. The coordinates of the left boundary of the current grid cell interface. To represent the first derivative operator with respect to spatial coordinate x, where Δx is the width of the grid cell, dx is the integration variable, representing integration along the x-axis of the spatial coordinate system.

[0078] Then calculate the nonlinear weights:

[0079] ;

[0080] ;

[0081] Then use The reconstructed value is obtained by weighted averaging. .

[0082] Back projection: Transforming the reconstructed feature space flux back into physical space.

[0083] ;

[0084] The five eigenvalues ​​are categorized according to their physical meaning and reconstructed in a low-dimensional subspace, thereby significantly reducing the amount of computation.

[0085] The five eigenvalues ​​are divided into two main categories:

[0086] Convection mode: transport of fluid, entropy, and eddies, propagating at the local fluid velocity u, including entropy waves. and two shear waves , .

[0087] The first dimensionality-reduced eigenvector can be represented as... ,in This represents the eigenvectors of the fluid, entropy, and vortex corresponding to the full-dimensional left eigenvector matrix L.

[0088] Acoustic mode: Pressure disturbance, propagating at velocities u±c, is the primary carrier of shock waves and expansion waves, and includes left-lateral acoustic waves. and right-handed sound waves .

[0089] The second dimension-reduced feature vector can be represented as:

[0090] ;

[0091] θ is the local covariant velocity, and c is the speed of sound. and These are the eigenvectors of the left and right sound waves in the full-dimensional left eigenvector matrix L, respectively.

[0092] Dimensionality reduction strategy: In five-dimensional projection, the five feature patterns need to be processed separately, and the left feature vector is a 5×5 matrix, requiring 25 multiplication and addition operations during the projection process. However, this method classifies the feature patterns into two categories, reducing the left feature variable from 5 dimensions to 2 dimensions, and the left feature variable becomes a 2×5 matrix, requiring only 10 multiplication and addition operations during the projection process.

[0093] The physical classification mechanism is based on wave patterns in fluids. The five-dimensional feature space corresponds to five independent feature patterns: a left-moving acoustic wave λ1, an entropy wave λ2, two shear waves λ3 and λ4, and a right-moving acoustic wave λ5. The physical mechanism of dimensionality reduction lies in aggregating physically closely related patterns: aggression waves and shear waves λ2, λ3, and λ4, which propagate at the same fluid velocity, are aggregated into a convection pattern; left and right acoustic waves λ1 and λ5, which propagate in opposite directions but have the same physical essence, are aggregated into an acoustic pattern. In this way, the complex dynamics of the five-dimensional space are approximately described by two dominant physical processes (convective transport and acoustic wave propagation), allowing for efficient reconstruction in the lower-dimensional subspace without losing key physical properties.

[0094] In step S4, L red Two representative vectors and The composition is a 2×5 matrix, which is projected to obtain... The projection process is as follows:

[0095] ;

[0096] right Reconstruction is performed to obtain nonlinear weighting factors. ;right Reconstruction is performed to obtain nonlinear weighting factors. Among them, for During reconstruction, the smoothness index of the convection mode is calculated by substituting it into the smoothness index formula. ;right During reconstruction, the smoothness index of the acoustic mode is calculated by substituting it into the smoothness index formula. .

[0097] In step S5, the two sets of nonlinear weighting factors obtained above are... and Expanded into a five-dimensional diagonal weight matrix:

[0098] .

[0099] In step S5, the nonlinear weighting factor ω k,a and ω k,c The expansion to the full-dimensional feature space is as follows: First, a 5×5 diagonal weight matrix is ​​constructed for each sub-template k (k=0,1,2). The diagonal elements of this matrix are filled with two sets of weights according to a specific rule: the first and fifth diagonal elements (corresponding to the acoustic modes) are filled with... The second, third, and fourth diagonal elements (corresponding to the flow pattern) are filled with... Then, using this weight matrix and the original full-dimensional left eigenvector matrix L and full-dimensional right eigenvector matrix R, the synthesis formula is obtained, and the numerical flux is calculated directly in the physical space, avoiding the need to perform complete projection and back-projection operations in the full-dimensional feature space.

[0100] This method spatially discretizes and splits the gas dynamics equations (Navier-Stokes equations); calculates the eigenvalues ​​and eigenvectors of the Jacobian matrix; the key is to construct a two-dimensional reduced-dimensional feature subspace composed of convective mode eigenvectors and acoustic mode eigenvectors; the physical flux is projected onto this two-dimensional subspace and nonlinearly reconstructed to obtain two sets of reconstruction weights; the weights are extended to full-dimensional space to construct a diagonal weight matrix; finally, the numerical flux in the physical space is calculated using the synthesis formula. The advantages of this invention are: by reducing the feature projection from five dimensions to two dimensions, the number of nonlinear reconstructions and matrix operations is significantly reduced, improving computational efficiency by 20%~50%, while perfectly maintaining the shock wave capture capability, high accuracy, and oscillation-free characteristics of the original format; it is highly versatile and easy to integrate. This method is particularly suitable for the refined simulation of complex flow structures such as shock waves, expansion waves, and shear layers in supersonic aircraft, and can accurately analyze the aerodynamic forces and aerothermal loads acting on the aircraft surface, providing reliable data support for the refined design of aircraft.

[0101] In this embodiment of the invention, the dimensionality-reduced feature subspace is composed of a first dimensionality-reduced feature vector l. c1 Second dimensionality reduction eigenvector l c2 Composition, the first dimensionality-reduced eigenvector l c1 The second dimensionality-reduced eigenvector, l, is a linear combination of eigenvectors related to entropy waves and shear waves, representing the convection process. c2 The acoustic process is represented by a linear combination of eigenvectors related to the left and right sound waves. The first reduced-dimensional eigenvector l c1 It is composed of a linear combination of eigenvectors related to entropy waves and shear waves.

[0102] In this embodiment of the invention, the second dimensionality reduction feature vector l c2 The construction method is as follows:

[0103] ;

[0104] Where θ is the local covariant velocity, c is the speed of sound, and l1 and l5 are the eigenvectors of the left and right sound waves in the full-dimensional left eigenvector matrix L, respectively.

[0105] In this embodiment of the invention, the nonlinear reconstruction format is WENO, TENO, or WENN format.

[0106] In this embodiment of the invention, the expression for the full-dimensional diagonal weight matrix is:

[0107] ;

[0108] in, It is a diagonal weight matrix. It is an identity matrix.

[0109] In this embodiment of the invention, the synthesis formula is specifically as follows:

[0110] ;

[0111] in, Let q be the numerical flux in physical space. k Let r1 and r5 be the physical flux reconstruction values ​​on the k-th sub-template, r1 and r5 be the column vectors corresponding to the left and right acoustic waves in the full-dimensional right eigenvector matrix R, l1 and l5 be the eigenvectors corresponding to the left and right acoustic waves in the full-dimensional left eigenvector matrix L, i represent the i-th grid point, and represents the grid point index.

[0112] The effect of this method is as follows: Figures 2-3 , and as shown in Table 1.

[0113] Table 1 compares the computation time of the dimension reduction method of the present invention applied to three traditional numerical calculation methods (WENO, TENO and WENN) on the Lax shock tube problem.

[0114] Table 1 compares the computation time of the dimension reduction method and the traditional method for solving the Lax shock tube problem;

[0115] Traditional methods 323s 1046s 476s Dimensionality reduction method of the present invention 237s 472s 238s Improve efficiency 26.6% 54.9% 50.0%

[0116] As can be seen from Table 1, for the three numerical methods, the dimensionality reduction method of the present invention can significantly reduce the computation time compared with the traditional calculation method.

[0117] Figure 2 This diagram illustrates a comparison of the density distribution obtained from simulating the Lax shock tube problem using the dimensionality reduction method of this invention and traditional numerical calculation methods; from Figure 2It can be seen that the dimensionality reduction method of this invention and the traditional method produce consistent results in simulating the Lax shock tube, indicating that the accuracy of the method of this invention is not reduced. Specifically, to facilitate observation of the shock region, Figure 2 The shock wave region was magnified and placed... Figure 2 The middle part is used as a subgraph, and corresponding coordinate axes are added, that is... Figure 2 The sub-image in the image is a magnified view of a local area of ​​the shock wave region.

[0118] Figure 3 The diagram shows a comparison of density cloud maps obtained by simulating a two-dimensional Riemann problem using the traditional method (left half) and the method of the present invention (right half). It can be seen that the dimensionality reduction method of the present invention and the traditional method produce the same results in simulating the Riemann problem.

[0119] Based on the data in the table above and Figures 2-3 The density distribution of the present invention and the traditional method can maintain simulation accuracy and improve computational efficiency.

[0120] Example 2;

[0121] Based on Embodiment 1, Embodiment 2 of the present invention provides a shock wave capture system for aerodynamic characteristic analysis of aircraft, the system comprising:

[0122] The modeling and discretization unit is used to determine local regions that may generate shock waves or flow discontinuities based on the incoming flow conditions and shape parameters of the aircraft; to obtain the geometric parameters of the local regions and construct the corresponding local computational grid; and to spatially discretize the gas dynamics control equations describing the flow on the local computational grid to obtain a set of discretized equations, wherein the convective flux term in the set of discretized equations represents the interface numerical flux to be determined.

[0123] The feature analysis and transformation preparation unit is used to perform eigenvalue decomposition on the local flow Jacobian matrix on the grid cell interface where flux reconstruction is required, according to the requirements of the shock wave capture format, to obtain the eigenvalues ​​of the local flow Jacobian matrix, the full-dimensional left eigenvector matrix L, and the full-dimensional right eigenvector matrix R; wherein, the full-dimensional left eigenvector matrix L is used to project the variables in the physical space to the feature space, and the full-dimensional right eigenvector matrix R is used to inversely transform the variables in the feature space back to the physical space;

[0124] The dimensionality reduction eigenspace construction unit is used to construct the first dimensionality reduction eigenvector l using the eigenvectors in the full-dimensional left eigenvector matrix L. c1 Second dimensionality reduction eigenvector l c2 , where l c1 It is a characteristic vector representing the convection mode, composed of the characteristic vectors corresponding to entropy waves and shear waves. c2It is a weighted combination of the corresponding eigenvectors in the full-dimensional left eigenvector matrix L and the full-dimensional right eigenvector matrix R, representing the eigenvectors of the acoustic mode. c1 and l c2 This forms a two-dimensional reduced feature subspace, resulting in a reduced projection matrix L. red =[l c1 ;l c2 ];

[0125] The dimensionality reduction space reconstruction and weight extraction unit is used to convert the physical flux F at the grid cell interface into the dimensionality reduction projection matrix L. red Projecting onto a two-dimensional reduced feature subspace yields two reduced feature variables F. a and F c For the characteristic variable F a and F c Interpolation is performed using a nonlinear reconstruction scheme, and two sets of nonlinear weighting factors ω, corresponding to the convection mode and the acoustic mode respectively, are extracted during the reconstruction process. k,a and ω k,c And obtain the reconstructed feature variable values ​​on both sides of the interface in the two-dimensional reduced feature subspace;

[0126] The interface numerical flux synthesis unit is used to synthesize the nonlinear weighting factor ω. k,a and ω k,c Extending to the full-dimensional feature space, a full-dimensional diagonal weight matrix is ​​constructed; in the full-dimensional diagonal weight matrix, the values ​​corresponding to the first and fifth diagonal positions of the convection mode are assigned as ω. k,a The values ​​ω are assigned to the second to fourth diagonal positions corresponding to the acoustic mode. k,c If the full-dimensional diagonal weight matrix has other diagonal positions, then it is assigned a value of 1; based on the physical flux reconstruction value on each sub-template, the synthesis formula is derived by decomposing the full-dimensional diagonal weight matrix to calculate the numerical flux in the physical space.

[0127] The time-progression and flow field solution unit is used to process all the mesh cell interfaces that need to be processed by the feature analysis and transformation preparation unit, the dimension reduction feature subspace construction unit, the dimension reduction space reconstruction and weight extraction unit, and the interface numerical flux synthesis unit, and calculate the final physical interface numerical flux of all mesh cell interfaces in the entire field. The final physical interface numerical flux of all mesh cell interfaces in the entire field is substituted into the discretized equation system, and time integration is performed to advance the flow field evolution until a convergent or transient flow field solution is obtained, thereby obtaining the shock wave structure or aerodynamic / thermal load distribution of the aircraft.

[0128] In this embodiment of the invention, the two-dimensional reduced feature subspace is composed of a first reduced feature vector l. c1 Second dimensionality reduction eigenvector lc2 Composition, the first dimensionality-reduced eigenvector l c1 The second dimensionality-reduced eigenvector l is formed by a linear combination of eigenvectors related to entropy waves and shear waves. c2 It is composed of a linear combination of eigenvectors related to the left and right sound waves.

[0129] In this embodiment of the invention, the second dimensionality reduction feature vector l c2 The construction method is as follows:

[0130] ;

[0131] Where θ is the local covariant velocity, c is the speed of sound, and l1 and l5 are the row vectors corresponding to the left and right sound waves in the full-dimensional left eigenvector matrix L, respectively.

[0132] In this embodiment of the invention, the nonlinear reconstruction is based on the WENO, TENO, or WENN format.

[0133] In this embodiment of the invention, the system may further include an input unit for inputting data and a display unit for displaying shock wave structure or aerodynamic / thermal load distribution data of the aircraft.

[0134] Although preferred embodiments of the invention have been described, those skilled in the art, upon learning the basic inventive concept, can make other changes and modifications to these embodiments. Therefore, the appended claims are intended to be interpreted as including both the preferred embodiments and all changes and modifications falling within the scope of the invention.

[0135] Obviously, those skilled in the art can make various modifications and variations to this invention without departing from its spirit and scope. Therefore, if these modifications and variations fall within the scope of the claims of this invention and their equivalents, this invention also intends to include these modifications and variations.

Claims

1. A shock wave capture method for aerodynamic characteristic analysis of aircraft, characterized in that, The method includes: Step S1: Based on the incoming flow conditions and shape parameters of the aircraft, determine the local regions that may generate shock waves or flow discontinuities; obtain the geometric parameters of the local regions and construct the corresponding local computational grid; on the local computational grid, spatially discretize the gas dynamics control equations describing the flow to obtain a discretized set of equations, wherein the convective flux term in the discretized set of equations represents the interface numerical flux to be determined. Step S2: According to the requirements of the shock wave capture format, on the interface of the grid cell where flux reconstruction is required, perform eigenvalue decomposition on the local flow Jacobian matrix to obtain the eigenvalues ​​of the local flow Jacobian matrix, the full-dimensional left eigenvector matrix L, and the full-dimensional right eigenvector matrix R; wherein, the full-dimensional left eigenvector matrix L is used to project the variables in the physical space to the eigenspace, and the full-dimensional right eigenvector matrix R is used to inversely transform the variables in the eigenspace back to the physical space; Step S3: Construct the first dimension-reduced eigenvector l using the eigenvectors in the full-dimensional left eigenvector matrix L. c1 Second dimensionality reduction eigenvector l c2 , where l c1 It is a characteristic vector representing the convection mode, composed of the characteristic vectors corresponding to entropy waves and shear waves. c2 It is a weighted combination of the corresponding eigenvectors in the full-dimensional left eigenvector matrix L and the full-dimensional right eigenvector matrix R, representing the eigenvectors of the acoustic mode. c1 and l c2 This forms a two-dimensional reduced feature subspace, resulting in a reduced projection matrix L. red =[l c1 ;l c2 ]; Step S4: The physical flux F at the grid cell interface is transformed by the dimension reduction projection matrix L. red Projecting onto a two-dimensional reduced feature subspace yields two reduced feature variables F. a and F c For the characteristic variable F a and F c Interpolation is performed using a nonlinear reconstruction scheme, and two sets of nonlinear weighting factors ω, corresponding to the convection mode and the acoustic mode respectively, are extracted during the reconstruction process. k,a and ω k,c And obtain the reconstructed feature variable values ​​on both sides of the interface in the two-dimensional reduced feature subspace; Step S5: Adjust the nonlinear weighting factor ω k,a and ω k,c Extending to the full-dimensional feature space, a full-dimensional diagonal weight matrix is ​​constructed; in the full-dimensional diagonal weight matrix, the values ​​corresponding to the first and fifth diagonal positions of the convection mode are assigned as ω. k,a The values ​​ω are assigned to the second to fourth diagonal positions corresponding to the acoustic mode. k,c If the full-dimensional diagonal weight matrix has other diagonal positions, then it is assigned a value of 1; based on the physical flux reconstruction value on each sub-template, the synthesis formula is derived by decomposing the full-dimensional diagonal weight matrix to calculate the numerical flux in the physical space. Step S6: Repeat steps S2 to S5 for all grid cell interfaces that need to be processed until the final physical interface numerical flux of all grid cell interfaces in the entire field is calculated; substitute the final physical interface numerical flux of all grid cell interfaces in the entire field into the discretized equation set, perform time integration, and advance the flow field evolution until a convergent or transient flow field solution is obtained, thereby obtaining the shock wave structure or aerodynamic / thermal load distribution of the aircraft.

2. The shock wave capture method for aerodynamic characteristic analysis of aircraft according to claim 1, characterized in that, The two-dimensional reduced eigenspace is composed of the first reduced eigenvector l c1 Second dimensionality reduction eigenvector l c2 Composition, the first dimensionality-reduced eigenvector l c1 The second dimensionality-reduced eigenvector l is formed by a linear combination of eigenvectors related to entropy waves and shear waves. c2 It is composed of a linear combination of eigenvectors related to the left and right sound waves.

3. The shock wave capture method for aerodynamic characteristic analysis of aircraft according to claim 2, characterized in that, Second dimensionality reduction eigenvector l c2 The construction method is as follows: ； Where θ is the local covariant velocity, c is the speed of sound, and l1 and l5 are the row vectors corresponding to the left and right sound waves in the full-dimensional left eigenvector matrix L, respectively.

4. The shock wave capture method for aerodynamic characteristic analysis of aircraft according to claim 1, characterized in that, The nonlinear reconstruction format is WENO, TENO, or WENN.

5. The shock wave capture method for aerodynamic characteristic analysis of aircraft according to claim 1, characterized in that, The expression for the full-dimensional diagonal weight matrix is: ； in, It is a diagonal weight matrix. It is an identity matrix.

6. The shock wave capture method for aerodynamic characteristic analysis of aircraft according to claim 1, characterized in that, The specific synthesis formula is as follows: ； in, Let q be the numerical flux in physical space. k Let be the reconstructed value of the physical flux F on the k-th sub-template, r1 and r5 be the column vectors corresponding to the left and right sound waves in the full-dimensional right eigenvector matrix R, l1 and l5 be the row vectors corresponding to the left and right sound waves in the full-dimensional left eigenvector matrix L, i represent the i-th grid point, and represent the grid point index.

7. A shock wave capture system for aerodynamic characteristic analysis of aircraft, characterized in that, The system includes: The modeling and discretization unit is used to determine local regions that may generate shock waves or flow discontinuities based on the incoming flow conditions and shape parameters of the aircraft; to obtain the geometric parameters of the local regions and construct the corresponding local computational grid; and to spatially discretize the gas dynamics control equations describing the flow on the local computational grid to obtain a set of discretized equations, wherein the convective flux term in the set of discretized equations represents the interface numerical flux to be determined. The feature analysis and transformation preparation unit is used to perform eigenvalue decomposition on the local flow Jacobian matrix on the grid cell interface where flux reconstruction is required, according to the requirements of the shock wave capture format, to obtain the eigenvalues ​​of the local flow Jacobian matrix, the full-dimensional left eigenvector matrix L, and the full-dimensional right eigenvector matrix R; wherein, the full-dimensional left eigenvector matrix L is used to project the variables in the physical space to the feature space, and the full-dimensional right eigenvector matrix R is used to inversely transform the variables in the feature space back to the physical space; The dimensionality reduction eigenspace construction unit is used to construct the first dimensionality reduction eigenvector l using the eigenvectors in the full-dimensional left eigenvector matrix L. c1 Second dimensionality reduction eigenvector l c2 , where l c1 It is a characteristic vector representing the convection mode, composed of the characteristic vectors corresponding to entropy waves and shear waves. c2 It is a weighted combination of the corresponding eigenvectors in the full-dimensional left eigenvector matrix L and the full-dimensional right eigenvector matrix R, representing the eigenvectors of the acoustic mode. c1 and l c2 This forms a two-dimensional reduced feature subspace, resulting in a reduced projection matrix L. red =[l c1 ;l c2 ]; The dimensionality reduction space reconstruction and weight extraction unit is used to convert the physical flux F at the grid cell interface into the dimensionality reduction projection matrix L. red Projecting onto a two-dimensional reduced feature subspace yields two reduced feature variables F. a and F c For the characteristic variable F a and F c Interpolation is performed using a nonlinear reconstruction scheme, and two sets of nonlinear weighting factors ω, corresponding to the convection mode and the acoustic mode respectively, are extracted during the reconstruction process. k,a and ω k,c And obtain the reconstructed feature variable values ​​on both sides of the interface in the two-dimensional reduced feature subspace; The interface numerical flux synthesis unit is used to synthesize the nonlinear weighting factor ω. k,a and ω k,c Extending to the full-dimensional feature space, a full-dimensional diagonal weight matrix is ​​constructed; in the full-dimensional diagonal weight matrix, the values ​​corresponding to the first and fifth diagonal positions of the convection mode are assigned as ω. k,a The values ​​ω are assigned to the second to fourth diagonal positions corresponding to the acoustic mode. k,c If the full-dimensional diagonal weight matrix has other diagonal positions, then it is assigned a value of 1; based on the physical flux reconstruction value on each sub-template, the synthesis formula is derived by decomposing the full-dimensional diagonal weight matrix to calculate the numerical flux in the physical space. The time-progression and flow field solution unit is used to process all the mesh cell interfaces that need to be processed by the feature analysis and transformation preparation unit, the dimension reduction feature subspace construction unit, the dimension reduction space reconstruction and weight extraction unit, and the interface numerical flux synthesis unit, and calculate the final physical interface numerical flux of all mesh cell interfaces in the entire field. The final physical interface numerical flux of all mesh cell interfaces in the entire field is substituted into the discretized equation system, and time integration is performed to advance the flow field evolution until a convergent or transient flow field solution is obtained, thereby obtaining the shock wave structure or aerodynamic / thermal load distribution of the aircraft.

8. A shock wave capture system for aerodynamic characteristic analysis of an aircraft according to claim 7, characterized in that, The two-dimensional reduced eigenspace is composed of the first reduced eigenvector l c1 Second dimensionality reduction eigenvector l c2 Composition, the first dimensionality-reduced eigenvector l c1 The second dimensionality-reduced eigenvector l is formed by a linear combination of eigenvectors related to entropy waves and shear waves. c2 It is composed of a linear combination of eigenvectors related to the left and right sound waves.

9. A shock wave capture system for aerodynamic characteristic analysis of an aircraft according to claim 8, characterized in that, Second dimensionality reduction eigenvector l c2 The construction method is as follows: ； Where θ is the local covariant velocity, c is the speed of sound, and l1 and l5 are the row vectors corresponding to the left and right sound waves in the full-dimensional left eigenvector matrix L, respectively.

10. A shock wave capture system for aerodynamic characteristic analysis of an aircraft according to claim 7, characterized in that, The nonlinear reconstruction is based on the WENO, TENO, or WENN format.

Images