An h-adaptive moment method applied to electromagnetic scattering

Abstract

The application discloses an h-adaptive moment method applied to electromagnetic scattering, and comprises the following steps: an electromagnetic model for calculating scattering characteristics of a to-be-measured target is created, and initial mesh division under plane wave excitation is carried out; current coefficients of each grid element at present are solved by using a moment method based on RWG base functions, and then surface charge density at a super-convergent point in the corresponding grid element is calculated; local errors of the corresponding grid element are calculated based on the surface charge density at the super-convergent point in each grid element and a super-convergent point recovery method; global errors are calculated based on the local errors of each grid element; whether an iteration stopping condition is met is judged; if not, mesh refinement is carried out on some grid elements based on a B-W interpolation point algorithm, re-divided grids are obtained, and the current coefficient solving step is returned; if yes, iteration is stopped, and the current solved current coefficient is output. The application can guarantee calculation precision, robustness, reduce error estimation, reduce the number of mesh refinement and reduce the consumption of calculation resources.

Description

Technical Field

[0001] This invention belongs to the field of electromagnetic fields, and specifically relates to an h-adaptive moment method applied to electromagnetic scattering. Background Technology

[0002] The Method of Moments (MoM), a classic numerical analysis method, has been widely applied in electromagnetic calculations due to its high computational accuracy and broad applicability. For example, stealth and anti-stealth technologies require precise analysis of the electromagnetic scattering characteristics of targets. However, accurate mesh generation is a necessary prerequisite for high-precision calculations in the electromagnetic field. This is because if the mesh size is too large, it cannot accurately describe the fine structure of the target or the complex changes in the electromagnetic flow near the source, leading to numerical errors. Conversely, if the mesh size is too small, the number of meshes will increase non-linearly, resulting in a sharp increase in computational resources. Furthermore, determining the correctness of the mesh and re-meshing can both lead to significant time waste.

[0003] Currently, a scheme to reduce numerical solution error by controlling local mesh size, namely h-adaptive refinement, can effectively solve the above problems. This scheme includes an error estimator and a controller; the error estimator is used to calculate the mesh cell error, while the controller determines which meshes to refine and sets the iteration termination condition.

[0004] In 2018, SK Kim proposed an h-adaptive moment method based on current density discontinuity error estimation. This method uses the current density discontinuity between grid cells to calculate the error value of the grid cells. This error estimation method consumes very little time, but has low robustness.

[0005] In contrast, other commonly used error estimation methods are typically residual error estimation, which estimates the error by changing the test function or boundary conditions in the selection process. While this method is highly robust, it requires computations comparable to matrix filling and extensive use of numerical integration, resulting in a long error estimation time. Other error estimation methods with shorter computation times, such as discontinuity error estimation, often fail. Therefore, reducing error estimation time while maintaining computational accuracy and robustness is a current technical challenge. Furthermore, minimizing mesh refinement to reduce computational resource consumption is also a pressing technical problem in this field. Summary of the Invention

[0006] To address the aforementioned problems in the prior art, this invention provides an h-adaptive moment method applied to electromagnetic scattering. The technical problem to be solved by this invention is achieved through the following technical solution:

[0007] An electromagnetic model is created to calculate the scattering characteristics of the target under test. An initial mesh is generated for the electromagnetic model under plane wave excitation to obtain a triangular mesh on the surface of the target under test.

[0008] The method of moments based on RWG basis functions is used to solve the current coefficient of each grid cell in the current division, and the surface charge density at the superconvergence point in the corresponding grid cell is calculated using the solved current coefficient.

[0009] Based on the surface charge density at the superconvergence point in each grid cell and the superconvergence point recovery method, the local error of the corresponding grid cell is calculated.

[0010] The global error is calculated based on the local error of each grid cell;

[0011] Based on the global error, determine whether the iteration stopping condition is met;

[0012] If not, select some mesh elements to refine the mesh based on the BW interpolation algorithm to obtain a re-divided mesh; and return to the step of using the method of moments based on RWG basis functions to solve the current coefficient of each mesh element currently divided, and using the solved current coefficient to calculate the surface charge density at the superconvergence point in the corresponding mesh element.

[0013] If so, stop the iteration and output the current coefficients of each grid cell that have been solved.

[0014] In one embodiment of the present invention, the initial mesh generation of the electromagnetic model under plane wave excitation includes:

[0015] The electromagnetic model is initially meshed with a size of λ4; where λ represents the wavelength of the plane wave.

[0016] In one embodiment of the present invention, the step of using the method of moments based on RWG basis functions to solve for the current coefficients of each currently divided grid cell, and using the solved current coefficients to calculate the surface charge density at the superconvergence point in the corresponding grid cell, includes:

[0017] The equivalent electromagnetic current of the target surface under test is discretized using the method of moments based on RWG basis functions and the currently divided grid. The basis functions generated after grid discretization are used to form a system of linear equations. The current coefficient of each grid cell is obtained by solving the system of linear equations.

[0018] Using the formulas for calculating the current coefficient and surface charge density of each grid cell, the surface charge density at the superconvergence point in the corresponding grid cell is calculated.

[0019] In one embodiment of the present invention, the calculation of the local error of the corresponding grid cell based on the surface charge density at the superconvergence point in each grid cell and the superconvergence point recovery method includes:

[0020] For each vertex belonging to each mesh cell, the recovered charge density value of the vertex is obtained by using the surface charge density at the superconvergence point in all mesh cells containing that vertex.

[0021] For each grid cell, the recovered charge density values ​​obtained from the three vertices of the grid cell are linearly fitted using the least squares method to construct the linear recovered charge density function of the grid cell.

[0022] To ensure the properties of the RWG basis function with respect to the charge on the positive and negative triangles in the triangle pair, the modified recovery charge density function of the grid cell is obtained by adding the linear recovery charge density function of the grid cell to a constant.

[0023] The local error of the grid cell is calculated based on the corrected recovery charge density function of the grid cell and the surface charge density at the original superconvergence point.

[0024] In one embodiment of the present invention, obtaining the recovered charge density value of the vertex using the surface charge density at the superconverged points in all mesh cells containing the vertex includes:

[0025] The surface charge density at the superconvergence point in all mesh cells containing the vertex is averaged, and the resulting average value is determined as the recovered charge density value of the vertex.

[0026] In one embodiment of the present invention, the formula used to calculate the local error of the mesh cell includes:

[0027]

[0028] in, ρ represents the local error of the i-th grid cell; i ' ,rec ρ represents the corrected recovery charge density function of the i-th grid cell; i S represents the surface charge density at the superconvergence point in the i-th grid cell; i This represents the i-th grid cell.

[0029] In one embodiment of the present invention, the calculation of the global error based on the local error of each grid cell includes:

[0030] For each grid cell, the local error of that grid cell is weighted by using the area of ​​its triangle as its own weight, and the weighted local error of that grid cell is obtained; and the sum of the weighted local errors of all grid cells is obtained, and the sum of the weighted local errors is obtained.

[0031] The weighted local error is normalized using the sum of the triangle areas of all grid cells;

[0032] The square root of the normalized result is used to obtain the global error.

[0033] In one embodiment of the present invention, the step of determining whether the iteration stopping condition is met based on the global error includes:

[0034] When the global error decreases to below 1 / Q of the preset initial global error, or the number of iterations reaches the preset maximum number of iterations, the iteration stopping condition is determined to be met; otherwise, the iteration stopping condition is determined not to be met; where Q is a natural number greater than 1.

[0035] In one embodiment of the present invention, the selection of partial mesh elements for mesh refinement based on the BW interpolation algorithm includes:

[0036] The weighted local error threshold is obtained by multiplying the maximum weighted local error in all grid cells by a preset percentage.

[0037] Mesh cells with a weighted local error greater than the weighted local error threshold are selected as the selected mesh cells;

[0038] The selected grid cells are refined using a BW interpolation algorithm.

[0039] In one embodiment of the present invention, mesh refinement is performed on any selected mesh cell using a BW interpolation algorithm, including:

[0040] In the current Delaunay triangulation M i The search found the point p to be inserted. i The base unit; where the current Delaunay triangulation M i This represents the currently defined triangular mesh;

[0041] Using the adjacency relationship of mesh elements, mark all mesh elements that do not conform to the circumcircle criterion at the adjacent mesh elements of the base element, and denote the set of such mesh elements as B. i ;

[0042] In the current Delaunay triangulation M i Delete B i The covered grid cells form cavity C. i , connect p i and C i The vertices on the boundary form a new triangulation B within the cavity. i+1 This yields the new Delaunay triangulation M. i+1 .

[0043] The beneficial effects of this invention are:

[0044] The h-adaptive moment method for electromagnetic scattering provided in this invention first creates an electromagnetic model for calculating the scattering characteristics of the target under test. An initial mesh is then created for the electromagnetic model under plane wave excitation, resulting in a triangular mesh on the target surface. Next, the moment method based on RWG basis functions is used to solve for the current coefficients of each mesh element. The surface charge density at the superconvergence point in the corresponding mesh element is then calculated using the solved current coefficients. Following this, the local error of the corresponding mesh element is calculated based on the surface charge density at the superconvergence point and the superconvergence point recovery method. Then, the global error is calculated based on the local errors of each mesh element. Finally, the iteration stopping condition is determined based on the global error. If not, a portion of the mesh elements is selected for mesh refinement using the BW interpolation algorithm, resulting in a re-divided mesh. The process then returns to the step of solving for the current coefficients using the moment method based on RWG basis functions. If yes, the iteration stops, and the current coefficients of each mesh element are output.

[0045] This invention aims to obtain accurate numerical results while controlling the number of mesh elements. Starting with an accurate and efficient error estimator and controller, it proposes an SPR error estimator that combines the Superconvergence Point Recovery (SPR) method with the characteristics of RWG basis functions. This estimator can efficiently and accurately calculate the error of mesh elements, exhibiting strong robustness and short computation time. Based on this, the h-adaptive moment method implemented using this error estimator achieves high computational accuracy, effectively controls computational resources, can adaptively refine the mesh without manual intervention, and automatically stops iteration after setting the stopping conditions. This effectively controls the number of mesh elements and yields accurate simulation results. Attached Figure Description

[0046] Figure 1 This is a schematic flowchart of an h-adaptive moment method for electromagnetic scattering provided in an embodiment of the present invention;

[0047] Figure 2 This is a schematic diagram of the triangular mesh and the superconvergence point of the mesh element in an embodiment of the present invention;

[0048] Figure 3 (a)~ Figure 3 (d) is a schematic diagram of the process of refining selected grid cells using the BW interpolation algorithm in an embodiment of the present invention;

[0049] Figure 4 This is a schematic diagram of the electromagnetic model of the cone and the plane wave excitation in the simulation experiment of this invention embodiment;

[0050] Figure 5 This is a graph showing the variation of the number of cone meshes and global error with adaptive iteration in the simulation experiment of this embodiment of the invention;

[0051] Figure 6 This is a diagram showing the initial mesh distribution of the cone in the simulation experiment of this invention.

[0052] Figure 7 This is a mesh distribution diagram of the cone after adaptive refinement in the simulation experiment of this embodiment of the invention;

[0053] Figure 8 This is a comparison of the bistatic RCS of the cone with adaptive refinement and the initial mesh in the simulation experiment of this embodiment of the invention;

[0054] Figure 9 This is a diagram showing the conventional mesh distribution of the cone in the simulation experiment of this invention.

[0055] Figure 10 This is a comparison of the RCS of the cone in the simulation experiment of the embodiment of the present invention, using adaptive refinement and conventional subdivision at two stations. Detailed Implementation

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

[0057] To facilitate understanding of the solutions of the embodiments of the present invention, the relevant technologies and inventive concepts of the embodiments of the present invention will be briefly introduced first.

[0058] Obtaining the radar cross section (RCS) of a target is of significant engineering importance for applications such as stealth technology and anti-stealth detection technology. However, obtaining data through experiments requires consideration of issues related to site, funding, and technology, making implementation difficult. Electromagnetic simulation, as an efficient, convenient, and practical method, is increasingly favored in obtaining the RCS of a target. This invention first studies existing algorithms such as the Method of Moments (MMT). Among them, the MMT based on RWG basis functions discretizes the equivalent electromagnetic flow of the target surface using a triangular mesh. Inaccurate meshing can introduce errors during the solution process. For example, some curved surfaces may be replaced by planar meshes, causing model deformation. In some fine structures, due to drastic changes in electromagnetic flow, the requirements for the mesh are also very high. To quantify this error and control the mesh size based on it, and to obtain accurate numerical results while controlling the number of meshes, this invention proposes an h-adaptive MMT based on the Superconvergent Patch Recovery (SPR) error estimation method for calculating the electromagnetic scattering characteristics of the target. This is based on how to implement an accurate and efficient error estimator and controller. Specifically, in this embodiment of the invention, the superconvergence point recovery method is introduced into the error estimation of the method of moments, and the error estimation method is optimized by combining the characteristics of the RWG basis function, and the mesh is refined by combining the fast BW interpolation algorithm.

[0059] like Figure 1 As shown in the embodiment of the present invention, an h-adaptive moment method for electromagnetic scattering may include the following steps:

[0060] S1, Create an electromagnetic model for calculating the scattering characteristics of the target under test, perform initial meshing of the electromagnetic model under plane wave excitation, and obtain a triangular mesh divided on the surface of the target under test;

[0061] In the embodiments of the present invention, the target to be tested can be an electromagnetic scattering object with an arbitrary geometric structure made of an ideal electrical conductor (PEC), such as an airplane, ship, missile, etc.

[0062] To calculate the scattering characteristics of the target under test, a geometric model corresponding to the design requirements needs to be constructed based on the target's structure, size, and boundary continuity. The outer surface of this geometric model is then set to a Perfect Electric Conductor (PEC) material to create an electromagnetic model of the target. Next, a plane wave excitation is applied to the electromagnetic model, followed by initial mesh generation.

[0063] In one optional implementation, the initial mesh generation of the electromagnetic model under plane wave excitation includes:

[0064] The electromagnetic model is initially meshed with a size of λ / 4.

[0065] Wherein, λ represents the wavelength of the plane wave.

[0066] After initial mesh generation, a triangular mesh is generated on the surface of the target object, with each triangle corresponding to a mesh cell.

[0067] S2, using the method of moments based on RWG basis functions, solve for the current coefficient of each grid cell in the current division, and use the solved current coefficient to calculate the surface charge density at the superconvergence point in the corresponding grid cell;

[0068] In one optional implementation, S2 may include:

[0069] S21, the equivalent electromagnetic current of the target surface under test is discretized using the method of moments based on RWG basis functions and the currently divided grid. The basis functions generated after grid discretization are used to form a linear equation system. The current coefficient of each grid cell is obtained by solving the linear equation system.

[0070] For details on this step, please refer to the relevant technical explanations:

[0071] [I n ] = [z mn ] -1 ·[g m ]

[0072] Among them, I n Z is the current coefficient. mn Let g be the impedance matrix. m This is the activation matrix.

[0073] Finally, the current density on the metal surface, characterized by the RWG basis functions, is:

[0074]

[0075] Among them, I n N is the current coefficient. s The number of basis functions. These are basis functions.

[0076] It is understandable that the surface current density can be obtained by solving for the current coefficient. Based on existing formulas, the radar cross section of the target, the gain of the antenna, the directivity, and the S-parameters of microwave devices can be calculated from the surface current density.

[0077] S22, using the current coefficient and surface charge density calculation formula of each grid cell, calculate the surface charge density at the superconvergence point in the corresponding grid cell.

[0078] The superconvergence point recovery method points out that the centroid of a triangle has superconvergence because the centroid is farthest from the edge with discontinuities. Therefore, the electromagnetic variables at the superconvergence point, i.e. the centroid, have the highest calculation accuracy. In other words, the electromagnetic variables at the centroid of a triangle can be used for error analysis. Therefore, this embodiment of the invention introduces the superconvergence point recovery method into the error estimation of the method of moments. For each grid cell, the surface charge density at the centroid of the grid cell is calculated using the surface charge density calculation formula and the current coefficient of the grid cell.

[0079] The formula for calculating surface charge density is as follows:

[0080]

[0081] in, Let ω be the surface current density and ω be the angular frequency.

[0082] S3, based on the surface charge density at the superconvergence point in each grid cell and the superconvergence point recovery method, calculate the local error of the corresponding grid cell;

[0083] This step is implemented by the SPR error estimator of this embodiment of the invention. In one optional implementation, S3 may include the following steps:

[0084] S31, For each vertex belonging to each mesh cell, the recovered charge density value of the vertex is obtained by using the surface charge density at the superconvergence point in all mesh cells containing the vertex.

[0085] Specifically, the recovered charge density value of the vertex is obtained by utilizing the surface charge density at the superconvergence points in all mesh cells containing the vertex, including:

[0086] The surface charge density at the superconvergence point in all mesh cells containing the vertex is averaged, and the resulting average value is determined as the recovered charge density value of the vertex.

[0087] Please see details. Figure 2 The diagram shows a triangular mesh and the superconvergence point of each mesh element; five triangles are given as an example, with the centroid of each triangle denoted by ρ1, ρ2, ρ3, ρ4, and ρ5, and the common vertex of the five triangles denoted by K. This step analyzes a vertex on a mesh element, such as... Figure 2 As shown, the average value obtained by averaging the surface charge density at the superconvergence points of all mesh cells containing the vertex is the recovered charge density value of that vertex.

[0088]

[0089] Where, ρ k This represents the recovered charge density value at that vertex; ρ i The surface charge density at the superconvergence point of all mesh cells containing the vertex is represented by h; h represents the number of all mesh cells containing the vertex.

[0090] Understandably, through the SPR method, each vertex of the mesh cell corresponding to the triangle obtains a restored charge density value.

[0091] S32, For each grid cell, the recovered charge density values ​​obtained from the three vertices of the grid cell are linearly fitted using the least squares method to construct the linear recovered charge density function of the grid cell;

[0092] The least squares fitting process in this step can be represented as follows:

[0093] 1) Let ρ rec =Pa;

[0094] Where, ρ rec This represents the linear recovery charge density function of the grid cell; P represents the coordinate vector; a represents the coefficient vector;

[0095] P=[1,x,y,z]; a=[a1,a2,a3,a4] T

[0096] 2) Fitting formula using the least squares method

[0097]

[0098] Here, (xi,yi,zi) represents the three coordinate values ​​of the corresponding point.

[0099] It can be seen that when F(a) reaches its minimum value:

[0100]

[0101] Solving for:

[0102] a = A -1 b

[0103] in:

[0104]

[0105]

[0106] Then substitute a into ρ rec The linear recovery charge density function ρ can be obtained from pa. rec .

[0107] S33, to ensure the characteristics of the RWG basis function with respect to the charge on the positive and negative triangles in the triangle pair, the modified recovery charge density function of the grid cell is obtained by adding the linear recovery charge density function of the grid cell to a constant;

[0108] According to the RWG basis function property, for any pair of triangles, the charge on the positive and negative triangles is constant, and its sum is zero. To ensure this property, this embodiment of the invention requires adding a constant to the constructed linear recoverable charge density function for each mesh cell, such that the recovered charge density at the center point of the triangle is the same as the original surface charge density at that point. Based on this objective, a constant can be determined, and the sum of the linear recoverable charge density function of the mesh cell and this constant is determined as the modified recoverable charge density function of the mesh cell.

[0109] S34. Based on the corrected recovery charge density function of the grid cell and the surface charge density at the original superconvergence point, calculate the local error of the grid cell.

[0110] Specifically, the formulas used to calculate the local error of this mesh element include:

[0111]

[0112] in, ρ represents the local error of the i-th grid cell; i ' ,rec ρ represents the corrected recovery charge density function of the i-th grid cell; i S represents the surface charge density at the superconvergence point in the i-th grid cell; i This represents the i-th grid cell.

[0113] For each grid cell, a local error can be calculated by following the steps described above.

[0114] This invention improves the error estimator's processing procedure by utilizing the characteristic of the RWG basis function method of moments (SPR) to discretize the target using a triangular mesh. The superconvergent point recovery (SPR) method is then introduced into the error calculation of the SPR method. Since the recovery method obtains a higher-order, more accurate numerical result based on the surface charge density of adjacent meshes, using the recovered solution as a reference solution for error estimation is more accurate than directly comparing the discontinuities in the surface charge density of mesh cells, thus exhibiting strong robustness. Furthermore, calculating the error of each mesh cell only requires a simple scalar multiplication of the surface charge density of adjacent meshes and involves only one numerical integration, resulting in short computation time. Therefore, this error estimator combines the SPR method with the characteristics of the RWG basis function, enabling efficient and accurate calculation of mesh cell errors with strong robustness and short computation time. It balances accuracy and efficiency, reduces the amount of mesh refinement, and saves time spent on error estimation.

[0115] S4, calculate the global error based on the local error of each grid cell;

[0116] This step is implemented by the controller. In one optional implementation, the calculation of the global error based on the local errors of each grid cell includes:

[0117] S41, for each grid cell, the local error of the grid cell is weighted by using the area of ​​its triangle as its own weight, and the weighted local error of the grid cell is obtained; and the sum of the weighted local errors of all grid cells is obtained.

[0118] S42, normalize the weighted local error using the sum of the triangle areas of all grid cells;

[0119] S43. Calculate the square root of the normalized result to obtain the global error.

[0120] The above steps are expressed by the following formula:

[0121]

[0122] Among them, e global Indicates global error; A represents the local error of the i-th grid cell; i The area of ​​the triangle in the i-th grid cell is represented by N; N represents the number of grid cells.

[0123] S5, determine whether the iteration stopping condition is met based on the global error;

[0124] This step may include:

[0125] When the global error is reduced to below 1 / Q of the preset initial global error, or the number of iterations reaches the preset maximum number of iterations, the iteration stopping condition is determined to be met; otherwise, the iteration stopping condition is determined not to be met; that is, when neither of the above two conditions is met, the iteration stopping condition is determined not to be met.

[0126] Where Q is a natural number greater than 1. The specific value can be set as needed; for example, in one possible implementation, Q is 4.

[0127] Similarly, the maximum number of preset iterations can be set as needed, such as 10 times, etc.

[0128] If not, proceed to step S6, select a portion of the mesh elements for mesh refinement based on the BW interpolation algorithm to obtain a re-divided mesh; and return to the step of using the method of moments based on RWG basis functions to solve the current coefficient of each currently divided mesh element, and using the solved current coefficient to calculate the surface charge density at the superconvergence point in the corresponding mesh element; that is, return to step S2.

[0129] Specifically, the mesh refinement based on the BW interpolation algorithm for the selected mesh elements in step S6, resulting in a re-divided mesh, refers to refining the mesh in large error regions based on the local errors of each mesh element, which is implemented by the controller.

[0130] Specifically, S6 may include:

[0131] S61, calculate the weighted local error threshold by multiplying the maximum weighted local error in all grid cells by a preset percentage;

[0132] The weighted local error threshold obtained in this step can be expressed as:

[0133]

[0134] Where u% represents the preset percentage, the specific value can be set reasonably as needed, such as 70%, and u is a natural number greater than 50; This represents the largest weighted local error among all mesh cells. The weighted local error used here is to prevent excessive differences in mesh area during the refinement process, which would lead to a deterioration in the overall mesh quality. At the same time, it ensures that mesh cells with large errors that are not included in this refinement meet the criteria for large error regions in the next iteration.

[0135] S62, determine the grid cells whose weighted local error is greater than the weighted local error threshold as selected grid cells;

[0136] In this embodiment of the invention, the large error region, i.e., the selected grid cell, is defined as a grid cell whose weighted local error is greater than the weighted local error threshold. Taking a preset percentage as an example, it is expressed as follows:

[0137]

[0138] The above process can identify at least one selected grid cell.

[0139] S63 uses a BW interpolation algorithm to refine the mesh for each selected mesh cell.

[0140] The BW interpolation algorithm is based on the circumcircle criterion of Delaunay triangulation and belongs to the incremental interpolation method. Assume T... n It is a point set P n ={p1,p2,…,p n Delaunay triangulation of} for T n For any simplex t, find the radius and center of its circumsphere, respectively. t and q t In T n Insert a new point p in the middle n+1 Define all circumscribed balls as containing p n+1 The set of simple shapes is B.

[0141] B = {t|t∈T} n ,d(p n+1 ,q t ) <r t}

[0142] If p is guaranteed n+1 In P n If B is inside the convex hull, then B is not empty. Delete all simplexes inside B, and then... n An internal cavity C is formed, which theoretically guarantees that cavity C is star-shaped, and that all points on it are perpendicular to p. n+1 Visual. Connect p n+1 By combining the points on cavity C, we can obtain the point set P. n+1 ={p1,p2,…,p n ,p n+1 Delaunay triangulation T n+1 .

[0143] For a summary of the above theories, please refer to [link / reference]. Figure 3 (a)~ Figure 3 (d) Understand the processing procedure of the single-step BW incremental interpolation method. Specifically, the BW-based interpolation algorithm is used to refine any selected mesh cell, including the following steps:

[0144] ① In the current Delaunay triangulation Mi The search found the point p to be inserted. i The base unit;

[0145] Among them, the current Delaunay triangulation M i This represents the currently generated triangular mesh; please refer to [link / reference]. Figure 3 (a) where all triangles form the current Delaunay triangulation M i Insertion point p i for Figure 3 (b) The points represented by the small green squares. The base element is... Figure 3 (b) contains p i That triangle with the green side.

[0146] ② Using the adjacency relationship of mesh elements, mark all mesh elements that do not conform to the circumcircle criterion at the adjacent mesh elements of the base element, and denote the set of such mesh elements as B. i ;

[0147] Among them, B i for Figure 3 (b) shows the three triangles inside the black border.

[0148] ③ In the current Delaunay triangulation M i Delete B i The covered grid cells form cavity C. i , connect p i and C i The vertices on the boundary form a new triangulation B within the cavity. i+1 This yields the new Delaunay triangulation M. i+1 .

[0149] Among them, C i for Figure 3 (c) The cavity formed by the black border; B i+1 for Figure 3 The five triangles within the black border in (d); M i+1 Depend on Figure 3 (d) is composed of all the triangles.

[0150] The entire process can be vividly described by the following formula:

[0151] M i+1 =M i -B i +B i+1

[0152] For details on the specific process, please refer to the relevant technical explanations; they will not be elaborated upon here.

[0153] The embodiments of the present invention use the BW interpolation algorithm for mesh refinement, which can meet the requirements of numerical simulation accuracy by continuously changing the mesh cells during the adaptive process. It can achieve adaptive mesh refinement without human intervention and obtain high-precision simulation results.

[0154] If so, execute S7 to stop the iteration and output the current coefficients of each grid cell that have been solved.

[0155] Understandably, the current coefficients of each grid cell can be used to further realize related applications, such as obtaining the radar cross section of a target, the gain and directivity of an antenna, and the S-parameters of microwave devices, which will not be described in detail here.

[0156] The solution provided in this invention aims to obtain accurate numerical results while controlling the number of meshes. Starting with an accurate and efficient error estimator and controller, an SPR error estimator is proposed. This error estimator combines the superconvergence point recovery (SPR) method with the characteristics of RWG basis functions, enabling efficient and accurate calculation of mesh element errors. It is robust and has a short computation time. Based on this, the h-adaptive moment method implemented using this error estimator achieves high computational accuracy, effectively controls computational resources, allows for adaptive mesh refinement without manual intervention, and automatically stops iterations after setting the stopping conditions. This effectively controls the number of meshes and obtains accurate simulation results.

[0157] To verify the effectiveness of the method in the embodiments of the present invention, the following simulation experiment is used for illustration.

[0158] The simulation uses a cone with a radius of 100mm and a height of 500mm, made of PEC. A vertically polarized plane wave with a frequency of 15GHz is incident on the tip of the cone. The corresponding electromagnetic model and plane wave excitation are as follows: Figure 4 As shown.

[0159] The large error region, i.e., the selected grid cell, is defined as follows:

[0160]

[0161] An adaptive refinement of the cone was performed using a controller based on the SPR error estimator results. The initial mesh for adaptive refinement was set to a conventional mesh with a maximum size of 1 / 4λ, and the initial mesh size was 6480. After the sixth iteration, the global error decreased to less than 1 / 4 of the initial global error, satisfying the iteration stopping condition, and the adaptive refinement ended. (See also...) Figure 5 The graph shows the variation of the number of meshes and the global error of the cone as a result of adaptive iteration. Table 1 shows the number of meshes and the global error for each iteration.

[0162] Table 1. Mesh count and global error after each iteration of adaptive refinement for the cone.

[0163]

[0164] On the one hand, in order to reflect the accuracy of the error estimator, we analyze it from the perspective of changes in the global error. Figure 5 The changes in the number of grids and global error with adaptive iteration are shown in Table 1, which shows the number of grids and global error after each iteration of adaptive refinement. Figure 6 and Figure 7 This is a map showing the initial and adaptively refined mesh distribution, with close-up details of some areas. Figure 6 The number of grids is 6480; Figure 7 The number of grid cells is 15120. (By...) Figure 5 It can be seen that, under the labeling of the SPR error estimator, the global error decreased significantly with each adaptive iteration while increasing the number of grid cells by a small amount. This reflects the accuracy of the SPR error estimator in this example.

[0165] On the other hand, the simulation results of this model under different mesh settings are analyzed. The bistatic RCS (radar cross section) calculated using this model with 88,544 meshes in the commercial software FEKO is used as a reference result. Figure 8 It can be seen that the bistatic RCS of the model after adaptive refinement is closer to the reference result than the initial mesh, especially showing a significant improvement in the θ = -180° to 0° range. -180°≤θ scat ≤180°, f=15GHz. θ represents the angle between the observation direction and the +Z axis. scat This indicates the angle between the observation direction and the +X axis; here... -180°≤θ scat ≤180° represents the xoz plane.

[0166] Finally, to demonstrate the advantages of adaptive refinement compared to conventional subdivision, such as... Figure 9 The cone is subjected to a conventional mesh of 0.1λ, with a mesh size of 32836. Figure 10 This is a comparison between adaptive refinement and conventional bi-station RCS, by Figure 10 The results show that the three curves coincide, indicating that adaptive refinement, with a significantly smaller mesh size than conventional meshing, produces simulation results that are essentially the same as conventional meshing. Table 2 shows that adaptive refinement reduces the mesh size by approximately 55% and saves about 80% of memory consumption compared to conventional meshing. This demonstrates a significant reduction in computational resources while maintaining high computational accuracy.

[0167] Table 2 Comparison of computational resources for adaptive refinement and conventional meshing of cones.

[0168] Subdivision type Adaptive refinement Conventional partitioning Reduction percentage (%) Number of grids 15120 32836 54.0 Memory usage (MB) 7848 37017 78.8

[0169] In summary, the SPR error estimator proposed in this invention has a significant advantage in computation time compared to commonly used residual error estimators, although it is slightly longer than discontinuous error estimators. In terms of accuracy, it performs comparably to residual error estimators and is significantly superior to discontinuous error estimators. This invention can achieve the shortest possible error estimation time and the fewest possible mesh refinements while ensuring computational accuracy, thus greatly reducing simulation time and computational resource consumption.

[0170] The above description is merely a preferred embodiment of the present invention and is not intended to limit the scope of protection of the present invention. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention are included within the scope of protection of the present invention.

Claims

1. An h-adaptive moment method applied to electromagnetic scattering, characterized in that, include: An electromagnetic model is created to calculate the scattering characteristics of the target under test. An initial mesh is generated for the electromagnetic model under plane wave excitation to obtain a triangular mesh on the surface of the target under test. The method of moments based on RWG basis functions is used to solve the current coefficient of each grid cell in the current division, and the surface charge density at the superconvergence point in the corresponding grid cell is calculated using the solved current coefficient. Based on the surface charge density at the superconvergence point in each grid cell and the superconvergence point recovery method, the local error of the corresponding grid cell is calculated. The global error is calculated based on the local error of each grid cell; Based on the global error, determine whether the iteration stopping condition is met; If not, select some mesh elements to refine the mesh based on the BW interpolation algorithm to obtain a re-divided mesh; and return to the step of using the method of moments based on RWG basis functions to solve the current coefficient of each mesh element currently divided, and using the solved current coefficient to calculate the surface charge density at the superconvergence point in the corresponding mesh element. If so, stop the iteration and output the current coefficients of each mesh element that have been solved so far; Among them, based on the surface charge density at the superconvergence point in each grid cell and the superconvergence point recovery method, the local error of the corresponding grid cell is calculated, including: For each vertex belonging to each mesh cell, the surface charge density at the superconvergence point in all mesh cells containing that vertex is averaged, and the average value is determined as the recovered charge density value of that vertex. For each grid cell, the recovered charge density values ​​obtained from the three vertices of the grid cell are linearly fitted using the least squares method to construct the linear recovered charge density function of the grid cell. To ensure the properties of the RWG basis function with respect to the charge on the positive and negative triangles in the triangle pair, the modified recovery charge density function of the grid cell is obtained by adding the linear recovery charge density function of the grid cell to a constant. Based on the corrected recovered charge density function of the mesh element and the surface charge density at the original superconvergence point, the local error of the mesh element is calculated using the following formulas: in, Indicates the first Local error of each grid cell; Indicates the first The corrected recovery charge density function for each grid cell; Indicates the first Surface charge density at the superconvergence point in each grid cell; Indicates the first Each grid cell.

2. The h-adaptive moment method for electromagnetic scattering according to claim 1, characterized in that, The initial mesh generation of the electromagnetic model under plane wave excitation includes: The electromagnetic model is The initial mesh is generated based on the dimensions; where, This represents the wavelength of the plane wave.

3. The h-adaptive moment method for electromagnetic scattering according to claim 2, characterized in that, The method of moments based on RWG basis functions is used to solve for the current coefficients of each grid cell in the current division, and the surface charge density at the superconvergence point in the corresponding grid cell is calculated using the solved current coefficients. This includes: The equivalent electromagnetic current of the target surface under test is discretized using the method of moments based on RWG basis functions and the currently divided grid. The basis functions generated after grid discretization are used to form a system of linear equations. The current coefficient of each grid cell is obtained by solving the system of linear equations. Using the formulas for calculating the current coefficient and surface charge density of each grid cell, the surface charge density at the superconvergence point in the corresponding grid cell is calculated.

4. The h-adaptive moment method for electromagnetic scattering according to claim 3, characterized in that, The calculation of global error based on the local error of each grid cell includes: For each grid cell, the local error of that grid cell is weighted by using the area of ​​its triangle as its own weight, and the weighted local error of that grid cell is obtained; and the sum of the weighted local errors of all grid cells is obtained, and the sum of the weighted local errors is obtained. The weighted local error is normalized using the sum of the triangle areas of all grid cells; The square root of the normalized result is used to obtain the global error.

5. The h-adaptive moment method for electromagnetic scattering according to claim 4, characterized in that, The step of determining whether the iteration stopping condition is met based on the global error includes: When the global error is reduced to a preset initial global error If the iteration count reaches the preset maximum iteration count, then the iteration stopping condition is determined to be met; otherwise, the iteration stopping condition is determined not to be met. It is a natural number greater than 1.

6. The h-adaptive moment method for electromagnetic scattering according to claim 5, characterized in that, The selection of partial mesh elements for mesh refinement based on the BW interpolation algorithm includes: The weighted local error threshold is obtained by multiplying the maximum weighted local error in all grid cells by a preset percentage. Mesh cells with a weighted local error greater than the weighted local error threshold are selected as the selected mesh cells; The selected grid cells are refined using a BW interpolation algorithm.

7. The h-adaptive moment method for electromagnetic scattering according to claim 6, characterized in that, Mesh refinement is performed on any selected mesh cell using a BW interpolation algorithm, including: In the current Delaunay triangulation The search found the point to be inserted. The base unit; where the current Delaunay triangulation This represents the currently defined triangular mesh; Using the adjacency relationship of mesh elements, mark all mesh elements that do not conform to the circumcircle criterion at the adjacent mesh elements of the base element, and denote the set of such mesh elements as . ; In the current Delaunay triangulation Delete The covered grid cells form a cavity. ,connect and The vertices on the boundary form new triangulations within the cavity. This yields a new Delaunay triangulation. .

Images