A numerical calculation method of electromagnetics based on Chebyshev FDTD

Abstract

The present invention discloses an electromagnetic numerical calculation method based on Chebyshev FDTD, firstly, the time interval of numerical calculation is equally divided into a plurality of smaller time segments, the length of the time segment being greater than the support boundary of the excitation signal; starting from the first time segment, numerical calculation is performed on the current time segment using Chebyshev FDTD, and the final value of the obtained time domain result is used as the initial value of the next time segment, and the next time segment is continuously numerically calculated using Chebyshev FDTD, and so on, so as to obtain the time domain results of all time segments; finally, the time domain results of all time segments are synthesized into a complete time domain result within the time interval. Compared with the prior art, the present invention effectively improves the computational efficiency of electromagnetic numerical calculation and reduces the consumption of computing resources by means of segmented calculation.

Description

Technical Field

[0001] The present invention relates to an electromagnetic numerical calculation method, in particular to an electromagnetic numerical calculation method based on Chebyshev FDTD (Finite Difference Time Domain). Background Art

[0002] The finite-difference time-domain (FDTD) method is widely used in numerical calculations in the electromagnetic field. This method takes Maxwell's equations as the starting point, converts the two curl equations in the equations into a differential format, discretely samples the time-domain electromagnetic field in space and time, obtains an iterative equation set, and numerically simulates the propagation law of electromagnetic waves and the interaction between electromagnetic waves and scatterers. The traditional FDTD method is limited by the CFL (Courant-Friedrich-Lewy) stability condition, and the time step is limited by the minimum grid size. When solving electromagnetic model problems with fine structures, the time step needs to be set very small to ensure the accuracy of the calculation results. At the same time, if the calculation time interval remains unchanged, the number of time steps will increase accordingly, reducing the calculation efficiency.

[0003] In order to solve the above problems of the traditional FDTD method, the inventor of the present application proposed an FDTD algorithm based on Chebyshev polynomials (hereinafter referred to as Chebyshev FDTD). This method uses Chebyshev polynomial basis functions to realize the expansion and reconstruction of the signal, analyzes the differential relationship between the basis functions to construct the Chebyshev differential matrix to replace the time domain differential, eliminates the time variable, and realizes the conversion of Maxwell's equations from the time domain to the Chebyshev domain. In space, the discrete method of traditional FDTD is adopted to construct the implicit equations of the electric field and magnetic field in the Chebyshev domain, solve the expansion coefficients of the electric field and magnetic field in the Chebyshev domain, and then reconstruct the time domain results. Compared with the traditional FDTD method, the calculation accuracy is guaranteed while the calculation efficiency is improved, and unconditional stable calculation is realized. For details of Chebyshev FDTD, please refer to the literature Unconditionallystable FDTD method based on Chebyshev polynomials—Chebyshev (CS) FDTD (C. Li, Z.-Y. Huang, 2021 International Applied Computational Electromagnetics Society (ACES-China) Symposium, 2021, pp. 1-2, doi: 10.23919 / ACES-China52398.2021.9581637.).

[0004] Although Chebyshev FDTD can improve both calculation accuracy and efficiency compared to traditional FDTD, it needs to increase the order of basis functions to ensure the stability and accuracy of calculation results when facing electromagnetic problems with long calculation time and large support intervals. However, the increase in order will increase the consumption of calculation memory, thereby reducing calculation efficiency. Therefore, it is necessary to further improve it to improve calculation efficiency and reduce consumption of calculation resources. Summary of the invention

[0005] The technical problem to be solved by the present invention is to overcome the shortcomings of the existing Chebyshev FDTD and provide an electromagnetic numerical calculation method based on Chebyshev FDTD, which improves the computational efficiency of electromagnetic numerical calculation and reduces the consumption of computing resources through segmented calculation.

[0006] The present invention specifically adopts the following technical solutions to solve the above technical problems:

[0007] An electromagnetic numerical calculation method based on Chebyshev FDTD, firstly, the time interval of numerical calculation is equally divided into a plurality of smaller time segments, the length of the time segment being greater than the support boundary of the excitation signal; starting from the first time segment, numerical calculation is performed on the current time segment using Chebyshev FDTD, and the final value of the obtained time domain result is used as the initial value of the next time segment, and the numerical calculation is continued to be performed on the next time segment using Chebyshev FDTD, and so on, so as to obtain the time domain results of all time segments; finally, the time domain results of all time segments are synthesized into a complete time domain result within the time interval.

[0008] Preferably, the Chebyshev differential matrix α of the Chebyshev FDTD in each time segment is as follows:

[0009]

[0010] Where q is the polynomial order, q = 0, 1, 2, ..., T f is the length of the time segment.

[0011] Furthermore, the method specifically comprises the following steps:

[0012] Step 1: Select a suitable Chebyshev polynomial order according to the electromagnetic problem to be solved, establish a Chebyshev differential matrix, and divide the time interval of the numerical calculation into a plurality of smaller time segments, wherein the length of the time segment is greater than the support boundary of the excitation signal;

[0013] Step 2: Grid division, introduction of material parameters and excitation sources to initialize the following coefficient matrix equation and determine the right side b of the equation i The electric field expansion coefficient of the Chebyshev domain in the first time segment is solved by the chase method:

[0014] α l(i-1) E| i-1 +α m(i) E| i +α r(i+1) E| i+1 =b i

[0015]

[0016]

[0017]

[0018]

[0019] Among them, E| i-1 、E| i 、E| i+1 is the electric field expansion coefficient of the Chebyshev domain at the adjacent position to be determined, ε i , μ i+1 / 2 , is the relative permittivity, conductivity, magnetic permeability and magnetic permeability after discretization; E y | i , H z | i+1 / 2 , J y | i 、M z | i+1 / 2 are the discretized electric field component, magnetic field component, current source density, and magnetic current source density; is the initial value of the electric field and magnetic field after discretization; i represents the position number of the spatial discrete point, Δx is the interval between adjacent electric fields or adjacent magnetic fields;

[0020] Step 3: Substitute the obtained Chebyshev domain electric field expansion coefficient into the following formula to solve the Chebyshev domain magnetic field expansion coefficient:

[0021]

[0022] Among them, H z | i+1 / 2 is the desired Chebyshev field expansion coefficient, E y | i+1 、E y | i is the solved Chebyshev domain electric field expansion coefficient at the adjacent position, M z | i+1 / 2 is the density expansion coefficient of the magnetic flux source, is the initial value of the magnetic field in this time segment, μ i+1 / 2 α(:,1) is the known coefficient matrix related to the medium parameters and the differential matrix, Δx is the interval between adjacent electric fields;

[0023] Step 4: Use the Chebyshev domain expansion coefficients and basis functions to reconstruct the electromagnetic field time domain results, and take the final values ​​of the electric field and magnetic field time domain results as the initial values ​​for the next time segment calculation to update b i Item, and then solve the Chebyshev domain electric field expansion coefficient for the next time segment;

[0024] Step 5, repeating steps 3 and 4 until all electromagnetic field quantities of all time segments are calculated, and combining the time domain results of all time segments into a complete time domain result within the time interval.

[0025] Compared with the prior art, the technical solution of the present invention has the following beneficial effects:

[0026] 1. It inherits the advantages of the original Chebyshev FDTD algorithm and makes up for the deficiency of excessive memory consumption during long-term calculation.

[0027] 2. Reduced memory usage while improving computing efficiency. BRIEF DESCRIPTION OF THE DRAWINGS

[0028] Figure 1 It is a schematic diagram of the basic principle of segmented calculation;

[0029] Figure 2 It is a basic flow chart of the electromagnetic numerical calculation method of the present invention;

[0030] Figure 3 It is a simplified diagram of the configuration of the numerical calculation model of the one-dimensional resonant cavity;

[0031] Figure 4 It is a comparison of the relative errors of the electric field time domain waveforms at the measurement point between the method of the present invention and the traditional FDTD method. DETAILED DESCRIPTION

[0032] In order to overcome the shortcomings of the Chebyshev FDTD algorithm, the solution of the present invention is to introduce a segmented calculation method based on the algorithm.

[0033] The electromagnetic numerical calculation method based on Chebyshev FDTD proposed in the present invention is as follows:

[0034] First, the time interval of numerical calculation is equally divided into multiple smaller time segments, and the length of the time segment is greater than the support boundary of the excitation signal; starting from the first time segment, the current time segment is numerically calculated using Chebyshev FDTD, and the final value of the time domain result obtained is used as the initial value of the next time segment, and the next time segment is continuously numerically calculated using Chebyshev FDTD, and so on, so as to obtain the time domain results of all time segments; finally, the time domain results of all time segments are synthesized into a complete time domain result within the time interval.

[0035] In order to facilitate public understanding, the technical solution of the present invention is described in detail below with reference to the accompanying drawings:

[0036] The existing Chebyshev FDTD uses the differential relationship of Chebyshev polynomials to derive the relationship between the Chebyshev expansion coefficients of the signal and the Chebyshev expansion coefficients after the signal is differentiated, establishes the Chebyshev integral matrix, and then obtains the Chebyshev differential matrix α by inversion. The Chebyshev differential matrix α is applied to the time-domain Maxwell equations to solve the electromagnetic field.

[0037] In order to overcome the problem that the existing Chebyshev FDTD calculation is time-consuming and occupies more computing resources, the present invention improves it and introduces a segmented calculation method: first, the complete time interval is divided into several smaller time segments before calculation, and the initial conditions of the electric field and magnetic field are introduced in the constructed Chebyshev domain electromagnetic field solution equation, and then the expansion coefficients of the electric field and magnetic field of the first segment are calculated and the time domain results are reconstructed, and the field values ​​of the electric field and magnetic field at the last moment in the time interval are recorded and converted to the Chebyshev domain, and the calculation is continued as the initial conditions of the electric field and magnetic field in the next time segment until the calculation of the electric field and magnetic field in all time segments is completed, and then integrated into a complete time domain result. This method can effectively reduce the order of the required basis functions, reduce the memory requirements of the calculation, and improve the calculation efficiency.

[0038] The key problem in realizing the piecewise calculation is how to establish the differential matrix of Chebyshev polynomials and update the initial value conditions of Maxwell equations in the piecewise calculation, which is described in detail below.

[0039] The Chebyshev differential matrix is ​​obtained by inverting the integral matrix. The construction of the integral matrix refers to the relationship between the Chebyshev expansion coefficients of the signal and the Chebyshev expansion coefficients after the signal is differentiated:

[0040]

[0041] in u q is the expansion coefficient of the original signal, is the differential signal expansion coefficient, q is the polynomial order, q=0,1,2,…,Tf is the length of each time segment. According to formula (1), the Chebyshev integral matrix α is:

[0042]

[0043] The segmented calculation refers to the time interval required for calculation, such as [0,T n ], and divide it into [0,T 1 ]、[T 1 ,T 2 ]、[T 2 ,T 3 ], …, [T n-1 ,T n ], where T 1 >T m , T m is the support boundary of the excitation signal (if the excitation signal is g(t), when t>T m When g(t) = 0 or approaches 0 infinitely). First, calculate the 1 ], and obtain the final value E(T 1 ) and H(T 1 ), and then calculate the next interval [T 1 ,T 2 ] exceeds the support boundary, so the original excitation signal no longer works, and E(T 1 ) and H(T 1 ) as the initial value condition of this interval and the new stimulus to continue the calculation, and finally complete [0,T n The basic principle of all calculations on ] is as follows Figure 1 shown.

[0044] In the specific calculation, the one-dimensional uniform lossy medium time domain Maxwell equations including initial value conditions are first used as the calculation basis, which can be written as:

[0045]

[0046]

[0047] where ε(x), μ(x), σ e (x),σ m (x) is the relative permittivity, magnetic permeability, electrical conductivity, and magnetic permeability, E y (x,t),H z (x,t),J y (x,t),M z (x, t) is the electric field component, magnetic field component, current source density, magnetic current source density, E y (x,t 0 ), Hz (x,t 0 ) is the initial value of the electric field and magnetic field, and δ(t) represents the impulse function.

[0048] Substitute the Chebyshev differential matrix α into formula (3) and formula (4) instead of Eliminating the time variable t, we can get

[0049]

[0050]

[0051] Where E y (x), H z (x), J y (x), M z (x) is no longer related to time t and is transformed into an expansion coefficient under the Chebyshev basis function; u(t) is a step function, corresponding to the zero-order Chebyshev basis function, and its expansion coefficient in the Chebyshev domain is u q =[1,0,0,…,0] T , is a q×1-dimensional vector, then the expansion coefficient of δ(t) is αu q =α(:,1), α(:,1) represents the first column of the differential matrix α.

[0052] Discretizing formula (5) and formula (6) in space yields

[0053]

[0054]

[0055] in ε i , μ i+1 / 2 , is the relative permittivity, conductivity, magnetic permeability and magnetic permeability after discretization; E y | i , H z | i+1 / 2 , J y | i 、M z | i+1 / 2 are the discretized electric field component, magnetic field component, current source density, and magnetic current source density; is the initial value of the electric field and magnetic field after discretization; i represents the position number of the discrete point in space, and Δx is the interval between adjacent electric fields or adjacent magnetic fields.

[0056] Substituting formula (8) into formula (7) to eliminate the magnetic field component first, then putting the terms related to the electric field component on the left side of the equation and the other terms on the right side of the equation, we can get the coefficient matrix equation:

[0057] α l(i-1) E| i-1 +α m(i) E| i +α r(i+1) E| i+1 =b i (9)

[0058]

[0059] Among them, E| i-1 、E| i 、E| i+1 is the electric field expansion coefficient of the Chebyshev domain at the adjacent position to be determined, α l(i-1) , α m(i) , α r(i+1) is the coefficient matrix related to the differential matrix α and the medium parameters, b i With current source J y | i , magnetic current source M z | i+1 / 2 , initial value of electric field Initial value of magnetic field Related parameters.

[0060] The specific process of the electromagnetic numerical calculation method of the present invention is as follows: Figure 2 As shown, the following steps are included:

[0061] Step 1: Select a suitable Chebyshev polynomial order according to the electromagnetic problem to be solved, establish a Chebyshev differential matrix, and divide the time interval of the numerical calculation into a plurality of smaller time segments, wherein the length of the time segment is greater than the support boundary of the excitation signal;

[0062] Step 2: Grid division, introduction of material parameters and excitation sources to initialize the coefficient matrix equation in equation (9) and determine the right side b of the equation i Item, use the chasing method to solve the Chebyshev domain electric field expansion coefficient of the first time segment;

[0063] Step 3, substitute the obtained Chebyshev domain electric field expansion coefficient into equation (8) to solve the Chebyshev domain magnetic field expansion coefficient;

[0064] Step 4: Use the Chebyshev domain expansion coefficients and basis functions to reconstruct the electromagnetic field time domain results, and take the final values ​​of the electric field and magnetic field time domain results as the initial values ​​for the next time segment calculation to update b i Item, and then solve the Chebyshev domain electric field expansion coefficient for the next time segment;

[0065] Step 5, repeating steps 3 and 4 until all electromagnetic field quantities of all time segments are calculated, and combining the time domain results of all time segments into a complete time domain result within the time interval.

[0066] The technical effect of the method of the present invention is verified by a specific numerical example below. Figure 3 The figure shows the model for calculating the electric field component at the sampling point after the current source is introduced into the one-dimensional resonant cavity containing a lossy dielectric. The total length of the cavity is 10 mm, with PEC boundaries at both ends and the current source J y 2mm from the left boundary, sampling point E y Located in J y 1mm to the right, there is a lossy dielectric plate in the middle, with a thickness of 0.6mm, relative dielectric constant ε = 4, and conductivity σ e =0.1. The entire area is divided into 100 uniform grids with a grid size of 0.1 mm; the time segment is divided into 3 segments, each of which is 3 ns long; the order of the basis function is selected as 200; the excitation source is a Gaussian pulse: where t 0 represents the time shift, τ represents the pulse width, and t 0 =4.5τ.

[0067] The method of the present invention, the original Chebyshev FDTD method and the traditional FDTD method are respectively used to perform numerical calculations and the obtained results are compared. Figure 4 The method of the present invention and the traditional FDTD method at the sampling point E y The comparison of the electric field time domain waveform and the relative error. The figure shows that the calculation results of the two methods in each section are very consistent, and the relative error is below -50dB, so the calculation results of the method of the present invention are correct. The computing resources consumed by the three methods are shown in Table 1. It can be seen from Table 1 that the method of the present invention greatly reduces the memory usage and improves the computing efficiency compared with the original Chebyshev FDTD method.

[0068] Table 1 Comparison of computing resources

[0069] Δt Memory (Mb) CPU time(s) FDTD 0.3ps 0.87 8.8 Original Chebyshev FDTD \ 122.07 1.25 Method of the present invention \ 61.03 0.95

Claims

1. An electromagnetic numerical calculation method based on Chebyshev FDTD, characterized in that: First, the time interval of numerical calculation is equally divided into a plurality of smaller time segments, the length of the time segment being greater than the support boundary of the excitation signal; starting from the first time segment, numerical calculation is performed on the current time segment using Chebyshev FDTD, and the final value of the obtained time domain result is used as the initial value of the next time segment, and the numerical calculation is continued using Chebyshev FDTD on the next time segment, and so on, so as to obtain the time domain results of all time segments; finally, the time domain results of all time segments are synthesized into a complete time domain result within the time interval; The method specifically comprises the following steps: Step 1: Select a suitable Chebyshev polynomial order according to the electromagnetic problem to be solved, establish a Chebyshev differential matrix, and divide the time interval of the numerical calculation into a plurality of smaller time segments, wherein the length of the time segment is greater than the support boundary of the excitation signal; Step 2: Grid division, introduction of material parameters and excitation sources to initialize the following coefficient matrix equation and determine the right side b of the equation i The electric field expansion coefficient of the Chebyshev domain in the first time segment is solved by the chase method: a l(i-1) E| i-1 +a m(i) E| i +a r(i+1) E| i+1 =b i Among them, E| i-1 、E| i 、E| i+1 is the electric field expansion coefficient of the Chebyshev domain at the adjacent position to be determined, ε i , μ i+1 / 2 , is the relative permittivity, conductivity, magnetic permeability and magnetic permeability after discretization; E y | i , H z | i+1 / 2 , J y | i 、M z | i+1 / 2 are the discretized electric field component, magnetic field component, current source density, and magnetic current source density; is the initial value of the electric field and magnetic field after discretization; i represents the position number of the spatial discrete point, Δx is the interval between adjacent electric fields or adjacent magnetic fields; Step 3: Substitute the obtained Chebyshev domain electric field expansion coefficient into the following formula to solve the Chebyshev domain magnetic field expansion coefficient: Among them, H z | i+1 / 2 is the desired Chebyshev field expansion coefficient, E y | i+1 、E y | i is the solved Chebyshev domain electric field expansion coefficient at the adjacent position, M z | i+1 / 2 is the density expansion coefficient of the magnetic flux source, is the initial value of the magnetic field in this time segment, μ i+1 / 2 α(:,1) is the known coefficient matrix related to the medium parameters and the differential matrix, Δx is the interval between adjacent electric fields; Step 4: Use the Chebyshev domain expansion coefficients and basis functions to reconstruct the electromagnetic field time domain results, and take the final values ​​of the electric field and magnetic field time domain results as the initial values ​​for the next time segment calculation to update b i Item, and then solve the Chebyshev domain electric field expansion coefficient for the next time segment; Step 5, repeating steps 3 and 4 until all electromagnetic field quantities of all time segments are calculated, and combining the time domain results of all time segments into a complete time domain result within the time interval.

2. The electromagnetic numerical calculation method based on Chebyshev FDTD as claimed in claim 1, characterized in that: The Chebyshev differential matrix α of the Chebyshev FDTD in each time segment is as follows: Where q is the polynomial order, q = 0, 1, 2, L, T f is the length of the time segment.

Images