Enhanced moving least squares meshfree-newmark integration method for reducing total dispersion error of transient wave propagation

Abstract

The application discloses an enhanced mobile least square meshless-Newmark integral method for reducing total dispersion error of transient wave propagation, and belongs to the technical field of numerical analysis and total dispersion error prediction of transient wave propagation, and comprises the following contents: a wave propagation problem domain is discretized into a series of nodes; a standard mobile least square interpolation base function of each node is derived; a local numerical approximation of an enhanced least square algorithm is completed; a global interpolation format of a general scalar field based on the enhanced least square algorithm is constructed; a system matrix equation of the dispersion analysis of the transient wave propagation problem is obtained by adopting a Galerkin weighted residual method and a Newmark time integral method; and based on a given node distribution, the total dispersion error of different calculation methods along different wave propagation directions under different time integral step lengths is calculated. The application adopts the above method, reduces the total dispersion error caused by the coupling of spatial dispersion and time dispersion, and improves the solving precision.

Description

Technical Field

[0001] This invention relates to the technical field of transient wave propagation numerical analysis and total dispersion error prediction, and in particular to an enhanced moving least squares meshless-Newmark integral method for reducing the total dispersion error of transient wave propagation. Background Technology

[0002] Transient wave propagation is a core topic in wave dynamics and is widely used in fields such as mechanical vibration, underwater acoustics, seismic exploration, medical ultrasound, and sonar technology. Analytical solutions are only applicable to simplified models with homogeneous media and regular boundaries, and are difficult to handle real-world complex scenarios such as faults, folds, and human tissue. Numerical simulation has become the mainstream method for solving such problems.

[0003] Numerical methods transform wave partial differential equations into a system of algebraic equations by discretizing space and time, enabling flexible handling of complex geometries, inhomogeneous media, and nonlinear problems. However, the traditional finite element method has significant drawbacks: to ensure accuracy, at least six elements are required for each wavelength, and high-frequency wave simulation leads to a surge in the number of elements, resulting in high computational costs and long processing times, making it difficult to meet engineering requirements such as parameter inversion and multiphysics coupling.

[0004] Spatial dispersion error is a core problem in numerical simulation. Discretization makes the system "overly rigid," causing numerical wavenumbers and wave velocities to deviate from their true values, with errors amplified sharply in the high-frequency band, affecting wavefront arrival time, amplitude, and waveform prediction accuracy. Existing methods optimize stiffness through higher-order shape functions and reduced-order integrals, but a balance between accuracy and efficiency remains a challenge.

[0005] Meshless methods rely solely on discrete nodes to construct approximate functions, eliminating the need for mesh generation. This avoids mesh distortion, reduces system stiffness, and improves dispersion characteristics. Among these, the moving least squares meshless method is widely used, achieving 40%–60% lower spatial dispersion error compared to the linear finite element method at the same node density. However, this method is still affected by factors such as node distribution, support domain, and weighting function, and cannot completely eliminate spatial dispersion error. Further optimization of interpolation and discretization strategies is needed.

[0006] The numerical accuracy of transient wave propagation is determined by both spatial and temporal discretization. The recursive scheme used for time integration introduces truncation errors that accumulate with the step size, easily leading to numerical instability. Furthermore, spatial and temporal errors are strongly coupled, which is particularly prominent in broadband transient problems. Existing meshless moving least squares methods still have shortcomings in controlling total dispersion error and suppressing numerical anisotropy. Summary of the Invention

[0007] The purpose of this invention is to provide an enhanced moving least squares meshless Newmark integration method to reduce the total dispersion error of transient wave propagation and improve the solution accuracy.

[0008] To achieve the above objectives, this invention provides an enhanced moving least squares meshless-Newmark integration method for reducing the total dispersion error of transient wave propagation, comprising the following steps: S1. Discretize the wave propagation problem domain into a series of nodes, and define an influence domain for each node; S2, selection weight function and interpolation basis function, based on the moving least squares method to fit the local approximate interpolation function, the standard moving least squares interpolation basis function of each node is derived; S3. Select the Lagrange polynomial as the local reinforcement interpolation basis function, construct the node degrees of freedom at each node consisting of the local reinforcement interpolation basis function and the coefficients to be determined, and complete the local numerical approximation of the enhanced moving least squares meshless algorithm. S4. Based on the node interpolation basis functions of the standard moving least squares meshless algorithm obtained in step S2 and the enhanced local numerical approximation obtained in step S3, a global interpolation scheme for a general scalar field based on the enhanced least squares algorithm is constructed under the premise of satisfying the unit decomposition property. S5. Based on the global interpolation scheme of the general scalar field obtained in step S4, based on the enhanced least squares algorithm, and based on the partial differential control equation of the wave propagation problem, the Galerkin weighted residual method and Newmark time integration method are used to complete the time discretization and spatial discretization of the transient wave propagation problem, and the system matrix equation of the dispersion analysis of the transient wave propagation problem is obtained. S6. Based on the given node distribution, solve the system matrix equation of the transient wave propagation problem obtained in step S5, and calculate the total dispersion error of different calculation methods along different wave propagation directions under different time integration step sizes. S7. Based on the total dispersion error of the transient wave propagation problem in different wave propagation directions under different time integration step sizes obtained in step S6, plot the dispersion error as a function of the dimensionless wave number, compare the dispersion error results obtained by different calculation methods, analyze the coupling relationship between spatial dispersion error, time dispersion error and total dispersion error in the transient wave propagation problem, and complete the dispersion analysis of the wave propagation problem.

[0009] Preferably, in step S1, the influence domain is circular or rectangular, and the size of the influence domain is adapted to the average spacing of the node distribution.

[0010] Preferably, the support domain of the integration point is a circular or rectangular region centered on the integration point, and the number of nodes participating in the interpolation is determined by the support domain of the integration point during numerical integration. The node influence domain is the area that each node can influence during numerical integration. The size of the node influence domain is closely related to the integration point support domain.

[0011] Preferably, step S2 specifically includes the following steps: S21. The weight function of nodes within the influence domain is non-zero, and the weight function of nodes outside the influence domain is zero. The weight functions are shown below: ; in, The radius of the influence region; S22. For the transient wave propagation problem in a scalar field, based on the selected nodes within the influence domain and the standard moving least squares meshless technique, the initial form of the global interpolation scheme is constructed as follows: ; in, For the physical variables of transient wave propagation problems, This represents the number of nodes participating in the interpolation. The position vector at the interpolation point, vector express There are unknown coefficients to be determined, and the vector is... This represents a vector composed of the interpolation basis functions used; S23. For one-dimensional, two-dimensional, and three-dimensional problems, the Lagrange linear basis functions used are expressed as follows: ; S24. Constructing the weighted squared error: ; in, The total number of nodes contained in the supporting domain for the node at coordinate x. The function value at the node with coordinate x; To minimize the weighted squared error, the weighted squared error must be about The extreme points of the weighted squared difference function with respect to The reciprocal of is zero: ; Therefore, we can conclude that: ; in, This represents an unknown node variable corresponding to n nodes in a supporting domain, represented by a matrix. Specifically, it is expressed as follows: ; matrix Specifically, it is expressed as follows: ; Undetermined coefficient matrix for: ; matrix , and Substituting the initial form of the global interpolation format, we get: ; in, For each node constructed, use the standard moving least squares interpolation basis function.

[0012] Preferably, step S3 specifically includes the following steps: S31. Enhanced interpolation scheme built on enhanced least squares algorithm: ; in, To increase the number of nodes involved in the least squares interpolation algorithm, To enhance the local numerical approximation of the least squares algorithm, and For each node, construct a local numerical approximation corresponding to the local enhancement interpolation basis function and the coefficients to be determined; S32. A local enhancement interpolation basis function is constructed using Lagrange polynomial basis functions and dimensionless relative coordinate values. The formula for the local enhancement interpolation basis function is: ; in, , and These are dimensionless relative coordinate values. This represents the average spacing between nodes.

[0013] Preferably, the global interpolation format in the enhanced least squares algorithm obtained in step S4 is: .

[0014] Preferably, step S5 specifically includes the following steps: S51. Based on the principle of virtual displacement, the partial differential control equations of the transient sound propagation problem are transformed into integral form, as shown below: ; in, For the Laplace operator, For the speed of sound propagation, For sound pressure, Indicates sound pressure Regarding time The second derivative; S52. Based on integration by parts and the Gaussian divergence theorem, the integral form of the partial differential control equations for the transient sound propagation problem is expressed as: ; in, The unit vector of the outward normal to the boundary of the computational domain along the sound propagation space; S53. Using the Galerkin weighted residual method, the sound pressure interpolation scheme for the transient acoustic scattering problem of underwater targets based on standard quadrilateral element grids is substituted into the formula in S52, resulting in the matrix form of the integral equation: ; in, for The sound pressure level of the unknown node corresponding to each node. The mass matrix for the transient acoustic scattering problem is... For the reason A vector composed of interpolation shape functions at each node. For the interpolation shape function corresponding to the node, This represents the total number of nodes within the spatial computational domain. Let be the smooth stiffness matrix for the transient acoustic scattering problem. For the transient sound scattering problem, the smoothed sound pressure gradient matrix is... The force matrix for a transient acoustic scattering problem. Let be the normal velocity vector of the acoustic particles at the nodes along the boundary of the spatial computational domain; S54. Based on the global interpolation scheme in the enhanced least squares algorithm obtained in step S4, the matrix governing equation for the transient wave propagation problem is written as follows: ; The solution to the matrix governing equation is: ; in, The imaginary unit, The wavenumber is calculated using numerical methods. Let be the unit vector in the direction of wave propagation. Let be the position vector of the field point along the wave propagation direction. The angular frequency is calculated using numerical methods. To and The magnitude vector corresponding to the coefficients of each unknown node. For time variation in numerical computation; S55. In the Newmark time integration method, the sound pressure vector at each node... , and its first and second derivative vectors with respect to time and Represented as: ; in, It is the time integration step size. and These are the parameters of the numerical damping in the Newmark time-integration algorithm under two conditions.

[0015] Preferably, in step S6, the total numerical dispersion error of the numerical solution to the transient wave propagation problem is: ; in, For accurate sound propagation speed, For accurate wavenumbers, For accurate angular frequency, The velocity of sound is obtained through numerical calculation.

[0016] Therefore, the present invention employs the enhanced moving least squares meshless-Newmark integration method described above for reducing the total dispersion error of transient wave propagation, which has the following advantages: (1) In this invention, high-order interpolation can be achieved without increasing the number of nodes. Compared with the standard moving least squares meshless method, it can significantly reduce the total dispersion error of transient wave propagation problems and improve the accuracy of numerical calculation.

[0017] (2) In this invention, the “numerical anisotropy” phenomenon in the standard moving least squares meshless method can be effectively suppressed, making the calculation results in different wave propagation directions more uniform and stable.

[0018] (3) In this invention, the total dispersion error will continue to decrease as the time integration step size decreases. The calculation accuracy can be flexibly controlled by adjusting the step size to adapt to transient wave propagation simulation analysis with different accuracy requirements in engineering.

[0019] The technical solution of the present invention will be further described in detail below with reference to the accompanying drawings and embodiments. Attached Figure Description

[0020] Figure 1 This is a schematic diagram of the enhanced moving least squares meshless-Newmark integral method for reducing the total dispersion error of transient wave propagation according to the present invention. Figure 2 This is a schematic diagram of the support domain of the integration point in the enhanced moving least squares meshless algorithm in this embodiment of the invention; Figure 3 This is a schematic diagram of the field influence domain in the enhanced moving least squares meshless algorithm in this embodiment of the invention; Figure 4 A schematic diagram of field point distribution in the standard moving least squares meshless algorithm and the enhanced moving least squares meshless algorithm for transient wave propagation dispersion analysis provided in embodiments of the present invention; Figure 5 The curves showing the variation of dispersion error with dimensionless wavenumber in different wave propagation directions, calculated by the standard moving least squares meshless-Newmark time integration algorithm provided in this embodiment of the invention. Figure 5(a) indicates the direction of wave propagation. The curve of change, Figure 5 (b) is the direction of wave propagation. The curve of change, Figure 5 (c) indicates the direction of wave propagation. The curve of change, Figure 5 (d) indicates the direction of wave propagation. The curve of change; Figure 6 The curves showing the variation of the total dispersion error with the dimensionless wavenumber in different wave propagation directions, calculated by the enhanced moving least squares meshless-Newmark time integration algorithm provided in this embodiment of the invention, are shown. Figure 6 (a) indicates the direction of wave propagation. The curve of change, Figure 6 (b) is the direction of wave propagation. The curve of change, Figure 6 (c) indicates the direction of wave propagation. The curve of change, Figure 6 (d) indicates the direction of wave propagation. The curve showing the change. Detailed Implementation

[0021] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, not all embodiments. The components of the embodiments of the present invention described and shown in the accompanying drawings can generally be arranged and designed in various different configurations. Specific model specifications need to be selected and determined according to the actual specifications of the device, etc. The specific selection calculation method adopts existing technology in the art, and therefore will not be described in detail.

[0022] Example like Figure 1 As shown, this invention provides an enhanced moving least squares meshless-Newmark integration method to reduce the total dispersion error of transient wave propagation, comprising the following steps: S1. Discretize the wave propagation problem domain into a series of nodes, and define an influence region for each node. The influence region is circular or rectangular (e.g., ...). Figure 2 and Figure 3 The size of the influence domain is adapted to the average spacing of the node distribution; in this embodiment, it is circular. S2, selection weight function and interpolation basis function, based on the moving least squares method to fit the local approximate interpolation function, the standard moving least squares interpolation basis function of each node is derived; The weight function value of nodes within this region is non-zero, while the weight function value of nodes outside this region is zero. A weight function of the following form can generally be used: (1); in, This indicates the radius of the influence region.

[0023] For any scalar field, numerical interpolation is performed using the standard moving least squares meshless technique. It can be represented as: (2); in, The vector represents the position vector at the interpolation point. This represents m unknown coefficients to be determined, and m is also the number of field points involved in the interpolation. The number of field points involved in the interpolation is determined by the support region of the integration point during numerical integration. The support region of the integration point is defined as a circular or rectangular region centered on the integration point (e.g., ...). Figure 3 As shown), this invention uses a rectangular region, with the center of this rectangular support region being the integration point. The distance from the center of the rectangular region to its side length is h (where h is the average spacing of the field points). The area that each field point can influence during numerical integration is called the field point influence domain. The size of the field point influence domain is closely related to the integration point support domain. In this invention, the influence domain of each field point is also a rectangular region, with the center of this rectangular field point influence domain being the field point itself. The distance from the center of the rectangular field point influence domain to its side length is h. (Vector) This represents a vector composed of the interpolation basis functions used.

[0024] For one-dimensional, two-dimensional, and three-dimensional problems, the Lagrange linear basis functions used can be expressed as: (3); The weighted squared error format is constructed as follows: (4); In the formula, n represents the total number of nodes contained in the supporting domain of the node at coordinate x. This represents the function value at the node with coordinate x.

[0025] In order to minimize the error scheme defined in formula (4), the weighted squared error must be about The extreme point, namely: (5); Therefore, we can conclude that: (6); in, This represents an unknown node variable that supports n nodes in the domain. and The related matrix has the following form: (7); (8); This embodiment uses a two-dimensional problem as an example to illustrate the matrix. and For the specific expression, the Lagrange linear basis function used for interpolation in the two-dimensional problem is: The result obtained in this way The matrix is ​​a 3x3 square matrix with the following form: (9); The results obtained The matrix is ​​a 3×n matrix and has the following form: (10); Solving equation (6) yields the matrix of undetermined coefficients. : (11); matrix , and Substituting into equation (2), we get: (12); in, It is the node interpolation function constructed based on the standard moving least squares algorithm.

[0026] S3. Select the Lagrange polynomial as the local reinforcement interpolation basis function, construct the node degrees of freedom at each node consisting of the local reinforcement interpolation basis function and the coefficients to be determined, and complete the local numerical approximation of the enhanced moving least squares meshless algorithm. Equation (12) represents the interpolation scheme of the standard moving least squares meshless algorithm. In the enhanced least squares algorithm, its enhanced interpolation scheme can be expressed as: (13); in, This represents the number of field points involved in the interpolation using the enhanced least squares algorithm. The second term in the equation represents the local numerical approximation of the enhanced least squares algorithm. and This represents the local reinforcement interpolation basis function and the coefficients to be determined corresponding to the local numerical approximation constructed at each node.

[0027] In general, Lagrange polynomial basis functions can still be used as local enhancement interpolation basis functions, that is: (14) To ensure sufficient stability of the obtained system matrix equations, dimensionless relative coordinate values ​​are often used to construct local enhancement interpolation basis functions, i.e.: (15) In the formula, , and These are dimensionless relative coordinate values. This represents the average spacing between field points.

[0028] S4. Based on the nodal interpolation basis functions of the standard moving least squares meshless algorithm obtained in step S2 and the enhanced local numerical approximation obtained in step S3, a global interpolation scheme for general scalar fields based on the enhanced least squares algorithm is constructed under the premise of satisfying the unit decomposition property. The final global interpolation scheme in the enhanced least squares algorithm can be expressed as: (16); S5. Based on the global interpolation scheme of the general scalar field obtained in step S4, based on the enhanced least squares algorithm, and based on the partial differential control equation of the wave propagation problem, the Galerkin weighted residual method and Newmark time integration method are used to complete the time discretization and spatial discretization of the transient wave propagation problem, and the system matrix equation of the dispersion analysis of the transient wave propagation problem is obtained. For the wave propagation problem in an ideal, stationary, homogeneous fluid medium, the governing equation for small-amplitude waves can be expressed as: (17); in, Represents the Laplace operator. Represents variables (such as sound pressure) in the problem domain of wave propagation. Indicates the speed of sound wave propagation. Indicates time.

[0029] Based on the principle of virtual displacement, the integral form ("weak form") of the partial differential governing equations for the transient sound propagation problem can be expressed as: (18); in, This represents any possible virtual displacement of the system (here referring to sound pressure). Indicates sound pressure Regarding time The second derivative of .

[0030] Based on integration by parts and Gaussian divergence theorem, the integral form of the partial differential governing equations for the transient sound propagation problem described above can be expressed as: (19); In the formula, This represents the outward normal unit vector along the boundary of the computational domain of sound propagation space.

[0031] By employing the Galerkin weighted residual method and substituting the sound pressure interpolation scheme for the transient acoustic scattering problem of underwater targets based on standard quadrilateral element grids into the above equation, we can obtain the matrix form of the above integral equation: (20); in, yes The sound pressure level of the unknown node corresponding to each node. It is the mass matrix of the transient acoustic scattering problem. It is by A vector composed of interpolation shape functions at each node. It is the interpolation form function corresponding to the node. It is the total number of nodes within the spatial computation domain. It is the smooth stiffness matrix of the transient acoustic scattering problem. It is the smoothed sound pressure gradient matrix of the transient sound scattering problem. It is the force matrix of the transient acoustic scattering problem. It is the normal velocity vector of the acoustic particles at the node along the boundary of the spatial computational domain.

[0032] When performing dispersion analysis on transient wave propagation problems, external loads and boundary conditions are generally not considered. Based on the global interpolation scheme in the previously established enhanced least squares algorithm, the matrix governing equations for the transient wave propagation problem can be written as follows: (twenty one); The solutions to the above matrix equations generally have the following form: (twenty two); in, The imaginary unit, This represents the wavenumber calculated using numerical methods. Represents the unit vector in the direction of wave propagation. This represents the position vector of a point in the field along the direction of wave propagation. This represents the angular frequency calculated using numerical methods. Indicates and The magnitude vector corresponding to the coefficients of each unknown node. This represents the time variable in numerical computation.

[0033] (twenty three); In the above formula, the magnitude vector corresponding to each node is exactly the same without considering the boundary conditions, that is, the vector... This will appear repeatedly in the above formula.

[0034] for Figure 4 Given the field point distribution, substituting equation (22) into equation (21), we get: (twenty four) in, and Representation and mass matrix and stiffness matrix The relevant coefficient matrix, Indicates the accurate wave number, It is a time Sum of angular frequencies The relevant magnitude vector.

[0035] (25); in, Is and Representing the stiffness matrix The Middle row and number The elements of the column.

[0036] (26); in, Is and Representing the mass matrix respectively The Middle row and number The elements of the column.

[0037] In the Newmark time integration method, the sound pressure vector at each node , and its first and second derivative vectors with respect to time and It can be represented as: (27); in, It is the time integration step size. and These are the parameters of the numerical damping in the Newmark time-integral algorithm under two conditions, when , At that time, the Newmark time integration algorithm does not have spurious numerical damping, and without numerical damping, numerical errors caused by numerical damping can be avoided.

[0038] From equation (24), we can obtain that , and The dynamic equilibrium equations at equal times are: (28); Based on equations (27) and (28), we can obtain: (29); in, For accurate angular frequency, parameters , and All are about parameters The function.

[0039] The formula Substituting into equation (29), we get: (30); From the above formula, we can obtain: (31); The total numerical dispersion error of the numerical solution to the transient wave propagation problem can be expressed as: (32); In the formula, , and For accurate sound propagation speed, accurate wavenumber, and accurate angular frequency; The velocity of sound is obtained through numerical calculation.

[0040] because: (33); in, This represents the average distance between nodes spatially discretized in the acoustic propagation computational domain. This indicates a dimensionless time step.

[0041] From equations (32) and (33), we can obtain: (34); As can be seen from equation (31), because the parameter and All are about parameters The function, therefore the parameters Also about parameters The function is such that equation (34) can be written as: (35); in, Represents a parameter The function.

[0042] For the transient partial differential equation of sound propagation as shown in equation (17), considering the time harmonic characteristics of the small amplitude wave propagation problem, the sound pressure distribution in the wave propagation problem domain can be expressed as: (36); in, This represents the spatial distribution of sound pressure amplitude.

[0043] Substituting equation (36) into equation (17), we obtain the following small-amplitude wave equation for the frequency case in an ideal fluid medium: (37); in, Represents wave number and .

[0044] Based on the global interpolation scheme in the previously established enhanced least squares algorithm, the matrix form of the above partial differential equation can be expressed as: (38); in, and Representing the acoustic stiffness matrix and mass matrix, This represents the unknown node vector to be determined.

[0045] If at each node in the enhanced least squares algorithm there is If there are unknown node coefficients, then the vector It can be represented as: (39); for Figure 4 Given the field point distribution, substituting equation (39) into equation (38), we get: (40); The condition for the existence of a non-zero solution in the above equation is that the determinant of its coefficient matrix is ​​zero, that is: (41); It can be seen from the above formula It is actually a matrix The eigenvalues, namely: (42); It can be seen from equations (25) and (26) that the matrix and All are numerical wavenumbers The function, therefore for any given exact wavenumber The corresponding numerical wavenumber can be calculated using the above formula. Furthermore, the spatial discretization error of the transient wave propagation problem can be calculated. The spatial discretization error of the transient wave propagation problem was obtained. Then, the total numerical dispersion error in the transient wave propagation problem can be calculated further according to equations (31) and (35). This allows us to determine the magnitude of dispersion error when using different numerical algorithms for wave propagation analysis, thus completing the dispersion error analysis.

[0046] S6. Based on the given node distribution, solve the system matrix equation of the transient wave propagation problem obtained in step S5, and calculate the total dispersion error of different calculation methods along different wave propagation directions under different time integration step sizes. S7. Based on the total dispersion error of the transient wave propagation problem in different wave propagation directions under different time integration step sizes obtained in step S6, plot the dispersion error as a function of the dimensionless wave number, compare the dispersion error results obtained by different calculation methods, analyze the coupling relationship between spatial dispersion error, time dispersion error and total dispersion error in the transient wave propagation problem, and complete the dispersion analysis of the wave propagation problem.

[0047] To compare and verify the ability of the enhanced moving least squares meshless-Newmark time integration algorithm presented in this invention to control the total dispersion error of transient wave propagation, the dispersion error along different wave propagation directions calculated by the standard moving least squares meshless algorithm is also calculated here.

[0048] Figure 5 The values ​​presented are the dispersion errors (including spatial discretization errors) in different wave propagation directions calculated using the standard moving least squares meshless-Newmark time integration algorithm. and total dispersion error As the dimensionless wavenumber changes, Figure 5 (a) indicates the direction of wave propagation. The curve of change, Figure 5 (b) is the direction of wave propagation. The curve of change, Figure 5 (c) indicates the direction of wave propagation. The curve of change, Figure 5 (d) indicates the direction of wave propagation. The curve of change; Figure 6 The given values ​​are the total dispersion error (including spatial discretization error) in different wave propagation directions calculated by the enhanced moving least squares meshless-Newmark time integration algorithm. and total dispersion error As the dimensionless wavenumber changes, Figure 6 (a) indicates the direction of wave propagation. The curve of change, Figure 6 (b) is the direction of wave propagation. The curve of change, Figure 6 (c) indicates the direction of wave propagation. The curve of change, Figure 6 (d) indicates the direction of wave propagation. The analysis results show that when using the standard moving least squares meshless-Newmark time integration algorithm, the dispersion error in each propagation direction is large and rapidly amplifies with the increase of the dimensionless wavenumber. Furthermore, there are significant differences in the total dispersion error in different directions, reflecting the obvious "numerical anisotropy" characteristic of this method, which affects the stability and consistency of the calculation results. In contrast, the enhanced moving least squares meshless-Newmark time integration algorithm proposed in this invention shows significant advantages under the same conditions. This method not only significantly reduces the total dispersion error but also effectively suppresses the numerical anisotropy phenomenon, making the calculation results along different wave propagation directions tend to be consistent, and significantly improving the symmetry and reliability of the numerical solution. Simultaneously, as the time integration step size decreases, the total dispersion error continues to decrease, further verifying the convergence and accuracy advantages of this algorithm.

[0049] Therefore, this invention employs an enhanced moving least squares meshless Newmark integral method to reduce the total dispersion error of transient wave propagation. While ensuring computational efficiency, it improves the solution accuracy and significantly optimizes the dispersion characteristics and anisotropy control capability of transient wave propagation problems. It provides a high-precision and high-stability numerical solution method for wave propagation simulation in complex engineering scenarios, and has broad application prospects and important engineering promotion value.

[0050] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and not to limit them. Although the present invention has been described in detail with reference to preferred embodiments, those skilled in the art should understand that modifications or equivalent substitutions can still be made to the technical solutions of the present invention, and these modifications or equivalent substitutions cannot cause the modified technical solutions to deviate from the spirit and scope of the technical solutions of the present invention.

Claims

1. An enhanced moving least squares meshless-Newmark integration method for reducing the total dispersion error of transient wave propagation, characterized by: Includes the following steps: S1. Discretize the wave propagation problem domain into a series of nodes, and define an influence domain for each node; S2, selection weight function and interpolation basis function, based on the moving least squares method to fit the local approximate interpolation function, the standard moving least squares interpolation basis function of each node is derived; S3. Select the Lagrange polynomial as the local reinforcement interpolation basis function, construct the node degrees of freedom at each node consisting of the local reinforcement interpolation basis function and the coefficients to be determined, and complete the local numerical approximation of the enhanced moving least squares meshless algorithm. S4. Based on the node interpolation basis functions of the standard moving least squares meshless algorithm obtained in step S2 and the enhanced local numerical approximation obtained in step S3, a global interpolation scheme for a general scalar field based on the enhanced least squares algorithm is constructed under the premise of satisfying the unit decomposition property. S5. Based on the global interpolation scheme of the general scalar field obtained in step S4, based on the enhanced least squares algorithm, and based on the partial differential control equation of the wave propagation problem, the Galerkin weighted residual method and Newmark time integration method are used to complete the time discretization and spatial discretization of the transient wave propagation problem, and the system matrix equation of the dispersion analysis of the transient wave propagation problem is obtained. S6. Based on the given node distribution, solve the system matrix equation of the transient wave propagation problem obtained in step S5, and calculate the total dispersion error of different calculation methods along different wave propagation directions under different time integration step sizes. S7. Based on the total dispersion error of the transient wave propagation problem in different wave propagation directions under different time integration step sizes obtained in step S6, plot the dispersion error as a function of the dimensionless wave number, compare the dispersion error results obtained by different calculation methods, analyze the coupling relationship between spatial dispersion error, time dispersion error and total dispersion error in the transient wave propagation problem, and complete the dispersion analysis of the wave propagation problem.

2. The enhanced moving least squares meshless-Newmark integral method for reducing the total dispersion error of transient wave propagation according to claim 1, characterized in that: In step S1, the influence domain is circular or rectangular, and the size of the influence domain is adapted to the average spacing of the node distribution.

3. The enhanced moving least squares meshless-Newmark integral method for reducing the total dispersion error of transient wave propagation according to claim 2, characterized in that: The support domain of the integration point is a circular or rectangular region centered on the integration point. The number of nodes participating in the interpolation is determined by the support domain of the integration point during numerical integration. The node influence domain is the area that each node can influence during numerical integration. The size of the node influence domain is closely related to the integration point support domain.

4. The enhanced moving least squares meshless-Newmark integral method for reducing the total dispersion error of transient wave propagation according to claim 3, characterized in that: Step S2 specifically includes the following steps: S21. The weight function of nodes within the influence domain is non-zero, and the weight function of nodes outside the influence domain is zero. The weight functions are shown below: ； in, The radius of the influence region; S22. For the transient wave propagation problem in a scalar field, based on the selected nodes within the influence domain and the standard moving least squares meshless technique, the initial form of the global interpolation scheme is constructed as follows: ； in, For the physical variables of transient wave propagation problems, This represents the number of nodes participating in the interpolation. The position vector at the interpolation point, vector express There are unknown coefficients to be determined, and the vector is... This represents a vector composed of the interpolation basis functions used. S23. For one-dimensional, two-dimensional, and three-dimensional problems, the Lagrange linear basis functions used are expressed as follows: ； S24. Constructing the weighted squared error: ； in, The total number of nodes contained in the supporting domain for the node at coordinate x. The function value at the node with coordinate x; To minimize the weighted squared error, the weighted squared error must be about The extreme points of the weighted squared difference function with respect to The reciprocal of is zero: ； Therefore, we can conclude that: ； in, This represents an unknown node variable corresponding to n nodes in a supporting domain, represented by a matrix. Specifically, it is expressed as follows: ； matrix Specifically, it is expressed as follows: ； Undetermined coefficient matrix for: ； matrix , and Substituting the initial form of the global interpolation format, we get: ； in, For each node constructed, use the standard moving least squares interpolation basis function.

5. The enhanced moving least squares meshless-Newmark integral method for reducing the total dispersion error of transient wave propagation according to claim 4, characterized in that: Step S3 specifically includes the following steps: S31. Enhanced interpolation scheme built on enhanced least squares algorithm: ； in, To increase the number of nodes involved in the least squares interpolation algorithm, To enhance the local numerical approximation of the least squares algorithm, and For each node, construct a local numerical approximation corresponding to the local enhancement interpolation basis function and the coefficients to be determined; S32. A local enhancement interpolation basis function is constructed using Lagrange polynomial basis functions and dimensionless relative coordinate values. The formula for the local enhancement interpolation basis function is: ； in, , and These are dimensionless relative coordinate values. This represents the average spacing between nodes.

6. The enhanced moving least squares meshless-Newmark integral method for reducing the total dispersion error of transient wave propagation according to claim 5, characterized in that: The global interpolation format obtained in step S4 of the enhanced least squares algorithm is as follows: 。 7. The enhanced moving least squares meshless-Newmark integral method for reducing the total dispersion error of transient wave propagation according to claim 6, characterized in that: Step S5 specifically includes the following steps: S51. Based on the principle of virtual displacement, the partial differential control equations of the transient sound propagation problem are transformed into integral form, as shown below: ； in, For the Laplace operator, For the speed of sound propagation, For sound pressure, Indicates sound pressure Regarding time The second derivative; S52. Based on integration by parts and the Gaussian divergence theorem, the integral form of the partial differential control equations for the transient sound propagation problem is expressed as: ； in, The unit vector of the outward normal to the boundary of the computational domain along the sound propagation space; S53. Using the Galerkin weighted residual method, the sound pressure interpolation scheme for the transient acoustic scattering problem of underwater targets based on standard quadrilateral element grids is substituted into the formula in S52, resulting in the matrix form of the integral equation: ； in, for The sound pressure level of the unknown node corresponding to each node. The mass matrix for the transient acoustic scattering problem is... For the reason A vector composed of interpolation shape functions at each node. For the interpolation shape function corresponding to the node, This represents the total number of nodes within the spatial computational domain. Let be the smooth stiffness matrix for the transient acoustic scattering problem. For the transient sound scattering problem, the smoothed sound pressure gradient matrix is... The force matrix for a transient acoustic scattering problem. Let be the normal velocity vector of the acoustic particles at the nodes along the boundary of the spatial computational domain; S54. Based on the global interpolation scheme in the enhanced least squares algorithm obtained in step S4, the matrix governing equation for the transient wave propagation problem is written as follows: ； The solution to the matrix governing equation is: ； in, The imaginary unit, The wavenumber is calculated using numerical methods. Let be the unit vector in the direction of wave propagation. Let be the position vector of the field point along the wave propagation direction. The angular frequency is calculated using numerical methods. To and The magnitude vector corresponding to the coefficients of each unknown node. For time variation in numerical computation; S55. In the Newmark time integration method, the sound pressure vector at each node... , and its first and second derivative vectors with respect to time and Represented as: ； in, It is the time integration step size. and These are the parameters of the numerical damping in the Newmark time-integration algorithm under two conditions.

8. The enhanced moving least squares meshless-Newmark integral method for reducing the total dispersion error of transient wave propagation according to claim 7, characterized in that: In step S6, the total numerical dispersion error of the numerical solution to the transient wave propagation problem is: ； in, For accurate sound propagation speed, For accurate wavenumbers, For accurate angular frequency, The velocity of sound is obtained through numerical calculation.

Images