A deep learning-based reflective metasurface phase prediction method

Abstract

The application discloses a kind of based on deep learning's reflective metasurface phase prediction method.The method includes: using 16×16 discrete coding symmetric structure to construct metasurface unit, improve the design space of metasurface, and generate the data set including structure image and corresponding reflection phase real part and imaginary part curve by Python-CST joint simulation;A kind of hybrid neural network model is constructed by local feature extraction module, global feature modeling module and regression output module, wherein local feature extraction is realized based on ResNet, and global feature modeling is completed by introducing position coding and frequency coding Transformer encoder;While designing the complex loss function including phase real part error, imaginary part error and physical constraint term, and the weight coefficient thereof is systematically optimized to balance prediction accuracy and physical rationality;Optimized loss function and set hyperparameter are used to train the network;Finally, the trained model is used to realize the rapid phase prediction of new metasurface structure.The application effectively solves the problems of phase jump processing difficulty in high degree of freedom metasurface design, and the lack of local and global feature collaborative modeling, realizes the rapid, high-precision and physically reasonable prediction of reflection phase in wide frequency band.

Description

Technical Field

[0001] This invention relates to the intersection of electromagnetic metasurface design and artificial intelligence technology, specifically to a method for rapid prediction of the reflection phase of a reflective metasurface unit structure based on a deep learning model. Background Technology

[0002] With the widespread application of electromagnetic metasurfaces in radar stealth, new antennas, electromagnetic shielding and other fields, their design complexity and performance requirements are increasing. Traditional metasurface design methods mainly rely on full-wave electromagnetic simulation software for parametric simulation and optimization, which has the following problems: (1) long design cycle, with a single simulation taking anywhere from several minutes to several hours; (2) the structural design space is limited by the parametric modeling method, making it difficult to achieve the design of high-degree-of-freedom, irregular structures; (3) the phase response curve will show abrupt changes near the resonant frequency, affecting the design stability; (4) there is a lack of rapid prediction tools that can simultaneously capture the local structure and global electromagnetic coupling relationship.

[0003] In recent years, deep learning technology has provided new ideas for the reverse design of metasurfaces. The electromagnetic metasurface design method based on deep learning can realize rapid prediction and reverse design by constructing a proxy model from structure to electromagnetic response. However, existing methods still have obvious limitations: (1) Most studies use parametric geometric structures (such as Jerusalem cross, ring, square ring, etc.), which have limited design freedom and are difficult to meet the requirements of complex electromagnetic functions; (2) When directly predicting the reflection phase curve, the phase change phenomenon that may occur in the frequency band is not handled well, which affects the stability of model training and prediction accuracy; (3) The network structure mostly uses pure convolutional neural network or pure Transformer, which is difficult to capture the local detailed features of the metasurface structure and the global electromagnetic coupling relationship at the same time, and the generalization ability is insufficient when dealing with broadband and high degree of freedom design tasks; (4) In broadband prediction tasks, the number of network parameters is large, the training difficulty is high, and the prediction accuracy is limited.

[0004] Therefore, there is an urgent need for a metasurface phase prediction method that can handle high-degree-of-freedom structures, stably predict broadband phase responses, and combine high accuracy and high efficiency. Summary of the Invention

[0005] The purpose of this invention is to overcome the shortcomings of existing technologies and provide a deep learning-based method for predicting the phase response of reflective metasurfaces. This method employs a high-degree-of-freedom discrete coding structure, an innovative data preprocessing approach, and a neural network architecture that combines local and global feature extraction capabilities. Furthermore, it customizes an attention mechanism to address the axisymmetric characteristics of the metasurface unit, thereby achieving rapid and accurate prediction of the phase response of reflective metasurface units.

[0006] To achieve the above objectives, the technical solution provided by the present invention is as follows:

[0007] A deep learning-based method for predicting the phase of a reflective metasurface includes the following steps:

[0008] Step S1: A 16×16 discrete coded symmetric structure is used as the metasurface unit. This unit consists of three layers: a top discrete metal patch, an intermediate dielectric layer, and a bottom metal ground plane.

[0009] Step S2 involves using a Python program in conjunction with CST electromagnetic simulation software to generate a large number of random structure samples and their corresponding reflection phase curves, thereby constructing a sufficient dataset.

[0010] Step S3: A phase prediction neural network model is constructed using a hybrid architecture that combines a convolutional neural network (CNN) and a Transformer encoder.

[0011] Specifically, a three-stage hybrid architecture is adopted.

[0012] The first stage is the local feature extraction module, implemented based on the ResNet-18 network. The input metasurface structure image first passes through a 7×7 convolutional layer and a max-pooling layer, and then performs deep feature extraction through four residual blocks. Each residual block contains two 3×3 convolutional layers, equipped with batch normalization and ReLU activation functions. The number of output channels for the four residual blocks are 64, 128, 256, and 512, respectively, with the spatial size gradually decreasing to 4×4. The local features extracted in this stage have clear physical correspondences: the features extracted by the shallow layers of the network (the first two residual blocks) mainly correspond to the geometric features of a single metal patch and its neighboring regions, such as the presence or absence of the patch, edge direction, and neighborhood density; the features extracted by the deep layers of the network (the last two residual blocks) correspond to a larger range of periodic structural patterns and spatial distribution rules, such as the regional distribution of metal content and the degree of symmetry. This feature extraction process from shallow to deep is consistent with the physical law of "near-field coupling decays with distance" in electromagnetic field theory, enabling the network to understand metasurface structures at multiple levels from micro to macro.

[0013] The second stage is the global feature modeling module. The 4×4×512 feature tensor output by ResNet is reshaped into a 16×512 sequence. It's important to note that each position in the 16×512 sequence (a total of 16 spatial points) corresponds to a 4×4 sub-region within the metasurface unit (i.e., a 4×4 metal patch block in the original 16×16 structure). The distribution pattern of the metal patches within this sub-region determines the local electromagnetic coupling characteristics. Therefore, reshaping the feature map into a sequence essentially transforms the local electromagnetic response units of the metasurface unit into sequential elements. The physical meaning of Transformer's self-attention mechanism here is to calculate the electromagnetic coupling strength between any two local blocks. Since electromagnetic waves propagate within the metasurface unit, there is electromagnetic interaction (including near-field coupling and far-field interference) between any two metal patches, this interaction mathematically manifests as a non-linear global dependency. Traditional pure CNN models are limited by a fixed receptive field, making it difficult to directly model the coupling between blocks that are far apart. In contrast, Transformer uses a self-attention weight matrix to explicitly calculate the electromagnetic influence coefficient of each block on all other blocks, thereby achieving direct modeling of global electromagnetic coupling relationships.

[0014] To preserve spatial location information, a two-dimensional learnable location code is added to the sequence. Simultaneously, to inject frequency information into the network, a frequency encoding vector is added, representing the normalized frequency values ​​of 20 equally spaced frequency points within the target frequency band. The processed sequence is then input into a module consisting of four Transformer encoder layers. Each Transformer encoder layer includes a multi-head self-attention mechanism and a feedforward neural network, with 8 attention heads and a 2048-dimensional hidden layer in the feedforward network.

[0015] Based on this, and considering the forced axisymmetry of the metasurface unit used in this invention, the self-attention mechanism of the standard Transformer is customized and improved, proposing a symmetry-aware self-attention module. Specifically, for the 16 sequence elements corresponding to 4×4 spatial positions, their symmetry equivalence classes are defined according to their spatial coordinates (i,j). Since the structure is symmetric about both the x-axis and y-axis, any position (i,j) is physically equivalent to its symmetric positions (15-i,j), (i,15-j), and (15-i,15-j). Based on this, when calculating the attention weights, the following condition is forcibly satisfied: A (i,j)→(p,q) = A (i’,j’)→(p’,q’)Here, (i',j') and (p',q') are the positions corresponding to (i,j) and (p,q) after symmetric transformation, respectively. The implementation is as follows: Before calculating Query and Key, symmetric pooling is performed on the input features. The features of each 4×4 block and its three symmetric blocks are averaged to obtain 8 independent "symmetric invariants" (since the 16 positions can be divided into 8 groups according to symmetry). Then, self-attention is calculated in the 8×8 compressed space, and finally, the attention weights are symmetrically mapped back to the original 16×16 attention matrix. The advantages of this design are: fewer parameters, reducing the computational complexity of the attention matrix from O(16²) to O(8²); physical constraint embedding, forcing the network to learn the symmetric relationships that inevitably exist under symmetric structures, avoiding the network learning this rule from the data itself; improved generalization ability, with symmetry as a hard constraint, making the model more robust to small structural perturbations during testing.

[0016] The third stage is the regression output module. The sequence output by the Transformer is reconstructed back into feature map form, and then downsampled through a three-layer convolutional network. Each convolutional layer is followed by batch normalization and ReLU activation. Normalizing the input data for each batch in the network improves training speed and alleviates the vanishing gradient problem, thus enhancing model training stability. ReLU is chosen as the activation function because its gradient is always 1 in the positive interval, effectively mitigating the vanishing gradient problem in deep networks and making the model easier to train. Next, global average pooling is performed to convert the feature map into feature vectors. Finally, two independent fully connected layers output 20 predicted points for the real and imaginary parts of the reflection phase, respectively.

[0017] Step S4: The optimizer selects the Adam optimizer and sets appropriate hyperparameters such as learning rate and weight decay to train and optimize the network. Specifically, the training process of the phase prediction network is as follows: input a binarized hypersurface grayscale image, the prediction network propagates forward while predicting the discrete values ​​on the real and imaginary curves of the reflection phase, calculates the loss between the prediction network and the discrete values ​​on the real and imaginary curves of the actual reflection phase, calculates the gradient of the loss and performs gradient backpropagation to update the weights, and iterates the above process until the loss function converges.

[0018] Step S5: Perform phase prediction using the trained model. For a new, unseen binary image of a metasurface structure, input it into the trained phase prediction network to obtain the corresponding prediction curves of the real and imaginary parts of the reflection phase through forward propagation. Using inverse trigonometric functions, a complete reflection phase curve can be further synthesized.

[0019] Preferably, step S2 further includes step S21, creating a zero matrix of size 8×8. .

[0020] Step S22: Each element in the matrix is ​​independently and randomly assigned a value. A uniform distribution is used, setting an element to 1 with a probability of 0.5 and setting it to 0 with a probability of 0. Mathematically, this is expressed as:

[0021] in

[0022] Step S23, to obtain a fully axisymmetric structure, for Perform two mirror symmetry operations:

[0023] First, perform horizontal mirror symmetry to generate an intermediate matrix. :

[0024]

[0025] Then, perform vertical mirror symmetry to generate the final matrix.

[0026]

[0027] After generation In the matrix, the element distribution satisfies both central and axial symmetry, which can effectively suppress cross-polarization components.

[0028] Step S24: Import the random structure from step S23 into CST software for reflection phase curve calculation. A frequency domain solver is used during simulation. Boundary conditions in the X and Y directions are set as periodic unit boundaries, and the boundary condition in the Z direction is set as an open boundary with added space. The electromagnetic wave is set as a vertically incident plane wave with a frequency range of 8-12 GHz.

[0029] Step S25: For the simulation-obtained reflection coefficients S11, to eliminate the adverse effects of potential phase jumps in the reflection phase curve on neural network training, the reflection coefficients are converted into curves representing the real and imaginary parts of the reflection phase as training labels. The conversion relationship between reflection phase, amplitude, and the real and imaginary parts of the reflection phase is as follows:

[0030]

[0031]

[0032] Step S26: Sample at fixed intervals within the target frequency band and perform max-min normalization on the real and imaginary data.

[0033] Step S27: Convert each 16×16 binary matrix into a 128×128 pixel grayscale image, where 0 values ​​correspond to black (RGB: 0,0,0) and 1 values ​​correspond to white (RGB: 255,255,255).

[0034] Step S28: A total of 40,000 valid data sets are generated, and the dataset is split into training, validation, and test sets in a 7:2:1 ratio. Each data set contains: a 16×16 binary structure image (as input), and the corresponding 81-dimensional phase real vector and 81-dimensional phase imaginary vector (as labels).

[0035] Preferably, step S3 further includes step S31, where the loss function is trained using the following composite loss function:

[0036]

[0037] in Represents the phase vector, subscript The subscript represents the predicted value of the phase prediction network. λ1, λ2, and λ3 represent the true label values ​​of the data, and λ3 are the corresponding weighting coefficients. The first term of the loss function represents the L1 loss between the predicted and true values ​​of the real part of the reflection phase; the second term represents the L1 loss between the predicted and true values ​​of the imaginary part of the reflection phase; and the third term represents the requirement that the predicted real and imaginary parts satisfy a trigonometric mathematical relationship. By introducing a constraint term based on trigonometric identities, the predicted real and imaginary parts are forced to satisfy a basic mathematical relationship, which can improve the physical rationality of the prediction results.

[0038] Preferably, step S3 further includes step S32, which systematically optimizes λ1, λ2, and λ3 to balance prediction accuracy and physical rationality. The optimization method is as follows: (1) Conduct a univariate scanning experiment to preliminarily determine the effective range of each coefficient; (2) Conduct a fine grid search within the effective range to evaluate the performance of the validation set under different coefficient combinations; (3) Construct a response surface model based on the grid search results to locate the theoretical optimal solution region; (4) Finally, through multi-objective trade-off analysis, achieve the best balance between prediction accuracy and physical constraint satisfaction to determine the optimal weight coefficient combination.

[0039] Preferably, step S4 further includes step S41, where training employs mini-batch gradient descent with a batch size of 128. The optimizer uses the Adam algorithm, setting appropriate parameters such as the initial learning rate, the exponential decay rate of the first moment estimate, the exponential decay rate of the second moment estimate, and the weight decay coefficient. When the amount of data is small, the model is prone to overfitting; weight decay regularization can be used in this case. The learning rate is scheduled using a cosine annealing strategy.

[0040] This method can complete predictions at an extremely fast speed (milliseconds), which is a significant improvement in efficiency compared to traditional simulation methods (minutes to hours).

[0041] The beneficial effects of this invention are as follows: (1) By adopting a discrete coding structure, the design freedom of metasurfaces is greatly expanded, providing a foundation for the design of complex functional metasurfaces. (2) By predicting the real and imaginary parts of the phase instead of directly predicting the phase angle, the training difficulties caused by phase jumps are effectively avoided, and the stability and accuracy of the model are improved. (3) The proposed CNN-Transformer hybrid network can simultaneously utilize the local feature extraction and inductive bias of CNN and the global relation modeling ability of Transformer, significantly enhancing the model's learning and generalization ability for the electromagnetic response of high-degree-of-freedom structures. (4) For the forced axisymmetric characteristics of metasurface units, a symmetry-aware self-attention mechanism and a frequency-space joint attention mechanism are designed, further enhancing the model's ability to model structural priors and frequency dependencies, making the prediction results more in line with physical laws. (5) The entire method is efficient, and the trained model can achieve millisecond-level phase prediction, greatly accelerating the design, analysis and optimization cycle of metasurfaces. Attached Figure Description

[0042] To more clearly illustrate the technical solutions of the embodiments of the present invention, the accompanying drawings used in the embodiments will be briefly introduced below. It should be understood that the following drawings only show some embodiments of the present invention and should not be regarded as a limitation on the scope. For those skilled in the art, other related drawings can be obtained based on these drawings without creative effort.

[0043] Figure 1 This is a schematic diagram of the 16×16 discrete coded symmetric metasurface unit structure used in the embodiments of the present invention.

[0044] Figure 2 The overall architecture of the reflective metasurface phase prediction network model provided by this invention.

[0045] Figure 3 This is a graph showing the change of the loss function during the network training process in an embodiment of the present invention.

[0046] Figure 4 This is a visualization example of three sets of phase prediction results on the test set in an embodiment of the present invention.

[0047] Figure 5 This is a scatter plot of the absolute phase error of 100 samples in the test set in this embodiment of the invention.

[0048] Figure 6 This is an analytical histogram of the absolute error of phase prediction for all samples in the test set in this embodiment of the invention. Detailed Implementation

[0049] To make the objectives, technical solutions, and advantages of this invention clearer, the invention will be further described in detail below with reference to the accompanying drawings and specific embodiments. It should be noted that the following embodiments are only for explaining the invention and are not intended to limit the scope of protection of the invention.

[0050] Example

[0051] This embodiment provides a phase prediction method for reflective metasurfaces in the 8-12 GHz frequency band.

[0052] Step S1, please refer to Figure 1 This embodiment employs a three-layer reflective metasurface unit with a unit period p set to 16 mm. The top layer, an 8×8 random matrix, generates a 16×16 axisymmetric metal patch pattern through mirror symmetry. The material is copper, with a thickness of 0.035 mm. The dielectric layer is made of F4B with a thickness of 2 mm. The bottom layer is a PEC ground plane. In the figure, yellow areas indicate locations with metal blocks, while blue areas indicate locations without metal blocks.

[0053] Step S2: Through joint simulation using Python program and CST electromagnetic simulation software, a large number of random structure samples and their corresponding reflection phase curves are generated in batches, thereby constructing a dataset of 40,000 sets.

[0054] Specifically, step S2 also includes step S21, creating a zero matrix of size 8×8.

[0055] Specifically, step S2 further includes step S22, which involves independently assigning a random value to each element in the matrix. The assignment uses a uniform distribution, setting an element to 1 with a probability of 0.5 and setting it to 0 with a probability of 0. Mathematically, this can be expressed as:

[0056] in

[0057] Specifically, step S2 also includes step S23, which, in order to obtain a fully axisymmetric structure, involves... Perform two mirror symmetry operations:

[0058] First, perform horizontal mirror symmetry to generate an intermediate matrix. :

[0059]

[0060] Then, perform vertical mirror symmetry to generate the final matrix.

[0061]

[0062] After generation In the matrix, the element distribution satisfies both central and axial symmetry, which can effectively suppress cross-polarization components.

[0063] Specifically, step S2 also includes step S24, which imports the random structure from step S23 into the CST software for reflection phase curve calculation. During simulation, a frequency domain solver is used, with boundary conditions in the X and Y directions set as periodic unit boundaries, and boundary conditions in the Z direction set as open boundaries with added space; the electromagnetic wave is set as a vertically incident plane wave with a frequency range of 8-12 GHz.

[0064] Specifically, step S2 also includes step S25, which extracts the real and imaginary parts of the simulation-obtained reflection coefficient S11 to obtain the phase real and imaginary part curves. The real and imaginary part curves are smooth and continuous, which is more conducive to neural network processing. The conversion relationship between the reflection phase, amplitude, and the real and imaginary parts of the reflection phase is as follows:

[0065]

[0066]

[0067] Specifically, step S2 also includes step S26, which involves sampling 81 points at 0.05 GHz intervals within 8-12 GHz of each curve, and performing max-min normalization on the real and imaginary data.

[0068] Specifically, step S2 also includes step S27, which converts each 16×16 binary matrix into a 128×128 pixel grayscale image, where 0 values ​​correspond to black (RGB: 0,0,0) and 1 values ​​correspond to white (RGB: 255,255,255).

[0069] Specifically, step S2 also includes step S28, which generates a total of 40,000 sets of valid data and splits the dataset into training, validation, and test sets in a 7:2:1 ratio. Each set of data contains: a 16×16 binary structure image (as input), and the corresponding 81-dimensional phase real vector and 81-dimensional phase imaginary vector (as labels).

[0070] Step S3, please refer to Figure 2 A phase prediction network architecture was constructed, employing a three-stage hybrid architecture. Detailed configurations of each layer are shown in Table 1. Specifically:

[0071] The first stage of the network is a local feature extraction module, implemented based on the ResNet-18 network. The input metasurface structure image first passes through a 7×7 convolutional layer and a max-pooling layer, and then performs deep feature extraction through four residual blocks. Each residual block contains two 3×3 convolutional layers, equipped with batch normalization and ReLU activation functions. The number of output channels for the four residual blocks are 64, 128, 256, and 512, respectively, with the spatial size gradually decreasing to 4×4. The local features extracted in this stage have clear physical correspondences: the features extracted in the shallow layers of the network mainly correspond to the geometric features of a single metal patch and its neighboring regions, such as the presence or absence of the patch, edge direction, and neighborhood density; the features extracted in the deeper layers correspond to a larger range of periodic structural patterns and spatial distribution rules, such as the regional distribution of metal content and the degree of symmetry. This feature extraction process from shallow to deep is consistent with the physical law in electromagnetic field theory that "near-field coupling decays with distance."

[0072] The second stage is the global feature modeling module. The 4×4×512 feature tensor output by ResNet is reshaped into a 16×512 sequence. It's important to note that each position in the 16×512 sequence (a total of 16 spatial points) corresponds to a 4×4 sub-region within the metasurface unit (i.e., a 4×4 metal patch block in the original 16×16 structure). The distribution pattern of the metal patches within this sub-region determines the local electromagnetic coupling characteristics. Therefore, reshaping the feature map into a sequence essentially transforms the local electromagnetic response units of the metasurface unit into sequential elements. The physical meaning of Transformer's self-attention mechanism here is to calculate the electromagnetic coupling strength between any two local blocks. Since electromagnetic waves propagate within the metasurface unit, there is electromagnetic interaction (including near-field coupling and far-field interference) between any two metal patches, this interaction mathematically manifests as a non-linear global dependency. Traditional pure CNN models are limited by a fixed receptive field, making it difficult to directly model the coupling between blocks that are far apart. In contrast, Transformer uses a self-attention weight matrix to explicitly calculate the electromagnetic influence coefficient of each block on all other blocks, thereby achieving direct modeling of global electromagnetic coupling relationships.

[0073] To preserve spatial location information, a two-dimensional learnable location code is added to the sequence. Simultaneously, to inject frequency information into the network, a frequency encoding vector is added, representing the normalized frequency values ​​of 20 equally spaced frequency points within the target frequency band. The processed sequence is then input into a module consisting of four Transformer encoder layers. Each Transformer encoder layer includes a multi-head self-attention mechanism and a feedforward neural network, with 8 attention heads and a 2048-dimensional hidden layer in the feedforward network.

[0074] Based on this, and considering the forced axisymmetry of the metasurface unit used in this invention, the self-attention mechanism of the standard Transformer is customized and improved, proposing a symmetry-aware self-attention module. Specifically, for the 16 sequence elements corresponding to 4×4 spatial positions, their symmetry equivalence classes are defined according to their spatial coordinates (i,j). Since the structure is symmetric about both the x-axis and y-axis, any position (i,j) is physically equivalent to its symmetric positions (15-i,j), (i,15-j), and (15-i,15-j). Based on this, when calculating the attention weights, the following condition is forcibly satisfied: A (i,j)→(p,q) = A (i’,j’)→(p’,q’) Here, (i',j') and (p',q') are the positions corresponding to (i,j) and (p,q) after symmetric transformation, respectively. The implementation is as follows: Before calculating Query and Key, symmetric pooling is performed on the input features. The features of each 4×4 block and its three symmetric blocks are averaged to obtain 8 independent "symmetric invariants" (since the 16 positions can be divided into 8 groups according to symmetry). Then, self-attention is calculated in the 8×8 compressed space, and finally, the attention weights are symmetrically mapped back to the original 16×16 attention matrix. The advantages of this design are: fewer parameters, reducing the computational complexity of the attention matrix from O(16²) to O(8²); physical constraint embedding, forcing the network to learn the symmetric relationships that inevitably exist under symmetric structures, avoiding the network learning this rule from the data itself; improved generalization ability, with symmetry as a hard constraint, making the model more robust to small structural perturbations during testing.

[0075] The third stage is the regression output module. The sequence output by the Transformer is reconstructed back into feature map form, then downsampled through a three-layer convolutional network, with batch normalization and ReLU activation following each convolutional layer. Next, global average pooling is performed to convert the feature map into feature vectors. Finally, two independent fully connected layers output 20 predicted points for the real and imaginary parts of the reflection phase, respectively.

[0076] Table 1 Detailed configuration of each layer of the phase prediction network

[0077]

[0078] Specifically, step S3 further includes step S31, designing the loss function as a sum of three terms. The first term is the phase real part error term, calculating the L1 loss between the predicted and true phase real parts; the second term is the phase imaginary part error term, calculating the L1 loss between the predicted and true phase imaginary parts; and the third term is a physical constraint term, forcing the predicted real and imaginary parts to satisfy the basic trigonometric identities. The specific expression is as follows:

[0079]

[0080] Specifically, step S3 also includes step S32, which systematically optimizes λ1, λ2, and λ3 to balance prediction accuracy and physical rationality. The optimization method is achieved through the following systematic experimental steps: (1) conducting univariate scanning experiments to initially determine the effective range of each coefficient; (2) conducting a fine grid search within the effective range to evaluate the performance of the validation set under different coefficient combinations; (3) constructing a response surface model based on the grid search results to locate the theoretical optimal solution region; (4) finally determining the optimal weight coefficient combination by achieving the best balance between prediction accuracy and physical constraint satisfaction through multi-objective trade-off analysis.

[0081] Step S4 involves training and optimizing the network. Specifically, the training process for the phase prediction network involves inputting a binarized metasurface grayscale image, propagating the network forward to predict the discrete values ​​of the real and imaginary parts of the reflection phase on the curves, calculating the loss between this prediction and the discrete values ​​of the actual reflection phase on the curves, calculating the gradient of the loss, and performing backpropagation to update the weights. This process is iterated until the loss function converges. Training was performed on a Windows 11 operating system, with an Intel Core i7-13700F processor and an NVIDIA GeForce RTX 4070 graphics card. The software environment was based on the Anaconda platform, using Python version 3.6.13 and the PyTorch 1.10.0 deep learning framework.

[0082] Specifically, step S4 also includes step S41, where training uses mini-batch gradient descent with a batch size of 128. The optimizer chosen is the Adam algorithm, with parameters set as follows: initial learning rate 1×10⁻ 4 The first-order moment estimates the exponential decay rate β1 = 0.9, the second-order moment estimates the exponential decay rate β2 = 0.999, and the weight decay coefficient is 1 × 10⁻ 4 Considering the limited amount of data in this paper, and the potential for overfitting in the model, a weight decay regularization method is used, with the regularization coefficient set to 10⁻. 4 The learning rate is scheduled using a cosine annealing strategy, and the total number of training cycles is set to 500. Table 2 shows all the hyperparameter settings for the phase prediction network.

[0083] Table 2 Hyperparameter settings in phase prediction networks

[0084]

[0085] Step S5: Perform phase prediction using the trained model. For a new, unseen binary image of a metasurface structure, input it into the trained phase prediction network to obtain the corresponding prediction curves of the real and imaginary parts of the reflection phase through forward propagation. Using inverse trigonometric functions, a complete reflection phase curve can be further synthesized.

[0086] Please see Figure 3 As training progresses, the loss on the training set gradually decreases, converging around the 200th iteration. The loss on the validation set initially increases, indicating poor generalization ability. However, with further training, the loss on the validation set begins to decrease, also approaching convergence around the 200th iteration. Furthermore, after convergence, the average loss on the validation set is slightly lower than the average loss on the training set, demonstrating good generalization ability of the network.

[0087] The trained model was evaluated on an independent test set containing 4000 samples. See also... Figure 4 The predicted results (scatter points) of the three random test samples are in high agreement with the actual curve (solid line).

[0088] Please see Figure 5 and Figure 6 The error scatter plot of 100 random samples and the histogram of the error of all 4000 samples were calculated. The prediction error of the vast majority of samples was concentrated between 2° and 5°, which fully demonstrates the high accuracy and reliability of the method proposed in this invention.

[0089] In summary, this embodiment verifies the effectiveness of the method of the present invention. This method can rapidly and accurately predict the broadband reflection phase response of a given high-degree-of-freedom reflective metasurface structure within milliseconds, providing a powerful tool for efficient reverse engineering and performance analysis of metasurfaces.

[0090] The above detailed embodiments fully illustrate the technical solution, implementation steps, and expected effects of the present invention. Those skilled in the art can make appropriate adjustments and modifications to the present invention based on the content of this specification and in conjunction with specific design requirements; all such adjustments and modifications should be included within the protection scope of the present invention.

Claims

1. A deep learning-based method for predicting the phase of a reflective metasurface, characterized in that, Includes the following steps: S1. Generate a reflective metasurface dataset: Design a 16×16 discrete coded symmetric structure as a metasurface unit, obtain the reflection coefficient data corresponding to the metasurface unit through automated co-simulation using python-CST, and convert the reflection coefficient data into real part curve data and imaginary part curve data of the reflection phase to form a dataset containing metasurface structure image, real part curve data and imaginary curve data. S2. Construct and train a phase prediction neural network: Construct a hybrid network model that includes a convolutional neural network module and a Transformer encoder module, and train the hybrid network model using the dataset; S3. Using the trained hybrid network model, predict the input target metasurface unit structure image and output the corresponding prediction curves of the real and imaginary parts of the reflection phase.

2. The method according to claim 1, characterized in that, In step S1, the design of the 16×16 discrete coding symmetric structure specifically includes: First, generate an 8×8 binary random matrix M8×8, where the values ​​of the matrix elements are randomly selected between 0 and 1; Perform mirror symmetry operations on the matrix M8×8 in the horizontal and vertical directions respectively to obtain a fully axisymmetric 16×16 matrix M16×16, where the position with a value of 1 represents the presence of a metal patch.

3. The method according to claim 1, characterized in that, In step S1, converting the reflection coefficient data into real and imaginary part curve data of the reflection phase specifically involves: extracting the real part (Re(S)) of the complex reflection coefficient. 11 The imaginary part (Im(S)) is extracted from the real part curve data. 11 ()) is used as the data for the imaginary part of the curve.

4. The method according to claim 1, characterized in that, In step S2, the hybrid network model adopts a three-stage architecture: The first stage is a local feature extraction module based on residual networks, which is used to extract features from the input metasurface structure image; The second stage is a global electromagnetic coupling modeling module based on a Transformer encoder, which is used to serialize the features output from the first stage and calculate the electromagnetic coupling strength between any two local regions within the metasurface unit through a self-attention mechanism. The input of the module integrates position coding and frequency coding. The third stage is the regression output module, which decodes the features processed in the second stage and outputs the predicted values ​​of the real and imaginary parts of the reflection phase, respectively.

5. The method according to claim 4, characterized in that, The Transformer encoder module employs a symmetry-aware self-attention mechanism. Based on the axisymmetric structure prior of the metasurface unit, it divides the spatial position into equivalence classes according to symmetry, performs attention calculation on the compressed symmetric feature space, and restores the calculation result to the original spatial dimension through symmetry mapping.

6. The method according to claim 4, characterized in that, The frequency coding is used to inject frequency information within the target frequency band into the network, and the location coding is used to inject spatial location information of the image into the network.

7. The method according to claim 1, characterized in that, In step S2, the hybrid network model is trained using a composite loss function, which is: Where, φ p For the phase predicted by the network, φ t For the true phase, Let λ1 represent the L1 norm, and λ1, λ2, and λ3 be weighting coefficients determined through systematic optimization. The first term of the loss function represents the L1 loss between the predicted and true values ​​of the real part of the reflection phase; the second term represents the L1 loss between the predicted and true values ​​of the imaginary part of the reflection phase; and the third term represents the requirement that the predicted real and imaginary parts satisfy a trigonometric relationship. By introducing constraint terms based on trigonometric identities, the predicted real and imaginary parts are forced to satisfy a basic mathematical relationship, which improves the physical plausibility of the prediction results.

8. The method according to claim 7, characterized in that, The optimization method for the weight coefficients λ1, λ2, and λ3 includes: determining the effective range through univariate scanning experiments, performing grid search within the effective range, constructing a response surface model based on the search results to locate the theoretically optimal region, and determining the final coefficient combination through multi-objective trade-off analysis.

9. The method according to claim 1, characterized in that, In step S2, when training the hybrid network model, the optimizer is Adam, the learning rate is set to 1e-4, the weight decay coefficient is set to 1e-4, and the batch size is set to 128.

Images