Open-domain electromagnetic target solving method and system for unidirectional fine structure

Abstract

The invention provides an open-domain electromagnetic target solving method and system for a one-way fine structure, and the method comprises the steps: introducing a scaling factor of a PML into a Maxwell equation, so as to construct a time evolution matrix of an electromagnetic field component; a mixed explicit-implicit strategy is adopted to process a matrix about a spatial derivative; reconstructing the time evolution matrix to obtain a two-step time domain discrete equation of the method; by introducing a single-step Leapfrog format, a unified numerical calculation format combining a CPML boundary condition and a single-step leapfrog CDHIE-FDTD method is derived. The technical problems that numerical calculation efficiency is low, time domain intermediate variables exist, and the correctness and effectiveness of integrated solution for the complex open domain electromagnetic problem containing a fine structure need to be improved are solved.

Description

Technical Field

[0001] The present invention relates to the field of electromagnetic problem data processing, and in particular to a method and system for solving open-domain electromagnetic targets of single-direction fine structures using a CD-HIE-FDTD numerical algorithm based on a CFS-PML technology. Background Art

[0002] Finite-difference time-domain (FDTD) is an effective computational tool for electromagnetic problems and is widely used in antenna design, optical device simulation, bioelectromagnetic research, electromagnetic compatibility analysis, and electromagnetic protection. For details, see the following references: DY Heh, and ELTan, "Divergence-Preserving Alternating Direction Implicit Scheme for Multi-Pole Debye Dispersive Media," IEEE Microw. Wirel. Compon. Lett., vol. 24, no. 2, pp. 69-71, 2014; G. Xie, N. Feng, M. Fang, M. Tong, WEI Sha, and Z. Huang, "AStability Extensible Symplectic FDTD (4, 4) Method and Its Application in Subgrid Technique for Arbitrary Grid Ratios," IEEE Trans. Antennas Propag., vol.71, no.11, pp.9124-9129, 2023; G.Xie, M.Fang, Z.Huang,X.Wu, Tech., vol.71, no.3, pp.1009-1018, 2023; S.Liu, B.Zou, L.Zhang, and S.Ren, "Heterogeneous CPU+GPU-Accelerated FDTD for Scattering ProblemsWith Dynamic Load Balancing," IEEE Trans.Antennas Propag., vol.68, no.9, pp.6734-6742, 2020; J.Chen, J.Li, and QHLiu, "Analyzing Graphene-Based Absorber by Using the WCS-FDTD Method,” IEEE Trans.Microw.Theory Tech., vol.65, no.10, pp.3689-3696, 2017. With the rapid development of technology, the research on electromagnetic problems is gradually oriented towards complex multi-scale structure modeling and simulation. For details, please refer to the following literature: J. Chen, and J. Wang, "A WCS-PSTD Method for Solving Electromagnetic Problems Both With Fine and Electrically Large Structures," IEEE Trans. Antennas Propag., vol. 62, no. 5, pp. 2695-2701, 2014; Q. Lei, J. Wang, J. Wang, and L. Shi, "An Efficient Unconditionally Stable FDTD Method for Analyzing Composite Right Left Handed Transmission Line," IEEE Antennas Wirel. Propag. Lett., pp. 1-1, 2017; Y. Wu, W. Liang, C. Ma, M. Wang, H. Zheng, and E. Li, "Thin Absorber With Dispersive Medium Simulated by Conformal Hybrid Explicit Implicit FDTD Method," IEEE Antennas Wirel.Propag.Lett.,vol.23,no.8,pp.2461-2465,2024. However, the CFL stability condition makes it difficult to give full play to the computational advantages of the FDTD method. In order to break through the CFL condition, some unconditionally stable fully implicit methods such as alternating-direction-implicit FDTD (ADI-FDTD), local one-dimensional FDTD (LOD-FDTD), Crank-Nicolson FDTD (CN-FDTD) methods have been proposed and widely used in solving complex electromagnetic problems. For details, see the following literature: S.Liu,and ELTan,“Leakage-Free Leapfrog ADI-FDTD Method for Lossy Media,”IEEE Antennas Wirel.Propag.Lett.,vol.23, no.2, pp.848-852, 2024; J.Feng, M.Fang, G.Xie, K.Song, W.Chen, Z.Huang, and X.Wu, "AnEfficient FDTD Method Based on Subgridding Technique and One-Step LeapfrogADI-FDTD," IEEE Microw.Wirel.Tech.Lett., vol.33, no.4, pp.375-378, 2023; G.Xie, G.Hou, K.Xu, K.Niu, N.Feng, Y.Li, Z.Huang, and AZElsherbeni, "Extension of FLOD-FDTD Method for Multiterm Modified Lorentz Model," IEEE Trans.AntennasPropag., vol.71, no.11, pp.8830-8840, 2023; H.Liu, Trans.Electromagn.Compat., vol.64, no.3, pp.827-839, 2022; Y.Ji, S.Wang, Y.Yu, H.Luan, Y.Wang, Q.Gao, and X.Zhao, "3-D Modeling and Analysis of Small-Loop Source TDEM Method Based on CFS-PML-CN-FDTD Method," IEEE Trans.Geosci.Remote Sensing,vol.62,pp.1-14,2024. These methods can use a larger time step to improve the efficiency of numerical calculations. However, the fully implicit method usually needs to solve 6 implicit equations. These operations are very time-consuming and are one of the key limiting factors affecting the efficiency of the fully implicit FDTD method. In addition, the selection of the time step in the fully implicit method is not without restrictions. A larger time step often causes greater damage to the numerical dispersion characteristics. When applying it, it is necessary to balance the relationship between calculation accuracy and calculation efficiency. .

[0003] In actual electromagnetic modeling and simulation, many device structures have fine structures in only one or two directions, such as electromagnetic coupling analysis of slender gaps, simulation of dielectric patch antennas, and shielding effectiveness analysis of thin shielding covers. For details, see the following references: J. Chen, and J. Wang, "A WCS-PSTD Method for Solving Electromagnetic Problems Both With Fine and Electrically Large Structures," IEEE Trans. Antennas Propag., vol. 62, no. 5, pp. 2695-2701, 2014; J. Chen, N. Xu, A. Zhang, and J. Guo, "Using Dispersion HIE-FDTD Method to Simulate the Graphene-Based Polarizer," IEEE Trans. Antennas Propag., vol. 64, no. 7, pp. 3011-3017, 2016; H. Mai, J. Chen, and A. Zhang, "A Hybrid Algorithm Based on FDTD and HIE-FDTD Methods for Simulating Shielding Enclosure," IEEE Trans.Electromagn.Compat., vol.60, no.5, pp.1393-1399, 2018; Y.-D. Kong, Trans.Antennas Propag., 2023; Q.Lei, J.Wang, andL.Shi, "Modified HIE-FDTD Algorithm for 2-D Periodic Structures With Fine Meshin Either Direction," IEEE Antennas Wirel.Propag.Lett., vol.16, pp.1455-1459, 2017.To this end, some weakly conditionally stable semi-implicit methods such as weakly conditionally stable FDTD (WCS-FDTD) and hybrid implicit-explicit (HIE-FDTD) methods have been proposed. In these methods, the selection of the time step is only restricted by the grid size in one or two directions, which weakens the requirement of the grid size on the time step in the traditional FDTD method. For details, see the following references: J. Wang, J. Wang, B. Zhou, and C. Gao, “An Efficient 3-D HIE-FDTD Method With Weaker Stability Condition,” IEEE Trans. Antennas Propag., vol. 64, no. 3, pp. 998-1004, 2016, W. Liang, Y. Wu, Y. Zhang, M. Wang, C. Fan, H. Zheng, and E. Li, “Nonuniform Grid Single-Field WCS FDTD Method Applied to Multithread Simulation of MIMO Dielectric Resonator Antenna,”IEEE Antennas Wirel. Propag. Lett.,vol.22,no.10,pp.2487-2491,2023, W. Liang,K. Zhang,Y. Wu,M. Wang,C. Fan,H. Zheng,and E. Li,“Implementation of ADE-CFS-PML for theSingle Field WCS-FDTD Method,”IEEE Microw. Wirel. Tech. Lett.,vol.33,no.7,pp.951-954,2023. Compared with the fully implicit FDTD method, the semi-implicit FDTD method only needs to solve fewer implicit equations and has higher computational efficiency. At the same time, the numerical dispersion error of the semi-implicit method is generally better than that of the fully implicit FDTD method. Therefore, the semi-implicit FDTD method shows stronger application potential when simulating electromagnetic models with fine structures in only one or two directions.In particular, during the development and numerical optimization of the HIE-FDTD method, a leapfrog HIE-FDTD method using the Peaceman-Rachford format was proposed. Compared with the traditional HIE-FDTD method, this method has a single-step leapfrog calculation characteristic similar to that of FDTD. At the same time, it further weakens the stability condition of the HIE-FDTD method, so that the selection of the time step is only limited by the grid size in one direction, thereby improving the efficiency of numerical calculation

[21] .

[0004] However, whether it is the single-step or double-step semi-implicit HIE-FDTD method or the fully implicit ADI-FDTD and LOD-FDTD methods, these methods do not have the divergence property, that is, they exhibit virtual charge phenomena in the passive region, which often introduces additional numerical calculation errors to the implicit method. Based on this, the literature ELTan, "ALeapfrog Scheme for Complying-Divergence Implicit Finite-Difference Time-Domain Method," IEEE Antennas Wirel. Propag. Lett., vol. 20, no. 5, pp. 853-857, 2021 proposed an effective CDI-FDTD method, which retains the unconditional stability characteristics of the fully implicit method while also having the divergence property.Subsequently, related research work includes: the construction of convolutional perfectly matched layer technology based on CDI-FDTD method

[23] , the application of CDI-FDTD method in dispersive medium simulation, for details, see: G. Xie, G. Hou, N. Feng, K. Song, M. Fang, Y. Li, X. Wu, and Z. Huang, “Auxiliary differential equation (ADE) method based complying-divergence implicit FDTD method for simulating the general dispersive anisotropic material,” Opt. Express, vol. 31, no. 11, 2023, G.-D. Xie, G.-L. Hou, K.-K. Niu, N.-X. Feng, M. Fang, Y.-S. Li, and Z.-X. Huang, “RI-CDI-FDTD method and program implementation for electromagnetic characteristics simulation of lossy Debyedispersive medium,” Acta Phys.Sin.,vol.72,no.15,2023 and CDI-FDTD method with lower numerical dispersion error, see: K.Zhang,P.Zhao,Y.-F.Cheng,S.Chen,N.Xu,G.Wang,andH.Zheng,“Optimized Traditional and Leapfrog CDI-FDTD Schemes With ReducedDispersion and Numerical Properties Analysis,”IEEE Trans.Microw.Theory Tech.,vol.71,no.11,pp.4725-4737,2023 were proposed one after another. These works have well verified the application value of unconditionally stable CDI-FDTD in electromagnetic simulation calculations.

[0005] In a recent report, Liu modified the curl matrix of the CDI-FDTD method and obtained a numerical algorithm that satisfies the divergence property and can unify the WCS and HIE calculation formats, achieving efficient and accurate solutions to electromagnetic problems. For details, see: S. Liu, and ELTan, "Unified Leapfrog CDI-FDTD Methods With WeaklyConditional Stability and Hybrid Implicit-Explicit Schemes," IEEE Trans. Antennas Propag., vol. 72, no. 8, pp. 6652-6662, 2024. However, to the best of our knowledge, no one has yet introduced PML technology into this effective single-step frog leap CDHIE-FDTD method for numerical simulation of open-domain electromagnetic problems. In addition, for this single-step frog leap CDHIE-FDTD method with divergence properties, the above work did not expand on the specific process, but only gave a more direct and conclusive description.Since Berenger proposed the perfect matching layer (PML), PML has been widely used in open-domain electromagnetic simulation calculations due to its excellent absorption performance. P.Wu, J.Wang, H.Jiang, Y.Xie, and T.Natsuki, "Fundamental Weaker HIE-PML Algorithm WithSimple Implementation and Enlarged Time Step for Fine Details Along a SingleDirection," IEEE Trans.Microw.Theory Tech., vol.72, no.1, pp.74-90, 2024, P.Wu, Y.Xie, H.Jiang, and L.Niu, "Performance-Enhanced Complex Envelope ADI-PML forBandpass EM Simulation," IEEE Microw.Wirel.Compon.Lett., vol.30, no.8, pp.729-732, 2020, N.Feng, J.Wang, J.Zhu, Z.Liang, GPWang, and WTJoines, "SwitchableTruncations Between the 1st-and 2nd-Order DZT-CFS-UPMLs for Relevant FDTDProblems," IEEE Trans.Antennas Propag., vol.68, no.1, pp.360-365, 2020, HLJiang, JFZhang, and TJCui, "An Alternative PML Implementation for ArbitraryMedia," IEEE Trans. Antennas Propag., vol.67, no.1, pp.612-615, 2019.In addition, researchers have proposed several advanced PML formulas, including anisotropic PML, recursive integral PML, complex frequency shift (CFS) PML, etc. For details, please refer to the following literature: A. Giannopoulos, "Multipole Perfectly Matched Layer for Finite-Difference Time-Domain Electromagnetic Modeling," IEEE Trans. Antennas Propag., vol. 66, no. 6, pp. 2987-2995, 2018, J. Liu, X. Xiao, J. Tang, Z. Ren, X. Huang, and J. Zhang, "Fast 3-D Controlled-Source Electromagnetic Modeling Combining UPML and Rational Krylov Method," IEEE Geosci. Remote Sens. Lett., vol. 19, pp. 1-5, 2022, N. Feng, Y. Zhang, Q. Sun, J. Zhu, WT Joines, and QH Liu, "An Accurate 3-D CFS-PMLBased Crank–Nicolson FDTD Method and Its Applications in Low-FrequencySubsurface Sensing,” IEEE Trans.Antennas Propag., vol.66, no.6, pp.2967-2975, 2018.Among them, convolutional PML (CPML) is an effective implementation of CFS-PML, which has been proven to minimize numerical reflections, especially in absorbing evanescent waves. For details, see the following references: M.-L.Zhai, H.-L.Peng, Z.Chen, and W.-Y.Yin, “The Implementation of Convolutional Perfectly Matched Layer (CPML) for One-Step Leapfrog HIE-FDTD Method,” IEEE Antennas Wirel.Propag.Lett., vol.14, pp.1694-1697, 2015, AK Pandey, and AK Saxena, “Implementation of Higher-Order CPML on the HIE-FDTD Methods,” IEEE Trans.Electromagn.Compat., vol.65, no.6, pp.1839-1847, 2023.

[0006] In summary, the existing technology has technical problems such as low numerical calculation efficiency, the existence of intermediate variables in the time domain, and the correctness and effectiveness of the integrated solution to complex open-domain electromagnetic problems containing fine structures need to be improved. Summary of the invention

[0007] The technical problem to be solved by the present invention is: how to solve the technical problems in the prior art of low numerical calculation efficiency, the existence of intermediate variables in the time domain, and the need to improve the correctness and effectiveness of the integrated solution to complex open-domain electromagnetic problems containing fine structures.

[0008] The present invention adopts the following technical solution to solve the above technical problem: an open domain electromagnetic target solution method for unidirectional fine structure includes:

[0009] S1. The PML stretch factor is introduced into Maxwell's equations to construct the time evolution matrix of the electromagnetic field components.

[0010] S2. Use a hybrid explicit and implicit difference strategy to process the matrix of spatial derivatives in the time evolution matrix; reconstruct the time evolution matrix to obtain a two-step time domain discrete equation;

[0011] S3. Based on the two-step time domain discrete equation and the single-step Leapfrog format, a unified numerical calculation format is derived to combine the CPML boundary conditions with the single-step frog leap CDHIE-FDTD method, and the single-step Leapfrog format is applied to the preset sub-time step logic to perform CD-HIE-FDTD numerical solution.

[0012] Specifically, the preset sub-time step logic includes the following equations (16) to (19).

[0013] The present invention introduces CFS into the Maxwell equations and constructs a time evolution matrix, then performs time-space operator splitting on the matrix to construct a two-step iterative solution with divergence properties. In order to improve the efficiency of numerical calculations, the present invention adopts the time stepping strategy of the Leapfrog format to eliminate the intermediate variables in the time domain and obtain a single-step differential discretization solution. The present invention uses recursive convolution technology to process the convolution operator in the discrete formula to obtain the iterative formula of the electromagnetic field components and the auxiliary variables in the PML region. In addition, for the above-mentioned calculation framework, the analysis of stability, divergence properties and complexity is completed, which can effectively realize the simulation of solving complex open-domain electromagnetic problems.

[0014] In a more specific technical solution, S1 includes:

[0015] S11. Use the following logic to express the curl equation of the electromagnetic field:

[0016]

[0017] In the formula, S w (w=x,y,z) is the scaling factor of PML:

[0018]

[0019] S12. Convert equation (1), equation (2) and from frequency domain to time domain and write them in matrix form:

[0020]

[0021] Among them, u=[E x ,E y ,E z ,H x ,H y ,H z ] T , the spatial curl matrix

[0022] S13. Assuming that there is no current source and magnetic current source, use the following logic to obtain the exact solution of equation (4) from time t = 0 to time t = Δt:

[0023] u(Δt)=exp(RΔt)u(0) (7)

[0024] R in equation (7) can be split into two sub-matrices to obtain the time evolution matrix.

[0025] In a more specific technical solution, S2 includes:

[0026] S21. Since the traditional FDTD method uses explicit differences, the time step Δt needs to meet the following conditions:

[0027]

[0028] Where Δx, Δy, Δz are the spatial steps in the x, y and z directions respectively. This means that the simulation time step value is affected by the selection of the spatial step. Based on the alternating direction implicit iteration method proposed by Peaceman and Rachford, the present invention uses the difference strategy to solve the Maxwell curl equations. In addition, for electromagnetic targets with fine structures in a single direction, the alternating direction implicit iteration scheme in three directions is adopted, which not only requires solving more complex matrix equations, but also has a larger numerical dispersion error under the same simulation time step. Therefore, the present invention only performs implicit iteration in one direction, and uses an explicit difference strategy in other directions.

[0029] The sub-matrices R1 and R2 of the time evolution matrix are processed using a hybrid implicit and explicit difference strategy in different directions. Implicit differences are used in the y direction, and explicit differences are used in the x and z directions. The sub-matrices R1 and R2 are expressed using the following logic:

[0030]

[0031] in:

[0032]

[0033] In addition, the leapfrog HIE-FDTD method proposed in the literature (see the literature: J. Wang, B. Zhou, L. Shi, C. Gao, and B. Chen, "A Novel 3-D HIE-FDTD Method With One-Step Leapfrog Scheme," IEEE Trans. Microw. Theory Tech., vol. 62, no. 6, pp. 1275-1283, 2014) does not obey the Gaussian divergence theorem and introduces non-physical virtual charges in numerical calculations. Therefore, this scheme further studies the exact solution of the Maxwell equations in the time domain in (7).

[0034] S22. Express the exponential matrix exp(RΔt) in equation (7) as the following matrix multiplication form:

[0035]

[0036] The spatial curl operator R in the time evolution matrix in equation (11) can be divided into two sub-matrices, R1 and R2, by the hybrid explicit and implicit strategy in equation (8). At the same time, according to the arrangement of the sub-matrices in the above equation, there are 4! forms of matrix multiplication. Some of them cannot satisfy the dual relationship between the electric field and the magnetic field in the Maxwell equations after discretization. In addition, some do not satisfy the Gaussian divergence condition. After testing, the scheme that satisfies the divergence condition is given below, and it is verified below that the proposed scheme satisfies the divergence characteristics. Finally, the discrete equation that can accurately describe the Maxwell curl equation is obtained and used to solve the electromagnetic target.

[0037] S23. Reconstruct the time evolution matrix in formula (11) to obtain the following exponential matrix multiplication scheme, which satisfies the divergence property:

[0038]

[0039] Through Taylor series expansion, equation (13) can be expressed as:

[0040]

[0041] According to equations (8) and (14), the time-space discrete form of the electromagnetic field components is obtained:

[0042]

[0043] S24. Use the following logic to divide equation (15) into two parts:

[0044] as well as

[0045] S25. Define auxiliary variable Q = [Q ex ,Q ey ,Q ez ,Q hx ,Q hy ,Q hz ] T , dealing with the inversion problem of the operator matrix in two sub-time steps, expressing the two-step time domain discrete equation.

[0046] In a more specific technical solution, in S25, the first sub-time step is derived using the following logic:

[0047]

[0048] Using the following logic, the second sub-time step is derived:

[0049]

[0050]

[0051] According to the first sub-time step and the second sub-time step of equations (16) to (19), the two-step time domain discretization equation of the CDHIE-FDTD method is described.

[0052] In more specific technical solutions, S3 includes:

[0053] In equations (16)-(19), we need to solve the field values ​​of the electromagnetic field components E and H at time n+1 / 2 and n+1. This means that when we perform electromagnetic simulation, we must calculate twice to get the updated electromagnetic field value for a time period. This scheme is based on the leapfrog update scheme, which processes (16)-(19), eliminates the intermediate field values ​​of the electromagnetic and magnetic fields, and obtains a more efficient CD-HIE-FDTD method.

[0054] S31, updating the electric field component in the time period from n-1 / 2 to n+1 / 2;

[0055] Now, taking equation (19) back one time step and substituting it into equation (16), we obtain

[0056]

[0057] Substituting equation (17) into equations (18) and (20), we obtain

[0058]

[0059] Next, we move the field component in equation (19) back one time step and substitute it into equation (22), and we get

[0060]

[0061] According to equation (23), the updated equation of the electromagnetic field component in the time period from n to n+1 / 2 can be obtained:

[0062]

[0063] The iterative formula for the electric field components from n to n+1 / 2 can be obtained by equation (24):

[0064]

[0065] From equations (9) and (10), we can see that the matrices C and D in equation (25) contain convolution terms. Then we define the auxiliary variable F w ( or ), according to equation (3), we can get

[0066]

[0067] in

[0068] Next, the convolution term in equation (25) is processed by the convolution technique, and the iterative formula of the field component in the PML region can be obtained:

[0069]

[0070] Due to the different spatial distributions of the electric and magnetic fields in the Yee grid, the spatial distributions of the relevant PML parameters in the PML region are also different. Therefore, C is set e ,D e and C h ,D h , they represent the operation operators in the electric field region and the magnetic field region respectively.

[0071]

[0072]

[0073] ψ E =[ψ Exy ,ψ Exz ,ψ Eyz ,ψ Eyx ,ψ Ezx ,ψ Ezy ] T , ψ H =[ψ Hxy ,ψ Hxz ,ψ Hyz ,ψ Hyx ,ψ Hzx ,ψ Hzy ] T ,

[0074]

[0075] G e =diag(g ex ,g ex ,g ey ,g ey ,g ez ,g ez ),G h =diag(g hx ,g hx ,g hy ,g hy ,g hz ,g hz ),

[0076]

[0077] Where g and f are the coefficients of the updated recursive convolution term, as follows

[0078] In order to obtain the electric field component at time n–1 / 2, we substitute equation (17) into equation (21) and obtain

[0079]

[0080] Next, substitute equation (19) into equation (30) and go back one time step to obtain

[0081]

[0082] Similar to the processing method of equation (23), the iterative formula of the electric field component at time n–1 / 2 in the PML region is further obtained through equation (31):

[0083]

[0084] Finally, by subtracting equation (32) from equation (27), we obtain the iterative formula for the electric field component from time n–1 / 2 to n+1 / 2:

[0085]

[0086] S32. Update the magnetic field components in the time period from n to n+1.

[0087] In step S32, equation (21) is substituted into equation (19) to obtain:

[0088]

[0089] Substituting equation (17) into equation (34), we obtain:

[0090]

[0091] The same processing method as equation (23) is used. By further deducing equation (35), the time step equation of the magnetic field component at time n+1 is obtained:

[0092]

[0093] in,

[0094] In this embodiment, in order to obtain the magnetic field component at time n, equation (19) is now retreated by one sub-time step and substituted into equation (20) to obtain:

[0095]

[0096] Substituting equation (17) into equation (37), we obtain:

[0097]

[0098] The time step equation of the magnetic field at time n is further obtained by the above equation (38):

[0099]

[0100] In this embodiment, equation (36) is subtracted from equation (39) to obtain:

[0101]

[0102] Equation (40) is a time-stepping equation for updating the magnetic field components based on a single-step frog leapfrog scheme.

[0103] S33. Use the following logic to define auxiliary variables:

[0104] and

[0105] S34. Based on the CDHIE-FDTD method, the electromagnetic field components in the CPML region are iterated:

[0106]

[0107]

[0108] In the formula, ψ e =[ψ exy ,ψ exz ,ψ eyz ,ψ eyx ,ψ ezx ,ψ ezy ] T , ψ h =[ψ hxy ,ψ hxz ,ψ hyz ,ψ hyx ,ψ hzx ,ψ hzy ] T .

[0109] In a more specific technical solution, S31 includes:

[0110] S311. Substitute equation (19) in the second sub-time step back one time step and substitute it into equation (16) in the first sub-time step to obtain:

[0111]

[0112] Substituting equation (17) of the first time step into equations (18) and (20) of the second sub-time step, we obtain:

[0113]

[0114] Substituting the field component in equation (19) back one time step into equation (22), we get:

[0115]

[0116] S312. According to formula (23), the updated equation of the electromagnetic field component in the time period from n to n+1 / 2 can be obtained:

[0117]

[0118] According to formula (24), the iterative formula of the electric field component from n to n+1 / 2 is obtained:

[0119]

[0120] S313, define auxiliary variable F w ;

[0121] S314. Write the convolution term in equation (25) as a discrete circular convolution form to obtain the iterative formula of the field component in the CPML region:

[0122]

[0123] S315. Now substitute equation (17) into equation (21) to obtain the electric field component at time n–1 / 2:

[0124]

[0125] Substituting equation (19) into equation (30) and going back one time step, we get:

[0126]

[0127] Processing equation (31), we obtain the iterative formula of the electric field component at time n–1 / 2 in the CPML region:

[0128]

[0129] S316. Subtract equation (32) from equation (27) to obtain the iterative formula of the electric field component from time period n–1 / 2 to n+1 / 2, which is used as the time step equation for updating the electric field component based on the single-step frog leapfrog format:

[0130]

[0131] In a more specific technical solution, in S313, or w=x,y,z. According to formula (3), we get:

[0132]

[0133] In the formula,

[0134] In a more specific technical solution, in S314, C is used in formula (27) respectively. e ,D e and C h ,D h , represents the spatial operation operator of the electric field region and the magnetic field region;

[0135] The following logic is used to iterate the auxiliary variables in equation (27):

[0136]

[0137] In a more specific technical solution, an open-domain electromagnetic target solving system for unidirectional fine structures includes:

[0138] The time evolution matrix construction module is used to introduce the PML's stretching factor into Maxwell's equations to construct the time evolution matrix of the electromagnetic field components;

[0139] A two-step time domain discrete equation acquisition module is used to process the matrix of spatial derivatives in the time evolution matrix by adopting a hybrid explicit and implicit difference strategy; the time evolution matrix is ​​reconstructed to obtain a two-step time domain discrete equation. The two-step time domain discrete equation acquisition module is connected with the time evolution matrix construction module;

[0140] The CPML boundary condition and single-step frog leap combined numerical solution module is used to derive a unified numerical calculation format based on the two-step time domain discrete equation and the single-step Leapfrog format, so as to combine the CPML boundary condition and the single-step frog leap CDHIE-FDTD method, apply the single-step Leapfrog format to equations (16) to (19) to perform CD-HIE-FDTD numerical solution, and the CPML boundary condition and single-step frog leap combined numerical solution module is connected with the two-step time domain discrete equation acquisition module.

[0141] Compared with the prior art, the present invention has the following advantages:

[0142] The present invention introduces CFS into the Maxwell equations and constructs a time evolution matrix, then performs time-space operator splitting on the matrix to construct a two-step iterative solution with divergence properties. In order to improve the efficiency of numerical calculations, the present invention adopts the time stepping strategy of the Leapfrog format to eliminate the intermediate variables in the time domain and obtain a single-step differential discretization solution. The present invention uses recursive convolution technology to process the convolution operator in the discrete formula to obtain the iterative formula of the electromagnetic field components and the auxiliary variables in the PML region. In addition, for the above-mentioned calculation framework, the analysis of stability, divergence properties and complexity is completed, which can effectively realize the simulation of solving complex open-domain electromagnetic problems.

[0143] The present invention solves the technical problems existing in the prior art of low numerical calculation efficiency, existence of intermediate variables in the time domain, and the need to improve the correctness and effectiveness of the integrated solution to complex open-domain electromagnetic problems containing fine structures. BRIEF DESCRIPTION OF THE DRAWINGS

[0144] Figure 1 A schematic diagram of basic steps of an open-domain electromagnetic target solving method for a unidirectional fine structure according to Embodiment 1 of the present invention;

[0145] Figure 2a Schematic diagram of the analysis of the divergence value (dB) of electric flux density obtained by two methods in Example 2 of the present invention;

[0146] Figure 2b This is a schematic diagram of the divergence analysis of the calculation formula of the FDTD method of Example 2 of the present invention;

[0147] Figure 2c Schematic diagram of divergence analysis of the CD-HIE-FDTD method calculation formula (51) in Example 2 of the present invention;

[0148] Figure 2d Schematic diagram of divergence analysis of calculation formulas (49)-(51) of the CD-HIE-FDTD method in Example 2 of the present invention;

[0149] Figure 3 This is a schematic diagram of the calculation area division of the CDHIE-FDTD method according to Example 2 of the present invention;

[0150] Figure 4a Schematic diagram of the distribution of the modulus value of the growth factor ζ in the first layer CPML plane region in Example 2 of the present invention;

[0151] Figure 4b Schematic diagram of the distribution of the modulus value of the growth factor ζ in the plane region of the fourth layer CPML of Example 2 of the present invention;

[0152] Figure 4c Schematic diagram of the distribution of the modulus value of the growth factor ζ in the seventh CPML plane region of Example 2 of the present invention;

[0153] Figure 4d Schematic diagram of the distribution of the modulus value of the growth factor ζ in the plane region of the tenth layer of CPML according to Example 2 of the present invention;

[0154] Figure 5a Schematic diagram of the distribution of the modulus value of the growth factor ζ in the edge region of the first layer of CPML according to Embodiment 2 of the present invention;

[0155] Figure 5b Schematic diagram of the distribution of the modulus value of the growth factor ζ in the edge region of the fourth layer CPML of Example 2 of the present invention;

[0156] Figure 5c Schematic diagram of the distribution of the modulus value of the growth factor ζ in the edge region of the seventh layer of CPML according to Example 2 of the present invention;

[0157] Figure 5d Schematic diagram of the distribution of the modulus value of the growth factor ζ in the edge region of the tenth layer of CPML in Example 2 of the present invention;

[0158] Figure 6a This is a schematic diagram of the distribution of the modulus value of the growth factor ζ in the corner vertex area of ​​the first layer of CPML according to Embodiment 2 of the present invention;

[0159] Figure 6b Schematic diagram of the distribution of the modulus value of the growth factor ζ in the corner vertex area of ​​the fourth layer CPML of Example 2 of the present invention;

[0160] Figure 6c Schematic diagram of the distribution of the modulus value of the growth factor ζ in the corner vertex area of ​​the seventh layer CPML of Example 2 of the present invention;

[0161] Figure 6d Schematic diagram of the distribution of the modulus value of the growth factor ζ in the corner vertex area of ​​the tenth layer of CPML in Example 2 of the present invention;

[0162] Figure 7 A 3D geometric diagram of placing a metal sheet in the xz plane according to Embodiment 3 of the present invention;

[0163] Figure 8 This is a schematic diagram of the relative reflection error of CDHIE-FDTD under different CFLNs under the 10-layer CPML boundary condition of Example 3 of the present invention;

[0164] Fig. 9 This is a structural diagram of a FSS unit loaded with a medium of 0.25 mm thickness according to Example 3 of the present invention;

[0165] Fig.10 It is a schematic diagram comparing transmission coefficients calculated by the traditional FDTD method, the CD-HIE-FDTD method and the commercial software CST according to Example 3 of the present invention;

[0166] Fig.11 Schematic diagram of shielding effectiveness of Example 3 of the present invention;

[0167] Fig.12 A comparison diagram of the calculation results of the traditional FDTD method and the CD-HIE-FDTD method in Example 3 of the present invention;

[0168] Fig.13 Schematic diagram of the geometric structure of different slots in Example 3 of the present invention;

[0169] Fig.14 A comparison diagram of shielding effectiveness of different slots in Example 3 of the present invention;

[0170] Fig.15 This is a diagram of an inverted-F antenna model according to Embodiment 3 of the present invention;

[0171] Fig.16 A comparison diagram of sampled voltages calculated at the port by two methods according to Embodiment 3 of the present invention;

[0172] Fig.17 A comparison diagram of sampled currents calculated at the port by two methods according to Embodiment 3 of the present invention;

[0173] Fig.18 4 is a comparison diagram of S11 parameters of the inverted-F antenna according to Embodiment 3 of the present invention;

[0174] Fig.19 2.4 GHz, antenna pattern in the xy plane of Embodiment 3 of the present invention;

[0175] Fig. 20 This is the antenna radiation pattern of 5.8 GHz, xy plane of Example 3 of the present invention. DETAILED DESCRIPTION

[0176] In order to make the purpose, technical solution and advantages of the embodiments of the present invention clearer, the technical solution in the embodiments of the present invention will be clearly and completely described below in combination with the embodiments of the present invention. Obviously, the described embodiments are part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.

[0177] Example 1

[0178] like Figure 1 As shown, the present invention provides an open-domain electromagnetic target solution method for a unidirectional fine structure, comprising the following basic steps:

[0179] S1. The PML stretch factor is introduced into Maxwell's equations to construct the time evolution matrix of the electromagnetic field components.

[0180] In the linear, lossless space of this embodiment, the PML truncation calculation domain in the form of complex frequency domain shift tensor coefficient coordinate scaling is adopted, and the curl equation of the electromagnetic field can be expressed as follows:

[0181]

[0182] Among them, S w (w=x,y,z) is the scaling factor of PML:

[0183]

[0184] Convert the above formula (1) and formula (2) from frequency domain to time domain and write them in matrix form:

[0185]

[0186] Among them, u=[E x ,E y ,E z ,H x ,H y ,H z ] T , the spatial curl matrix

[0187] In this embodiment, assuming that there is no current source and magnetic current source, the exact solution of the above equation (4) from time t=0 to time t=Δt can be expressed as follows:

[0188] u(Δt)=exp(RΔt)u(0) (7)

[0189] The R in the above formula (7) can be split into two sub-matrices, namely, R=R1+R2. Since R can be combined with matrices in various forms, the matrix expressions of R1 and R2 are also diverse.

[0190] S2. A hybrid explicit and implicit strategy is used to process the matrix of spatial derivatives. By reconstructing the time evolution matrix, the double-step time domain discrete equation of this method is obtained.

[0191] In this embodiment, for the hybrid implicit and explicit difference strategies in different directions, the construction forms of matrices R1 and R2 are similar. If implicit difference is used in the y direction and explicit difference is used in the x and z directions, the matrix forms of R1 and R2 are as follows:

[0192]

[0193] in:

[0194]

[0195] In this embodiment, the exponential matrix exp(RΔt) in equation (7) can be expressed in the form of matrix multiplication as follows:

[0196]

[0197] In the above equation (11), exp(R1Δt / 2), exp(R2Δt / 2), exp(-R1Δt / 2) -1 and exp(-R2Δt / 2) -1 There are 24 (4×3×2×1) combinations of matrix multiplication, and different combinations often correspond to different numerical discrete forms of Maxwell equations. Among them, the time evolution matrix constructed by the Peaceman-Rachford scheme used in the traditional frog leaping format HIE-FDTD method is:

[0198]

[0199] However, the numerical discrete formula obtained by expanding the Maxwell equations according to formula (12) does not satisfy the divergence property.

[0200] Based on this, in this embodiment, the time evolution matrix in equation (11) is reconstructed to obtain the following exponential matrix multiplication scheme that satisfies the divergence property:

[0201]

[0202] In this embodiment, through Taylor series expansion, equation (13) can be expressed as follows:

[0203]

[0204] According to formula (8) and formula (14), the time-space discrete form of the electromagnetic field component can be obtained as follows:

[0205]

[0206] Equation (15) can be divided into two parts for calculation, namely:

[0207] as well as

[0208] In addition, in order to deal with the problem of inverting the operator matrix in two sub-time steps, an auxiliary variable Q = [Q ex ,Q ey ,Q ez ,Q hx ,Q hy ,Q hz ]T , the specific implementation is as follows:

[0209] In this embodiment, the first sub-time step is derived using the following formula:

[0210]

[0211] Using the following formula, the second sub-time step is derived:

[0212]

[0213] In this embodiment, the aforementioned equations (16)-(19) describe the double-step time domain discretization CDHIE-FDTD method.

[0214] S3. By introducing the single-step Leapfrog format, a unified numerical calculation format is derived and the single-step Leapfrog format is applied to the above equations (16) to (19).

[0215] In this embodiment, the aforementioned unified numerical calculation format can combine the CPML boundary condition and the single-step frog leap CDHIE-FDTD method.

[0216] In this embodiment, the electric field component is updated in the time period from n-1 / 2 to n+1 / 2;

[0217] In this embodiment, equation (19) is backed up by one time step and substituted into equation (16), resulting in:

[0218]

[0219] Substituting equation (17) into equations (18) and (20), we obtain:

[0220]

[0221] Substituting the field component in equation (19) back one time step and into equation (22) yields:

[0222]

[0223] According to equation (23), the updated equation of the electromagnetic field component in the time period from n to n+1 / 2 can be obtained:

[0224]

[0225] The iterative formula for the electric field components from n to n+1 / 2 can be obtained by equation (24):

[0226]

[0227] In this embodiment, it can be seen from equations (9) and (10) that the matrices C and D in equation (25) above contain convolution operators, so the auxiliary variable F is defined as w ( or w=x,y,z), according to equation (3), we can get:

[0228]

[0229] in,

[0230] Writing the convolution term in equation (25) as a discrete circular convolution form, we can obtain the iterative formula for the field component in the CPML region:

[0231]

[0232] Specifically, according to the basic electromagnetic theorem, in the Yee grid discretization scheme, the time-space configuration of the electric field and the magnetic field is different, resulting in different configurations of the corresponding PML parameters in the PML region. Therefore, C is used in formula (27) to represent e ,D e and C h ,D h , represents the spatial operation operator of the electric field region and the magnetic field region, and its specific form is as follows:

[0233]

[0234] At the same time, the iterative formula of the auxiliary variables in formula (27) is as follows:

[0235]

[0236] in:

[0237] ψ E =[ψ Exy ,ψ Exz ,ψ Eyz ,ψ Eyx ,ψ Ezx ,ψ Ezy ] T , ψ H =[ψ Hxy ,ψ Hxz ,ψ Hyz ,ψ Hyx ,ψ Hzx ,ψ Hzy ] T ,

[0238]

[0239] G e=diag(g ex ,g ex ,g ey ,g ey ,g ez ,g ez ),G h =diag(g hx ,g hx ,g hy ,g hy ,g hz ,g hz ),

[0240]

[0241] In this embodiment, for G e and F e The matrix element g in ew and f ew (w=x,y,z), whose forms are: G h and F h The matrix element g in hw and f hw The algebraic form of is similar.

[0242] In this embodiment, in order to obtain the electric field component at time n–1 / 2, equation (17) is substituted into equation (21) to obtain:

[0243]

[0244] Substituting equation (19) into equation (30) and backing off one time step, we obtain:

[0245]

[0246] In this embodiment, equation (31) is processed using a strategy similar to that used to process equation (23), and an iterative formula for the electric field component at time n–1 / 2 in the CPML region can be obtained:

[0247]

[0248] In this embodiment, the above equation (32) is subtracted from the above equation (27) to obtain the iterative formula of the electric field component in the time period from n–1 / 2 to n+1 / 2:

[0249]

[0250] Equation (33) is a time-stepping equation for updating the electric field components based on a single-step frog leapfrog scheme.

[0251] In this embodiment, the magnetic field component of the time period from n to n+1 is updated. Substituting equation (21) into equation (19), we obtain:

[0252]

[0253] Substituting equation (17) into equation (34), we obtain:

[0254]

[0255] The same processing method as equation (23) is used. By further deducing equation (35), the time step equation of the magnetic field component at time n+1 is obtained:

[0256]

[0257] in,

[0258] In this embodiment, in order to obtain the magnetic field component at time n, equation (19) is now retreated by one sub-time step and substituted into equation (20) to obtain:

[0259]

[0260] Substituting equation (17) into equation (37), we obtain:

[0261]

[0262] The time step equation of the magnetic field at time n is further obtained by the above equation (38):

[0263]

[0264] In this embodiment, equation (36) is subtracted from equation (39) to obtain:

[0265]

[0266] Equation (40) is a time-stepping equation for updating the magnetic field components based on a single-step leapfrog scheme.

[0267] In this embodiment, although equation (33) and equation (40) constitute the basic equations for solving the electromagnetic field components by the CDHIE-FDTD method, the matrix inversion problem in the equation is still very complex and time-consuming. Therefore, the following auxiliary variables are defined in this embodiment: and Finally, the iterative formula of the electromagnetic field components in the CPML region based on the CDHIE-FDTD method is as follows:

[0268]

[0269] Among them, ψ e =[ψ exy ,ψ exz ,ψ eyz ,ψ eyx ,ψ ezx ,ψ ezy ] T , ψ h =[ψ hxy ,ψ hxz ,ψ hyz ,ψ hyx ,ψ hzx ,ψ hzy ] T .

[0270] Example 2

[0271] In this embodiment, a numerical characteristic analysis of stability, divergence property and complexity is performed on an open-domain electromagnetic target solving method for a unidirectional fine structure provided by the present invention.

[0272] In the divergence analysis performed in this embodiment, the divergence condition of Gauss's theorem is a key feature of the Maxwell curl equation. When a current source is introduced into the Maxwell equations, it is crucial to maintain the divergence condition during the process of discretizing the equations through a numerical difference scheme. The divergence of the electric flux density of the Maxwell equations including the scaling factor can be expressed as:

[0273]

[0274] The above equation (47) shows that in the active region, the divergence of the electric displacement vector is equal to the charge density (if in the passive region, the divergence of the electric flux vector should be zero).

[0275] In this embodiment, the divergence property of the leapfrog CDHIE-FDTD method is numerically verified based on the open domain electromagnetic simulation space. The simulation space mainly includes two regions, the free space region and the CPML region.

[0276] In the CPML region, the divergence of the electric flux density can be expressed as follows:

[0277]

[0278] Substituting equation (26) into equation (48) and further processing using the recursive convolution method

[35] , the divergence of the electric flux density in the CPML region is:

[0279]

[0280] Among them, the numerical iteration formula of the auxiliary variable is:

[0281]

[0282] In free space, the expansion factor is 1, so formula (49) can be simplified to the following form:

[0283]

[0284] like Figure 2a to Figure 2d As shown, the divergence of the electric flux density solved by the traditional FDTD method and the CD-HIE-FDTD method is shown. Figure 2a and Figure 2c The calculation results are obtained by solving formula (51) using the traditional FDTD method and CD-HIE-FDTD. It can be found that in the CPML region, the divergence value of the electric flux density is greater than the divergence value in the free space region. This is because the constitutive relationship between the electric flux density and the electric field in the CPML region is different from that in the free space. Therefore, when solving the divergence of the electric flux density in the CPML region, it is necessary to use formulas (49) and (50). Figure 2b and Figure 2d The results are obtained by using formula (51) in free space and formulas (49) and (50) in the CPML region. It can be found that the traditional FDTD method and the CD-HIE-FDTD method can obtain similar calculation results to those in free space in the CPML region. In addition, according to Figure 2a to Figure 2d It can be seen that the results calculated by the traditional FDTD method and the CD-HIE-FDTD method have an obvious peak at the location where the excitation source exists, which is also consistent with the electromagnetic theory described by formula (47). This shows that when the CDHIE-FDTD method is used to discretize Maxwell's equations, it can have the same good divergence characteristics as the traditional FDTD method.

[0285] In the algorithm complexity analysis performed in this embodiment, the traditional HIE-FDTD method, see the literature: J. Chen, and J. Wang, "A three-dimensional semi-implicit FDTD scheme for calculation of shielding effectiveness of enclosure with thin slots," IEEE Trans. Electromagn. Compat., vol. 49, no. 2, pp. 354-360, 2007, the leapfrog HIE-FDTD method (see the literature: J. Wang, B. Zhou, L. Shi, C. Gao, and B. Chen, "A Novel 3-D HIE-FDTD Method With One-Step Leapfrog Scheme," IEEE Trans. Microw. Theory Tech., vol. 62, no. 6, pp. 1275-1283, 2014) and the numerical complexity of the leapfrog CDHIE-FDTD method are analyzed. During the numerical iteration, two field components in all three methods involve solving tridiagonal matrix equations (for the implicit case in the y direction, the leapfrog HIE-FDTD method is E x ,H x , the leapfrog CDHIE-FDTD method is e x ,h x , the traditional HIE-FDTD method is E x and E z). Table I lists the computational complexity comparison of the three methods. The results show that the CDHIE-FDTD method has fewer computational operators. This is because the right side of the equation of the CDHIE-FDTD method has a concise form of the Fundamental scheme. See the literature: ELTan, “From Time-Collocated to Leapfrog Fundamental Schemes for ADI and CDI FDTD Methods,” Axioms, vol. 11, no. 1, 2022. In addition, when the CDHIE-FDTD method introduced in equations (41)-(46) implements CPML, 36 A / S and 30 M / D calculation operators need to be added, which is consistent with the calculation operators required to be added when the frog leapfrog HIE-FDTD method introduced in the document M.-L.Zhai, H.-L.Peng, Z.Chen, and W.-Y.Yin, "The Implementation of Convolutional Perfectly Matched Layer (CPML) for One-Step Leapfrog HIE-FDTD Method," IEEE Antennas Wirel. Propag. Lett., vol.14, pp.1694-1697, 2015 implements CPML. Therefore, the CPML numerical implementation scheme based on the CDHIE-FDTD method proposed in the present invention has high computational efficiency when solving open domain electromagnetic problems.

[0286] Table I Computational complexity of different methods

[0287] Method A / S M / D Conventional-HIE-FDTD 42 16 LF-HIE-FDTD 52 18 LF-CDHIE-FDTD 40 16

[0288] like Figure 3 As shown, in the stability analysis performed in this embodiment, the von Neumann method is used to test the stability of the proposed CD-HIE-FDTD algorithm combined with the CPML technology. As shown in FIG2 , we divide the calculation area by selecting different values ​​for CFS.

[0289] (1) Simulation area: Set the scaling factor (S x =1, S y =1, S z =1) into equations (41)-(46), and through further simplification, the following matrix form is finally obtained:

[0290] M L u n+1 =M R u n (52)

[0291] in:

[0292]

[0293] D w represents the first-order spatial derivative operator. Since the first-order spatial derivative is approximated by central difference, D w =2jsin(k w Δw / 2) / Δw, k w Represents the wave number in the w direction.

[0294] In the case of a time-harmonic field, the solution of a three-dimensional plane wave at the nth moment is

[0295] u(x,y,z,n)=u0ζ n exp[j(k x x+k y y+k z z)] (53)

[0296] Among them, u0 represents the amplitude of the field component and ζ represents the growth factor.

[0297] In this embodiment, substituting equation (53) into equation (52) yields

[0298] (M L ζ-M R ) n =Mu n =0 (54)

[0299] In order for equation (54) to have a non-zero solution, the determinant of its coefficient matrix must be 0, that is, |M| = 0. Through algebraic processing, we can obtain

[0300]

[0301] By solving equation (55), we obtain ζ 1,2 =1, in

[0302] ζ 1,2 The magnitude is 1, for the growth factor ζ 3,4 And 5,6 , only if N 2 ≤M 2 When , the amplitude can be kept at 1, that is

[0303]

[0304] Since M>0, and Just need to satisfy as well as This ensures that the magnitude of the growth factor is 1. Further deduction yields Δt≤Δx / c0 and Δt≤Δz / c0. Therefore, the CFL stability condition of this method is:

[0305] Δt≤min[Δx / c0,Δz / c0] (57)

[0306] (2) CPML region: Equations (41)-(46) are described in the form of a matrix:

[0307]

[0308] where u PML =[E x ,E y ,E z H x ,H y ,H z ,ψ exy ,ψ exz ,ψ eyz ,ψ eyx ,ψ ezx ,ψ ezy ,ψ hxy ,ψ hxz ,ψ hyz ,ψ hyx ,ψ hzx ,ψ hzy ] T , substituting equation (53) into (58), we can obtain:

[0309]

[0310] Coefficient matrix M L and M R for:

[0311]

[0312]

[0313] P w =1-bdD ew D hw .

[0314] In this embodiment, since the expression of the matrix M in equation (59) is very complex, we analyze the modulus of the growth factor by numerically solving M. Here, Δx=5Δy=Δz=5mm, Δt=min[Δx / c0,Δz / c0]. The relevant parameters in the CPML technology are defined as follows:

[0315] Where d is the thickness of the PML layer, w0 is the position of the PML layer close to the computational domain, m = 4, κ max =7, α max =0.05. In addition, sin(k w Δw / 2) is still unknown. According to the symmetry of the sine function, k w The value range of Δw is [0,π]

[0316] like Figures 4a to 4d As shown, in this embodiment, the numerical characteristic analysis is performed on the CPML plane area. First, consider the plane area perpendicular to the x-axis in the CPML area. At this time, there is only the auxiliary variable ψ along the x-axis direction. eyx , ψ ezx and ψ hyx , ψ hzx Need to be updated. Therefore, the coordinate scaling factor should be chosen as S y =1 and S z = 1. Next, the CPML parameters are substituted into equation (59) to calculate the modulus of the growth factor ζ. Considering that the parameter values ​​in CPML vary with the number of layers, we selected four layers from the outermost layer to the innermost layer of the CPML region for stability analysis. Figures 4a to 4d The modulus of the growth factor ζ is shown. It can be found that all the moduli of the growth factor ζ are less than 1, indicating that the iterative solution of the field component in this CPML region is stable.

[0317] like Figure 5a to Figure 5d As shown, in this embodiment, numerical characteristic analysis is performed on the CPML edge region. In the CPML edge region, for the edge parallel to the z-axis, the coordinate scaling factor should be selected as and S z = 1. By substituting the above parameters into equation (59), the distribution of the growth factor ζ modulus is obtained as Figure 5a to Figure 5d It can be found that the modulus values ​​of all growth factors ζ are not greater than 1.

[0318] like Figures 6a to 6d As shown, in this embodiment, the numerical characteristic analysis PML is performed on the CPML corner vertex area. In the CPML corner vertex area, the coordinate scaling factor can be selected as as well as The growth factor ζ modulus distribution is as follows Figures 6a to 6d As shown, the modulus values ​​of all growth factors ζ are no greater than 1.

[0319] Through the above stability analysis, it can be found that the CPML technology based on the CD-HIE-FDTD format we constructed can achieve the stability condition that the time step is only related to the grid size in one direction.

[0320] Example 3

[0321] In this embodiment, a variety of examples are given to verify the correctness and effectiveness of the CD-HIE-FDTD numerical solution method based on the CFS-PML technology provided by the present invention.

[0322] like Figure 7 As shown, in this embodiment, in order to verify the effectiveness of CPML in the proposed CDHIE-FDTD method, an open domain electromagnetic problem containing a PEC sheet structure is simulated. As shown in Figure 6, the entire calculation domain contains 60×100×60 units, and the thickness of the CPML is 10 layers. The unit size is set to Δx=5Δy=Δz=5mm. Among them, the metal sheet is placed in the center of the calculation domain, which occupies 32×20 units on the xz plane. The excitation source is set at point A (31, 51, 31), and the time domain expression of the excitation source is

[0323]

[0324] Where t = nΔt, τ=120Δt / CFLN,t0=τ。 The detected field component is located at point B (49, 48, 31), 2 grid cells away from the CPML boundary. The parameters of CPML are: m=4,κ max =7, α max =0.05.

[0325] In order to verify the absorption performance of CPML, the relative reflection error function is used to quantify

[0326]

[0327] Among them, E cpml represents the electric field value obtained at the observation point under the condition of loading 10 layers of CPML, E ref It is the electric field value at the detection point obtained in a larger calculation space.

[0328] like Figure 8As shown in the figure, the relative reflection error of CPML under different CFLN conditions is shown. It can be found that when CFLN is 1, the relative reflection error remains below -65dB. When the maximum time step specified by the CDHIE-FDTD method is used (CFLN = 5.19), the relative reflection error is around -50dB, which verifies that the CPML technology based on the frog-leaping CDHIE-FDTD format has relatively stable and efficient absorption performance. At the same time, this also verifies in numerical form that the CPML technology constructed in this paper will not affect the stability conditions of the CDHIE-FDTD method.

[0329] like Fig. 9 As shown in the figure, the unit structure of the frequency selective surface is shown. The number of grids occupied by the computational domain is 20×120×20. Periodic boundaries are used in the x and z directions, and the unit period is T x =T z =10mm, and 10 layers of CPML are used for truncation in the y direction. The outer ring side length of the square metal patch is 8mm, and the inner ring side length is 6mm. The relative dielectric constant ε of the dielectric substrate r =2.2, thickness is 0.5mm. The simulation unit grid is set to Δx=5Δy=Δz=0.5mm, and the time step of the traditional FDTD method is set to The time step of the leapfrog CDHIE-FDTD method is set to Δt CDHIE-FDTD =Δx / c0. The Gaussian plane wave source is introduced by the total field / scattered field technique, and the Gaussian pulse form is Among them, τ=280Δy / c0, t0=1.2τ. Fig. 9 The transmission characteristic curves of the frequency selective surface calculated by the traditional FDTD method, CDHIE-FDTD method and the commercial software Computer Simulation Technology (CST) are shown. It can be seen from the figure that the proposed method is consistent with the results of the traditional FDTD method and CST calculation (the resonant frequency of FSS is 12.6GHz). In addition, the CPU calculation time of the FDTD method and the CDHIE-FDTD method is 128.704 seconds and 49.055 seconds, respectively. This shows that compared with the traditional FDTD method, the CDHIE-FDTD method has higher computational efficiency when solving electric field problems with fine structures.

[0330] See also Fig.10 , transmission coefficients calculated by the traditional FDTD method, the CD-HIE-FDTD method, and the commercial software CST; see Fig.11 In order to further verify the effectiveness of the proposed method, this section analyzes the thin slot shielding chassis. Fig.10The size of the metal shielding box is 30×12×30cm 3 The grid size of the computational area is set to Δx = 2.5 cm, Δy = 0.6 cm, Δz = 3 cm. It is assumed that all surfaces of the chassis are PEC with zero thickness and the size of the slot is 20 × 3 cm. 2 , located in the center of the chassis. A Gaussian plane wave is loaded in front of the slot. The time domain form of the Gaussian excitation source is the same as that of Example B, where τ = 0.637ns and t0 = 3τ. The calculation area is truncated using 10 layers of CPML.

[0331] like Fig.12 As shown, in this embodiment, the shielding effectiveness calculation results of the chassis obtained by the traditional FDTD method (CFLN = 1) and the CD-HIE-FDTD method (CFLN = 1, 3, 4.64, CFLN = 4.64 corresponds to the maximum discrete time step that the frog-leaping CD-HIE-FDTD method can take) are shown respectively. It can be found that the calculation results obtained by the CD-HIE-FDTD method with different time steps are consistent with the calculation results of the traditional FDTD method, which once again verifies the correctness of this method.

[0332] At the same time, Table II gives a comparison of the computational efficiency of the two methods. It can be seen that when the CD-HIE-FDTD method takes CFLN = 3, its simulation time is less than that of the traditional FDTD method. When the CD-HIE-FDTD method takes a larger time step, the computational efficiency is more obvious.

[0333] Table II CPU time calculated by different methods

[0334] Method Time step CPU time Traditional FDTD 5000 341.53s CD-HIE-FDTD (CFLN=1) 5000 862.23s CD-HIE-FDTD(CFLN=3) 1666 304.65s CD-HIE-FDTD(Max Time) 1145 207.87s

[0335] Next, the method of the present invention is used to calculate Fig.13 The four different slot positions of the chassis are shown in the figure to analyze the impact of the slot position on the chassis shielding performance. As shown in the figure, the gap area of ​​the four different slot positions is the same (60cm 2 ). Fig.14 The shielding effectiveness of four different slot chassis is demonstrated. The results show that although the gap area is the same, the effects of different slot designs on the shielding effectiveness of the chassis are significantly different. Therefore, the method proposed in this paper can provide an efficient and accurate simulation platform for designing and optimizing thin slot chassis with excellent shielding effectiveness.

[0336] Dual-band Printed F Antenna

[0337] To further evaluate the accuracy and effectiveness of the proposed method, we simulated the electromagnetic characteristics of a printed inverted-F antenna (PIFA) supporting dual-band operation, such as Fig.15As shown. The metal floor and the inverted F structure are etched on both sides of the dielectric substrate. The length of the dielectric substrate along the x and z directions is 40 mm, the thickness along the y direction is 0.787 mm, and the relative dielectric constant is 2.2. The entire calculation area occupies 140×59×140 Yee grids in the x, y and z directions, of which the x and z directions contain 10 layers of CPML, and the y direction contains 15 layers of CPML. The grid sizes are set to Δx=0.4mm, Δy=0.087mm, Δz=0.4mm respectively. The port is equipped with a 50Ω resistive voltage source and Gaussian pulse excitation is used to meet the accuracy requirement of 20 grids corresponding to the minimum wavelength. Fig.16 and 17 The time domain waveforms of the sampled voltage and current near the voltage source using the traditional FDTD method and the CD-HIE-FDTD method (CFLN = 4.79, which is the maximum discrete time step of the CD-HIE-FDTD method) are shown. It can be found that the calculation results of the two methods are consistent, where the CPU calculation time of the FDTD method and the CDHIE-FDTD method are 59.52min and 34.43min, respectively, which verifies the effectiveness of the proposed method. In addition, the time domain waveforms calculated by the two methods decay to zero after 2ns, further verifying the excellent absorption effect of CPML on electromagnetic waves in the CD-HIE-FDTD method. Finally, Figures 18 to 20 The S parameters of the traditional FDTD method and the CD-HIE-FDTD method (CFLN = 4.79) in the frequency domain, as well as the radiation patterns cut in the xy plane at 2.4 GHz and 5.8 GHz are shown respectively. It can be found that the calculation results of the proposed method are consistent with those of the traditional FDTD method, which further verifies the effectiveness of the proposed method.

[0338] In summary, the present invention introduces CFS into the Maxwell equations and constructs a time evolution matrix, then performs time-space operator splitting on the matrix to construct a two-step iterative solution with divergence properties. In order to improve the efficiency of numerical calculations, the present invention adopts the time stepping strategy of the Leapfrog format to eliminate the intermediate variables in the time domain and obtain a single-step differential discretization solution. The present invention uses recursive convolution technology to process the convolution operator in the discrete formula to obtain the iterative formula of the electromagnetic field components and the auxiliary variables in the PML region. In addition, for the above-mentioned calculation framework, the analysis of stability, divergence properties and complexity is completed, which can effectively realize the simulation of solving complex open-domain electromagnetic problems.

[0339] The present invention solves the technical problems existing in the prior art of low numerical calculation efficiency, existence of intermediate variables in the time domain, and the need to improve the correctness and effectiveness of the integrated solution to complex open-domain electromagnetic problems containing fine structures.

[0340] The above embodiments are only used to illustrate the technical solutions of the present invention, rather than to limit the same. Although the present invention has been described in detail with reference to the aforementioned embodiments, those skilled in the art should understand that the technical solutions described in the aforementioned embodiments may still be modified, or some of the technical features may be replaced by equivalents. However, these modifications or replacements do not deviate the essence of the corresponding technical solutions from the spirit and scope of the technical solutions of the embodiments of the present invention.

Claims

1. An open-domain electromagnetic target solution method for unidirectional fine structures, characterized in that: The method comprises: S1. The PML stretch factor is introduced into Maxwell's equations to construct the time evolution matrix of the electromagnetic field components. S2, using a hybrid explicit and implicit difference strategy to process the matrix of spatial derivatives in the time evolution matrix; reconstructing the time evolution matrix to obtain a two-step time domain discrete equation; S3. Based on the two-step time domain discrete equation and the single-step Leapfrog format, a unified numerical calculation format is derived to combine the CPML boundary conditions with the single-step frog leap CDHIE-FDTD method, and the single-step Leapfrog format is applied to the preset sub-time step logic to perform CD-HIE-FDTD numerical solution.

2. The open-domain electromagnetic target solving method for unidirectional fine structure according to claim 1 is characterized in that: The S1 includes: S11. Use the following logic to express the curl equation of the electromagnetic field: In the formula, S w (w=x,y,z) is the scaling factor of PML: S12, converting the formula (1), the formula (2) and from the frequency domain to the time domain and writing them in matrix form: Among them, u=[E x ,E y ,E z ,H x ,H y ,H z ] T , the spatial curl matrix S13. Assuming that there is no current source and magnetic current source, the following logic is used to obtain the exact solution of equation (4) from time t=0 to time t=Δt: u(Δt)=exp(RΔt)u(0)(7) The R in the formula (7) can be split into two sub-matrices to obtain the time evolution matrix.

3. The open-domain electromagnetic target solving method for unidirectional fine structure according to claim 1 is characterized in that: The S2 includes: S21, using the hybrid implicit and explicit difference strategy in different directions to process the sub-matrices R1 and R2 of the time evolution matrix; using implicit difference in the y direction and explicit difference in the x and z directions, and expressing the sub-matrices R1 and R2 using the following logic: in: S22. Express the exponential matrix exp(RΔt) in equation (7) as the following matrix multiplication form: S23, reconstructing the time evolution matrix in the formula (11) to obtain the following exponential matrix multiplication scheme, which satisfies the divergence property: By Taylor series expansion, the equation (13) is expressed as: According to the equation (8) and the equation (14), the time-space discrete form of the electromagnetic field component is obtained: S24. Using the following logic, the formula (15) is divided into two parts for calculation: as well as S25. Define auxiliary variable Q = [Q ex ,Q ey ,Q ez ,Q hx ,Q hy ,Q hz ] T , dealing with the inversion problem of the operator matrix in the two sub-time steps, expressing the two-step time domain discrete equation.

4. The open-domain electromagnetic target solving method for unidirectional fine structure according to claim 3 is characterized in that: In S25, the first sub-time step is derived using the following logic: Using the following logic, the second sub-time step is derived: According to the first sub-time step and the second sub-time step of equations (16) to (19), the two-step time domain discrete equation of the CDHIE-FDTD method is described.

5. The open-domain electromagnetic target solving method for unidirectional fine structure according to claim 1, characterized in that: The S3 includes: S31, updating the electric field component in the time period from n-1 / 2 to n+1 / 2; S32. Update the magnetic field components in the time period from n to n+1. In S32, equation (21) is substituted into equation (19) to obtain: Substituting equation (17) into equation (34), we obtain: The same processing method as equation (23) is used. By further deducing equation (35), the time step equation of the magnetic field component at time n+1 is obtained: in, Reverse equation (19) by one sub-time step and substitute it into equation (20), and we get: Substituting equation (17) into equation (37), we obtain: The time step equation of the magnetic field at time n is further obtained by the above equation (38): Subtracting equation (39) from equation (36) yields: Equation (40) is a time-stepping equation for updating the magnetic field components based on a single-step leapfrog scheme. S33. Use the following logic to define auxiliary variables: and S34. Based on the CDHIE-FDTD method, the electromagnetic field components in the CPML region are iterated: where ψ e = [ψ exy , ψ exz , ψ eyz , ψ eyx , ψ ezx , ψ ezy T , ψ h = [ψ hxy , ψ hxz , ψ hyz , ψ hyx , ψ hzx , ψ hzy T .​​ 6. The method for solving an open-domain electromagnetic target for a unidirectional fine structure according to claim 5, characterized in that: The S31 includes: S311. Back one time step from equation (19) in the second sub-time step and substitute it into equation (16) in the first sub-time step to obtain: Substituting equation (17) of the first time step into equation (18) of the second sub-time step and equation (20), we obtain: The field component in equation (19) is stepped back by one time step and substituted into equation (22), and we get: S312. According to formula (23), the updated equation of the electromagnetic field component in the time period from n to n+1 / 2 can be obtained: According to equation (24), the iterative formula of the electric field component from n to n+1 / 2 is obtained: S313, define auxiliary variable F w ; S314, write the convolution term in the formula (25) as a discrete circular convolution form to obtain the iterative formula of the field component in the CPML region: S315. Now substitute equation (17) into equation (21) to obtain the electric field component at time n–1 / 2: Substituting equation (19) into equation (30) and going back one time step, we get: Processing equation (31), the iterative formula of the electric field component at time n–1 / 2 in the CPML region is obtained: S316. Subtract the formula (32) from the formula (27) to obtain an iterative formula for the electric field component in the time period n–1 / 2 to n+1 / 2, which is used as a time step equation for updating the electric field component based on the single-step frog leapfrog format:

7. The method for solving an open-domain electromagnetic target for a unidirectional fine structure according to claim 6, characterized in that: In the S313, or w=x,y,z. According to formula (3), we get: In the formula, 8. The open-domain electromagnetic target solving method for unidirectional fine structure according to claim 6, characterized in that: In the S314, C is used in the formula (27) e ,D e and C h ,D h , represents the spatial operation operator of the electric field region and the magnetic field region; The auxiliary variables in equation (27) are iterated using the following logic:

9. An open-domain electromagnetic target solving system for unidirectional fine structures, characterized in that: The system comprises: The time evolution matrix construction module is used to introduce the PML's stretching factor into Maxwell's equations to construct the time evolution matrix of the electromagnetic field components; A two-step time domain discrete equation acquisition module is used to process the matrix of the space derivative in the time evolution matrix by adopting a hybrid explicit and implicit difference strategy; reconstruct the time evolution matrix to obtain a two-step time domain discrete equation, and the two-step time domain discrete equation acquisition module is connected to the time evolution matrix construction module; The CPML boundary condition and single-step frog leapfrog combined numerical solution module is used to derive a unified numerical calculation format based on the two-step time domain discrete equation and the single-step Leapfrog format, so as to combine the CPML boundary condition and the single-step frog leapfrog CDHIE-FDTD method, apply the single-step Leapfrog format to the equations (16) to (19) to perform CD-HIE-FDTD numerical solution, and the CPML boundary condition and single-step frog leapfrog combined numerical solution module is connected to the two-step time domain discrete equation acquisition module.

Images