A multi-feature adaptive mesh generation method based on GPU-accelerated multi-scale FDTD

By using a GPU-accelerated multi-scale FDTD method, a linear hierarchical bounding box acceleration structure is constructed and a two-layer feature extraction is performed to generate an efficient and robust non-uniform mesh. This solves the problem of balancing accuracy and efficiency in complex multi-scale electromagnetic simulations in existing technologies, and achieves optimized allocation of computing resources and improved simulation accuracy.

CN122287258APending Publication Date: 2026-06-26ANHUI UNIV
View PDF 0 Cites 0 Cited by

Patent Information

Authority / Receiving Office
CN · China
Patent Type
Applications(China)
Current Assignee / Owner
ANHUI UNIV
Filing Date
2026-05-08
Publication Date
2026-06-26

AI Technical Summary

Technical Problem

Existing mesh generation techniques struggle to balance accuracy and computational resource consumption when handling complex multi-scale electromagnetic simulations. Furthermore, traditional adaptive meshing techniques cannot intelligently identify key electromagnetic features, resulting in low computational efficiency and resource waste.

Method used

A GPU-accelerated multi-scale FDTD method is adopted. By constructing a linear hierarchical bounding box acceleration structure, two-layer feature extraction and adaptive mesh generation are performed. Combined with four-stage post-processing correction, an efficient and robust non-uniform mesh is generated.

Benefits of technology

It significantly reduces the number of grids and computation time, while improving simulation accuracy and efficiency. It can intelligently identify key electromagnetic regions and perform adaptive encryption, solving the simulation challenges of complex models.

✦ Generated by Eureka AI based on patent content.

Smart Images

  • Figure CN122287258A_ABST
    Figure CN122287258A_ABST
Patent Text Reader

Abstract

This invention discloses a multi-feature adaptive mesh generation method based on GPU-accelerated multi-scale FDTD. By constructing a two-layer feature representation framework and combining it with GPU-accelerated BVH ray tracing technology, this invention can intelligently identify key electromagnetic regions and achieve adaptive mesh refinement. Compared with existing uniform mesh methods, while maintaining high-fidelity simulation accuracy, it can significantly reduce the number of meshes by 33% to 51% and reduce computation time by 52% to 70.1%. Simultaneously, a four-stage spatial connectivity algorithm ensures the robustness and numerical stability of non-uniform meshes, effectively avoiding numerical reflections caused by mesh abrupt changes, and solving the problem of balancing efficiency and accuracy in electromagnetic simulation of complex industrial models.
Need to check novelty before this filing date? Find Prior Art

Description

Technical Field

[0001] This invention relates to the fields of computational electromagnetics and numerical simulation technology, specifically to a multi-feature adaptive mesh generation method based on GPU-accelerated multi-scale FDTD. Background Technology

[0002] The finite-difference time-domain (FDTD) method, as a widely used numerical method for electromagnetic field calculation, plays a crucial role in handling complex multi-scale problems such as aircraft stealth design and electromagnetic compatibility analysis of integrated circuit boards. The quality of its mesh generation directly determines the accuracy and efficiency of the simulation. However, existing mesh generation techniques have significant limitations when dealing with large-scale models containing fine geometric features (such as sharp edges, slits, and micro-vias).

[0003] First, traditional uniform mesh generation strategies often struggle to balance approximation accuracy with computational resource consumption when fitting complex geometries. To reduce step approximation errors and accurately capture the electromagnetic properties of fine structures, the global mesh size is typically set to an extremely small value. This leads to a cubic explosion in the number of meshes, placing extremely high demands on computer memory and computation time, sometimes exceeding the capabilities of conventional workstations. Second, most existing adaptive meshing techniques focus on identifying purely geometric features, lacking in-depth consideration of electromagnetic physics mechanisms. Many methods refine the mesh based solely on simple geometric edges or surface roughness, failing to intelligently identify critical areas prone to electromagnetic singularities, such as regions with strong curvature or sharp-effect regions. This results in valuable computational resources being wasted in electromagnetically insensitive areas, reducing overall computational efficiency.

[0004] Furthermore, with the increasing complexity of industrial models, the number of model patches often reaches millions. Traditional geometric analysis algorithms based on CPU serial computation are inefficient in processing such large-scale data, becoming a bottleneck restricting the simulation process. Simultaneously, in the process of generating non-uniform meshes, without effective spatial connectivity control, distorted meshes with excessively large aspect ratios or topological breaks are easily generated in the transition regions between meshes of different densities. This not only affects the CFL stability condition of the FDTD algorithm but may also lead to computational divergence or result distortion. Therefore, developing an efficient mesh generation method that can intelligently sense key electromagnetic features, utilize modern hardware acceleration technology, and ensure mesh robustness has become an urgent technical challenge in the field of computational electromagnetics. Summary of the Invention

[0005] To address the technical problems mentioned above, this invention provides a multi-feature adaptive mesh generation method based on GPU-accelerated multi-scale FDTD, comprising the following steps: S1. Obtain the three-dimensional geometric model data to be simulated; S2. Based on the 3D geometric model data, construct a linear hierarchical bounding box acceleration structure on the GPU. S3. Based on the linear hierarchical bounding box acceleration structure, the model is sampled using GPU parallel threads, and two-layer feature extraction is performed to generate geometric and structural feature information of the model surface; wherein, the two-layer feature extraction includes the first layer of continuous geometric descriptor calculation and the second layer of discrete structural feature classification. S4. Based on geometric and structural features, project and map the three-dimensional priority field onto... x , y , z A one-dimensional linear density distribution is formed on the three coordinate axes, and the grid lines are automatically refined in high-priority regions according to the cumulative distribution function to generate an initial non-uniform grid. S5. Perform four-stage post-processing correction on the generated initial non-uniform mesh to obtain the final non-uniform mesh.

[0006] Preferably, in step S2, when constructing the linear hierarchical bounding box acceleration structure, the spatial geometry is divided into a hierarchical bounding box tree using the linear hierarchical bounding box algorithm to reduce the complexity of intersection testing for spatial ray tracing.

[0007] Preferably, the calculation of the first layer of continuous geometric descriptors in S3 further includes: Based on the distance between adjacent points, the difference between adjacent unit normals, and adjacent angles in the 3D geometric model data, calculate the curvature, sharpness, and local scale of the model surface, respectively. Curvature is used to characterize the degree of bending of the surface, sharpness is used to identify edges and corners, and local scale is used to quantify the physical dimensions of the geometry.

[0008] Preferably, the second-layer discrete structural feature classification in S3 further includes: Based on curvature, sharpness, and local scale, geometric regions are classified into geometric singularities, curvature-dominated connecting regions, and planar regions, thereby establishing a mapping relationship between geometric features and electromagnetic field singularities.

[0009] Preferably, in the process of generating the initial non-uniform mesh, step S4 introduces a scaling constraint algorithm, which includes the following steps: During the generation of grid lines, the length of the next step is determined based on the step size of the current region, and the adjacent edge ratio of each step is limited to a preset growth rate. Transition grid lines are inserted by forcibly implementing hierarchical proportional constraints between adjacent grid cells.

[0010] Preferably, the four-stage post-processing correction in S5 includes: For blank areas in the model that lack geometric features, a two-way balanced advancement strategy is adopted for processing; The bidirectional balanced advancement strategy is as follows: First, blank intervals are automatically identified based on vertex density detection. Then, a grid is generated synchronously from both ends of the blank interval towards the center. The grid size is increased exponentially and smoothly using a preset maximum geometric growth rate. Finally, the remaining gaps are precisely processed through a center intersection algorithm.

[0011] Preferably, the process of forming a one-dimensional linear density distribution in S4 further includes: Each coordinate axis direction is initially uniformly divided. Within each initially uniformly divided small interval, the density importance is obtained based on the number of vertices in the interval. This density importance is then integrated and weighted with the geometric-feature importance obtained in step S3 to obtain the final importance. ψ .

[0012] Preferably, the process of automatically encrypting the grid lines in S4 further includes: Based on the final importance ψ ,by τ high and τ low As the boundary, the importance is determined. ψ Convert to local step size h m This is used to divide the area into a fine-grid region corresponding to high-priority features, a medium-density region corresponding to the transition region, and a coarse-discrete region corresponding to the smooth manifold.

[0013] Compared with the prior art, the beneficial effects of the present invention are as follows: This invention constructs a two-layer feature representation framework and combines it with GPU-accelerated BVH ray tracing technology to intelligently identify key electromagnetic regions and achieve adaptive mesh refinement. Compared with existing uniform mesh methods, it can significantly reduce the number of meshes by 33% to 51% and reduce computation time by 52% to 70.1% while maintaining high-fidelity simulation accuracy. Simultaneously, a four-stage spatial connectivity algorithm ensures the robustness and numerical stability of non-uniform meshes, effectively avoiding numerical reflections caused by abrupt mesh changes, and solving the problem of balancing efficiency and accuracy in electromagnetic simulation of complex industrial models. Attached Figure Description

[0014] To more clearly illustrate the technical solution of the present invention, the drawings used in the embodiments are briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0015] Figure 1 This is the overall framework of the adaptive mesh generation strategy based on electromagnetic information proposed in this invention; Figure 2 This invention presents a hierarchical classification and decomposition of the geometric features Q1(v) of the UAV model proposed in this invention. Figure 3 The structural features J1, J2, J3, and Q2(v) hierarchical structure of the UAV model proposed in this invention are shown in the following diagram. Figure 4 This invention presents an important-driven adaptive step-size non-uniform mesh line generation process. Figure 5 This is a visualization of the hierarchical structure classification of the UAV model based on δ proposed in this invention; where (a) is a top view and (b) is a side view. Figure 6 This is the initial STL model for generating non-uniform meshes for a 10m-class large aircraft proposed in this invention; Figure 7 Voxelization visualization of the adaptive non-uniform mesh for the 10m-class large aircraft proposed in this invention; Figure 8 This is a comparative schematic diagram of the bistatic radar cross section in the 75MHz band proposed in this invention; wherein, (a) is a comparison of the non-uniform grid (467×453×147 grid points) with the FEKO reference solution; (b) is a comparison of the uniform grid (598×578×184 grid points) with the FEKO reference solution; Figure 9 This paper compares the relative errors between the RCS results of uniform and non-uniform meshes for the large aircraft model proposed in this invention and the FEKO results. Figure 10 The smartphone antenna model with multi-port feeding configuration proposed in this invention: sub-millimeter structural details at the 0.1 mm scale, used for S-parameter verification; Figure 11 This is a schematic diagram of the port current and reflected voltage waveforms of the present invention, wherein (a) is a comparison diagram of the port current of the adaptive non-uniform grid (787×853×112) and the uniform grid (1039×1010×132), and (b) is a comparison diagram of the reflected voltage waveforms of the adaptive non-uniform grid and the uniform grid. Figure 12 This is a comparison chart of the reflection coefficient S11 results of the non-uniform mesh model, the uniform mesh model, and the CST in the smartphone antenna model proposed in this invention. Detailed Implementation

[0016] 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.

[0017] To make the above-mentioned objects, features and advantages of the present invention more apparent and understandable, the present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments.

[0018] Example like Figure 1 The diagram shown is a schematic representation of the method flow in this embodiment, and the steps include: S1. Obtain the three-dimensional geometric model data to be simulated.

[0019] First, read the 3D geometric model data to be simulated (such as an STL format file) to obtain the vertex coordinates, face normals, and topological connectivity of the model.

[0020] S2. Based on the 3D geometric model data, construct a linear hierarchical bounding box acceleration structure on the GPU.

[0021] A linear hierarchical bounding box (LBVH) acceleration structure is constructed on the GPU side. The LBVH algorithm is used to divide the spatial geometry into a hierarchical bounding box tree, which reduces the complexity of intersection testing for spatial ray tracing.

[0022] S3. Based on the linear hierarchical bounding box acceleration structure, the model is sampled using GPU parallel threads to perform two-layer feature extraction, thereby generating geometric and structural feature information of the model surface.

[0023] In this embodiment, to quantify the electromagnetic importance of a region, the present invention defines a two-layer feature index. The first layer is a continuous geometric descriptor calculation, which is described by calculating three key parameters for each vertex of the model surface: curvature (average curvature of points around each vertex) to describe the degree of surface curvature; sharpness (standard deviation of curvature around each vertex) to describe the degree of irregularity or anisotropy of the curvature, thereby identifying edges and corners; and local scale (average edge length of the triangle connecting the vertices) to describe the physical size of the geometric features. The basic geometric information of the physical size of the geometric body is as follows: Normalized weighted average of the normal vectors within a ring at a vertex: ; Distance between adjacent points: ; Difference between adjacent unit normal vectors: ; Adjacent angles: ; Curvature: ; Sharpness: ; Local scale: .

[0024] in, nf Represents a piece of dough f The unit normal vector, This represents the set of triangular faces connected to a vertex. vi Indicates the coordinates of the current vertex. vj express vi The coordinates of the first-order neighborhood vertex, This represents the change in the normal vector between adjacent vertices. N This represents the total number of vertices in the first-order neighborhood of a vertex. Represents vertices vi With neighboring vertices vj The angle between the normal vectors, This represents the average angle between the normal vectors in the neighborhood of a vertex.

[0025] The geometric feature, one of the judgment conditions, is obtained by integrating the three geometric judgment parameters—curvature, sharpness, and local scale—according to their respective weighting coefficients. Q 1( v ): (1) in v Representing each vertex, κ norm , σ norm and s norm These represent normalized curvature, sharpness, and local scale, respectively. Weighting coefficients. w κ , w σ and w s The relative contributions of curvature, abrupt change features, and feature size were adjusted, with the sharpness weight typically set slightly higher to prioritize singularities in fields that meet the Meckner edge condition. The inverse scaling term [1- s This ensures that smaller features trigger correspondingly higher refinement priorities. Specifically... Q 1 ( v The distribution heatmap of complex UAV models is as follows: Figure 2 As shown.

[0026] The second layer involves discrete structural feature classification, categorizing geometric structures into geometric singularities, curvature-dominated connections, and regions dominated by planar patches outside these two categories. This approach emphasizes important structures while minimizing less important regions, thereby increasing mesh discriminability and avoiding over-subdivision. This serves as the basis for subsequent mesh refinement. In summary, these three categories form a hierarchical refinement priority scheme. J 1→ J 2→ J 3. This scheme allocates computational resources accordingly based on electromagnetic complexity. To quantitatively encode these priorities, this embodiment defines the importance of vertex-level structural features using a classification indicator function. Q 2( v ): (2) Specific examples and structural features of the three structures Q 2( v Heat map as follows Figure 3 As shown.

[0027] Depend on Figure 2 , Figure 3 The comparison shows the structural importance score Q 2( v As an important macroscopic supplement to the continuous geometric importance score Q1(v), it can combine fine-scale surface variations with advanced structural topologies. Figure 3 As shown, Q 1( v This describes the complexity of the local terrain, such as fuselage curvature and truss density; while Q 2( v This acts as a semantic filter, distinguishing key areas from a macroscopic perspective. Specifically, on propeller blades, Q 1( v It captures the overall surface curvature, while Q 2( v This algorithm rigorously refines the sharp tail edges to prevent feature loss during voxelization. For cylindrical motor housings, this algorithm will... Q 1( v The curvature of ) and J The two categories are combined to ensure that it is not for J 1. Accurate modeling of surface flow fields even with excessive refinement due to discontinuities. This two-layer approach dynamically allocates computational resources to the critical regions required for high-fidelity FDTD simulations based on geometric sensitivity and structural importance.

[0028] For large-scale models, this embodiment employs a GPU parallel computing architecture to improve analysis efficiency. First, a linear hierarchical bounding box (LBVH) is constructed on the GPU, reducing the time complexity of spatial search from...O ( N Reduced to O(log N Then, using GPU threads to fire rays in parallel, the distances and intersections from each point on the model surface to the mesh lines are quickly calculated. This data will be used for subsequent voxelization and mesh adjustment.

[0029] S4. Based on geometric and structural features, project and map the three-dimensional priority field onto... x , y , z A one-dimensional linear density distribution is formed on the three coordinate axes, and the grid lines are automatically refined in high-priority regions according to the cumulative distribution function to generate an initial non-uniform grid.

[0030] Specifically, first, an initial uniform division is performed in each direction, such as... Figure 4 As shown.

[0031] Within each initially uniformly divided subinterval, a normalized density importance is obtained based on the number of vertices within the subinterval. This density importance is then combined with the geometric-feature importance obtained above and weighted to obtain the final importance score. (3) in, V m Indicates the first m The set of all vertices contained within a spatial interval.

[0032] Importance ψ Convert to local step size h m High and low thresholds are respectively τ high , τ low To represent this, the high and low values ​​are typically taken as 0.7 and 0.4, as shown in formula (4), and the allowable step size variation is limited to the dynamic range. h min , h max Inside. Normally, h min It is determined by the dimensions of the finest local geometric features, and h max It is then limited by the wavelength in free space ( λ The constraint is / 20). This divides the region into a fine-mesh region with high-priority features, a moderately refined region in the transition zone, and a coarse-discrete region with a smooth manifold, such as... Figure 2 As shown.

[0033] (4) Grid lines are obtained by... x minarrive x max It is generated in a progressive manner. The step size of the current region. h m The length of the next step is determined, and the adjacent edge ratio of each step is restricted to a certain value. ρ Within. Figure 2 The final mesh structure is shown. The black dashed line represents the initial uniform mesh baseline, while the solid red line shows the adaptive mesh generated through iterative steps. It can be seen that the final generated region includes fine mesh areas, coarse mesh areas, and transition matrix areas.

[0034] To maintain numerical stability and reduce human-induced reflections, strict hierarchical scaling constraints are enforced between adjacent grid cells. ρ (Typically 1.2 to 1.3). This constraint manifests as a smooth change in mesh spacing in transition regions, preventing abrupt mesh changes while allowing for aggressive refinement in electromagnetically critical regions. The resulting non-uniform mesh concentrates computational resources where electromagnetic field variations require finer spatial sampling, while maintaining a stable, YEE-compatible Cartesian mesh structure. Figure 5 As shown, the final non-uniform mesh of the unmanned aerial vehicle is generated based on its spatial importance field ψ, where the gradient corresponds to the local feature mesh size ( ).

[0035] S5. Perform four-stage post-processing correction on the generated initial non-uniform mesh to obtain the final non-uniform mesh.

[0036] In this embodiment, the four stages are as follows: mesh validity verification and cleanup, precise geometric feature fitting and correction, smooth transition optimization of adjacent mesh sizes, and final confirmation of global mesh compliance. The first stage traverses the initial mesh structure, removes invalid and redundant voxels, corrects topological anomalies and element errors, and ensures the integrity and reliability of the basic mesh topology. The second stage performs local mesh adaptation calibration on key geometric feature areas such as edges and surfaces based on model vertex curvature, sharpness, and comprehensive importance parameters, ensuring accurate restoration of the original structural details. The third stage significantly weakens the abrupt size changes at the boundary between coarse and fine meshes through adjacent mesh step size ratio constraints, constructing a smooth and uniform mesh transition zone to avoid distortion and artifacts in the numerical simulation process. The fourth stage performs an overall verification of the rationality of the division of the global refinement zone, transition zone, and coarsening zone, finalizes the coordinates of all nodes and element attributes, and outputs a final non-uniformly adapted mesh that balances computational accuracy, computational efficiency, and numerical stability.

[0037] Next, for blank areas without geometric features in the model, this embodiment adopts a two-way balanced advancement strategy: First, blank areas are automatically identified based on vertex density detection. Then, meshes are generated synchronously from both ends of the area (i.e., the boundaries of adjacent geometries) towards the center. The mesh size is increased exponentially and smoothly using a preset maximum geometric growth rate (ρ), realizing a gradual transition from fine meshes at the boundaries to sparse meshes at the center. Finally, the remaining gaps are precisely processed through the center intersection algorithm to ensure that the meshes in the entire region meet strict adjacent size ratio constraints, effectively eliminating numerical reflection errors caused by mesh abrupt changes.

[0038] Figure 6 This study uses a 10m-class large aircraft model as an example. Three-dimensional electromagnetic scattering analysis is performed on a geometrically complex commercial aircraft model (characteristic dimensions 75.4 × 72.3 × 23.1 meters) to verify the proposed method. The computational domain employs a six-sided perfectly matched layer PML absorbing boundary condition to simulate unbounded free-space propagation. Bistatic RCS calculations are performed using both the proposed adaptive non-uniform mesh and a traditional uniform Cartesian mesh, and the results are compared with the benchmark solution of the commercial method of moments.

[0039] The adaptive mesh uses 467×453×147 grid points (total grid count 31,097,997). The aircraft model meshing results are as follows. Figure 7 As shown, the discretization results are visualized using voxelization. This distribution indicates that regions with high geometric curvature are preferentially discretized, while aerodynamically smooth fuselage sections remain coarser. Compared to a uniform grid (598×578×184 grid points, total grid count 63,598,496), the adaptive grid reduces the number of grid points by 51% and the computation time by 70.1% (from 6629.4 seconds to 1978.8 seconds).

[0040] Dual-station RCS comparison ( Figure 8 The results show that both meshing strategies are in good agreement with the FEKO benchmark solution in the range θ∈[0°,180°]. The mean absolute error (MAE) of the adaptive non-uniform mesh is 1.802 dB, which is close to the accuracy of the uniform mesh (MAE=1.99 dB). This high consistency indicates that the adaptive mesh can accurately capture the singular field behavior near the sharp edges and complex connections of the aircraft wing (these regions are crucial for high-frequency scattering). Despite the significant reduction in mesh size, non-uniform discretization still achieves better accuracy. Figure 9The relative error curves show that the adaptive mesh performs better in the low-frequency scattering range, indicating that for non-uniform mesh generation of 10-meter-class commercial aircraft structures, even compared with uniform meshes with more than twice the mesh density, it can still maintain superior detail resolution. This verifies that the method can achieve discriminative refinement of fine / coarse structural features without sacrificing accuracy, thereby significantly saving computational memory and establishing an optimal balance between geometric fidelity and computational efficiency for large-scale aerospace electromagnetic characterization applications.

[0041] Next, we will conduct S-parameter verification for a 0.1mm millimeter-scale multi-port antenna array (including fine structures such as feeding structures, microstrip discontinuities, and impedance matching networks). The model is as follows: Figure 10 As shown. The aim is to verify the applicability of the proposed method in the simulation of microscale electronic components, and the results are based on the commercial software CSTMicrowaveStudio.

[0042] In the simulation, the adaptive non-uniform mesh used 787×853×112 grid points (a total of 75,206,072 points) to adaptively refine highly complex regions, while coarsening was applied to regions with uniform substrates. The traditional uniform mesh required 1039×1010×132 grid points (a total of 138,519,480 points). Both methods underwent time-domain current and reflected voltage waveform tests, as well as S-parameter analysis of the reflection coefficient |S11| in the 0-7GHz frequency band. Results showed that the adaptive mesh reduced the number of grid points by 45.7% compared to the uniform mesh, and the computation time decreased from 17264.8 seconds to 8286.8 seconds, a reduction of 52.1%. Regarding accuracy, the port current and reflected voltage waveform results are as follows: Figure 11 As shown, it closely matches the uniform grid and can accurately capture features such as peaks and oscillations. A comparison of the relationship between S11 magnitude and frequency for the three benchmarks relative to the CST benchmark is shown in the figure below. Figure 12 As shown, the mean absolute error (MAE) of the adaptive mesh is 0.43 dB and the root mean square error (RMSE) is 0.58 dB, which are 69.5% and 71.8% higher than those of the uniform mesh (MAE=1.41 dB, RMSE=2.06 dB), respectively. It can accurately reproduce the characteristics of the primary resonance near 4 GHz and the secondary resonance of 5-6 GHz, fully demonstrating the advantages of this method in balancing efficiency and accuracy in microscale electromagnetic simulation.

[0043] The embodiments described above are merely preferred embodiments of the present invention and are not intended to limit the scope of the present invention. Various modifications and improvements made by those skilled in the art to the technical solutions of the present invention without departing from the spirit of the present invention should fall within the protection scope defined by the claims of the present invention.

Claims

1. A multi-feature adaptive mesh generation method based on GPU-accelerated multi-scale FDTD, characterized in that, Includes the following steps: S1. Obtain the three-dimensional geometric model data to be simulated; S2. Based on the 3D geometric model data, construct a linear hierarchical bounding box acceleration structure on the GPU. S3. Based on the linear hierarchical bounding box acceleration structure, the model is sampled using GPU parallel threads, and two-layer feature extraction is performed to generate geometric and structural feature information of the model surface; wherein, the two-layer feature extraction includes the calculation of the first layer of continuous geometric descriptors and the classification of the second layer of discrete structural features. S4. Based on geometric and structural features, project and map the three-dimensional priority field onto... x , y , z A one-dimensional linear density distribution is formed on the three coordinate axes, and the grid lines are automatically refined in high-priority regions according to the cumulative distribution function to generate an initial non-uniform grid. S5. Perform four-stage post-processing correction on the generated initial non-uniform mesh to obtain the final non-uniform mesh.

2. The multi-feature adaptive mesh generation method based on GPU-accelerated multi-scale FDTD according to claim 1, characterized in that, In S2, when constructing the linear hierarchical bounding box acceleration structure, the linear hierarchical bounding box algorithm is used to divide the spatial geometry into a hierarchical bounding box tree to reduce the intersection testing complexity of spatial ray tracing.

3. The multi-feature adaptive mesh generation method based on GPU-accelerated multi-scale FDTD according to claim 1, characterized in that, The calculation of the first-level continuous geometric descriptor in S3 further includes: Based on the distance between adjacent points, the difference between adjacent unit normals, and adjacent angles in the 3D geometric model data, calculate the curvature, sharpness, and local scale of the model surface, respectively. Curvature is used to characterize the degree of bending of the surface, sharpness is used to identify edges and corners, and local scale is used to quantify the physical dimensions of the geometry.

4. The multi-feature adaptive mesh generation method based on GPU-accelerated multi-scale FDTD according to claim 3, characterized in that, The second-level discrete structure feature classification in S3 further includes: Based on curvature, sharpness, and local scale, geometric regions are classified into geometric singularities, curvature-dominated connecting regions, and planar regions, thereby establishing a mapping relationship between geometric features and electromagnetic field singularities.

5. The multi-feature adaptive mesh generation method based on GPU-accelerated multi-scale FDTD according to claim 1, characterized in that, In the process of generating the initial non-uniform mesh, S4 introduces a scaling constraint algorithm, the steps of which include: During the generation of grid lines, the length of the next step is determined based on the step size of the current region, and the adjacent edge ratio of each step is limited to a preset growth rate. Transition grid lines are inserted by forcibly implementing hierarchical proportional constraints between adjacent grid cells.

6. The multi-feature adaptive mesh generation method based on GPU-accelerated multi-scale FDTD according to claim 1, characterized in that, The four-stage post-processing correction in S5 includes: For blank areas in the model that lack geometric features, a two-way balanced advancement strategy is adopted for processing; The bidirectional balanced advancement strategy is as follows: First, blank intervals are automatically identified based on vertex density detection. Then, a grid is generated synchronously from both ends of the blank interval towards the center. The grid size is increased exponentially and smoothly using a preset maximum geometric growth rate. Finally, the remaining gaps are precisely processed through a center intersection algorithm.

7. The multi-feature adaptive mesh generation method based on GPU-accelerated multi-scale FDTD according to claim 1, characterized in that, The process of forming a one-dimensional linear density distribution in S4 further includes: Each coordinate axis direction is initially uniformly divided. Within each initially uniformly divided small interval, the density importance is obtained based on the number of vertices in the interval. This density importance is then integrated and weighted with the geometric-feature importance obtained in step S3 to obtain the final importance. ψ .

8. The multi-feature adaptive mesh generation method based on GPU-accelerated multi-scale FDTD according to claim 7, characterized in that, The process of automatically encrypting grid lines in S4 further includes: Based on the final importance ψ ,by τ high and τ low As the boundary, assign importance ψ Convert to local step size h m This is used to divide the area into a fine-grid region corresponding to high-priority features, a medium-density region corresponding to the transition region, and a coarse-discrete region corresponding to the smooth manifold.