Implicit finite volume calculation method, equipment, medium and system for time-varying electromagnetic field

Abstract

The present invention discloses a time-varying electromagnetic field implicit finite volume calculation method, equipment, medium and system, which belongs to the field of numerical solution of time-domain electromagnetic finite volume method, including: iterative solution of Maxwell equations time-varying electromagnetic field by implicit double time step method based on time iterative advancement of conserved electromagnetic field and spatial flux residual; the outer layer of the simulation model is a physical time step cycle until the calculation converges; the inner layer of the simulation model is a virtual time step sub-iteration cycle until the sub-iteration converges; in each virtual time sub-iteration process, each grid block and each grid unit is sequentially subjected to spatial flux calculation and implicit iterative solution calculation, and the next level virtual time sub-iteration step number of conserved electromagnetic field values ​​is updated. The present invention can solve large-scale electromagnetic scattering problems of arbitrarily complex shapes, high-frequency electrical large-size targets.

Description

Technical Field

[0001] The present invention relates to the field of numerical solution of time-domain electromagnetic finite volume method, and more specifically, to a time-varying electromagnetic field implicit finite volume calculation method, equipment, medium and system. Background Art

[0002] The time domain method can compatibly, dynamically and intuitively observe and simulate the interaction process between electromagnetic waves and objects. The electromagnetic waveforms that can be processed include both continuous waves and broadband pulse waveforms. It is widely used in target radiation characteristics, antenna design and electromagnetic compatibility.

[0003] The time-varying electromagnetic field in time-domain calculation satisfies the time-domain Maxwell equations and has the same hyperbolic characteristics as the Euler equations in fluid mechanics. Therefore, computational fluid dynamics (CFD) technology can be applied to electromagnetic field calculations. Among them, the finite difference time domain (FDTD) method in 1966 is the most well-known, and the finite volume time domain (FVTD) method appeared in the 1990s.

[0004] The traditional finite-difference time-domain method has the biggest flaw because it is an explicit method. Due to stability limitations, the calculation must use the global minimum time step, which greatly increases the amount of time-domain electromagnetic field calculations. In contrast, the implicit method can relax the stability restrictions on the calculation step, but the time accuracy decreases and matrix inversion calculations are required.

[0005] ADI (Alternating Direction Implicit, ADI) can overcome the limitation of CFL (Courant FriedrichsLewy, CFL) number, so the time step selection does not depend on the Courant stability condition and is determined by the accuracy. The step size can be doubled and the calculation time can be greatly reduced. The alternating direction implicit time-domain difference (ADI-FDTD) method, which appeared in 1999, is artificially unconditionally stable. The traditional FDTD method uses Cartesian grids, which has step errors that affect the accuracy. FVTD using body-fitting curve coordinates does not have this problem. However, the implicit time-domain finite volume electromagnetic field algorithm of alternating direction (ADI) that can improve the efficiency of time advancement is not well developed at present. For this reason, the present invention proposes a new time-advancing finite volume algorithm based on the local orthogonal coordinate system, which is dimension-by-dimension (two-dimensional and three-dimensional coupling is converted into a series of one-dimensional calculations) to fill this deficiency.

[0006] In summary, it is necessary to develop an efficient direct solution to the Maxwell equations in the time domain, which can accurately simulate complex-shaped targets and adopt an implicit time-marching algorithm to relax the restrictions of very small grid scales on iterative physical time steps while maintaining high numerical accuracy, thereby improving computational performance and efficiency. Summary of the invention

[0007] The purpose of the present invention is to overcome the shortcomings of the prior art and provide a time-varying electromagnetic field implicit finite volume calculation method, equipment, medium and system, so as to realize an efficient time-domain electromagnetic numerical method for solving large-scale electromagnetic scattering problems of complex-shaped, high-frequency, and electrically large-sized conductor targets, and to improve the calculation efficiency while ensuring high precision of the time and space calculation formats.

[0008] The object of the present invention is achieved through the following solutions:

[0009] An implicit finite volume calculation method for time-varying electromagnetic fields, comprising:

[0010] The time-varying electromagnetic field of Maxwell's equations is iteratively solved using an implicit dual-time-step method based on time iterative advancement of the conservative electromagnetic field and spatial flux residuals; the outer layer of the simulation model is a physical time step loop until the calculation converges; the inner layer of the simulation model is a virtual time step sub-iteration loop until the sub-iteration converges; in each virtual time sub-iteration process, the spatial flux calculation and implicit iterative solution calculation are performed on each grid block and each grid unit in turn, and the conservative electromagnetic field value of the next level of virtual time sub-iteration step is updated.

[0011] Furthermore, the time-varying electromagnetic field of the Maxwell equations is iteratively solved in an implicit dual-time-step method based on time iterative advancement of the conserved electromagnetic field and spatial flux residual; the outer layer of the simulation model is a physical time step loop until the calculation converges; the inner layer of the simulation model is a virtual time step sub-iteration loop until the sub-iteration converges; in each virtual time sub-iteration process, the spatial flux calculation and the implicit iterative solution calculation are performed on each grid block and each grid unit in turn, and the conserved electromagnetic field value of the next virtual time sub-iteration step is updated, which specifically includes the sub-steps:

[0012] Among them, the process of physical time step cycle and virtual time step iteration cycle is:

[0013] Add virtual time derivative , the passive time-domain Maxwell equations to be solved in the form of scattered field are modified to:

[0014] ; (1)

[0015] ; (2)

[0016] in, is the conservation variable of the scattered electromagnetic field, It's virtual time. is the normalized physical time, , , is the electromagnetic flux in rectangular coordinates , , Quantity, is the real scattered magnetic induction intensity vector, is the real-valued electric displacement vector of the scattered field, including the subscript , , The scalars are the rectangular coordinates of the corresponding vectors , , Quantity;

[0017] when This system of equations approaches the original system of equations, with constant virtual time The subiterations are:

[0018] ; (3)

[0019] in, yes The value of the electromagnetic field conservation variable after virtual time iteration, yes The superscript * represents the value of the conserved variable corresponding to the virtual time sub-iteration. is the flux residual, is the number of physical time steps; is the physical time step, It is The electromagnetic conservation variables at the physical time step, It is Electromagnetic conservation variables at the physical time step; for The corresponding intermediate state flux residual; the physical time derivative uses the backward second-order difference;

[0020] The constant subiteration (3) uses a fully implicit algorithm and performs Virtual time step flux linearization expansion:

[0021] ; (4)

[0022] ;

[0023] ;

[0024] in, is the virtual time step; , , Corresponding to the curvilinear coordinate system Directional difference operator; , , Corresponding to the curvilinear coordinate system Directional electromagnetic flux Jacobian coefficient matrix; is the difference of conservation variables of scattering fields in adjacent physical times; It is Physical time step spatial flux residual, It is the total spatial flux residual of the previous physical time step plus the double time step correction; It is Grid unit Electromagnetic conservation variables for virtual time steps; It is Grid unit Electromagnetic conservation variables for virtual time steps; , , They are the corresponding curvilinear coordinate systems Directional electromagnetic flux;

[0025] The calculation process of spatial flux is as follows:

[0026] Steger-Warming splitting is used to calculate the interface flux of grid cells;

[0027] The fully implicit time iteration process of the local orthogonal system alternating direction is as follows:

[0028] ; (5)

[0029] in, is the identity matrix, is the time proportional coefficient, and , superscript is the virtual time step number, rhs is the spatial flux residual within the virtual time step;

[0030] The approximation factor AF is decomposed as follows:

[0031] ; (6)

[0032] Where h is the virtual time step in the alternating direction, ;

[0033] A local orthogonal coordinate system XYZ is established in each grid integration unit. The three directions of the local coordinate system are expressed in vector form. , , Represented by, thus transforming the global xyz rectangular coordinate system into a local orthogonal system XYZ:

[0034] ;

[0035] ;

[0036] ;

[0037] ;

[0038] ;

[0039] in, , , Is a vector The x, y, and z components of , , Is a vector The x, y, and z components of , , Is a vector The x, y, and z components of

[0040] In the passive case, Maxwell's equations are expressed in this local coordinate system as:

[0041] ; (7)

[0042] in, , , They are The x, y, and z components of

[0043] The Jacobian coefficient matrix corresponding to the flux in the local coordinate system is:

[0044] ;

[0045] The local coordinate system corresponds to the similarity matrix:

[0046] ; (8)

[0047] ; (9)

[0048] in, is the magnetic susceptibility, is the dielectric constant;

[0049] For formula (6) to calculate the alternating direction ADI, when ,Right now , the volume integral of the differential governing equations for each grid cell is obtained by applying Gauss's theorem:

[0050] ; (10)

[0051] in, is the local grid cell volume, is the local time step, is the physical time step, Corresponding curvilinear coordinate system Directional difference operator, is the mesh interface area in the corresponding direction, is the local electromagnetic flux Jacobian coefficient matrix, is the difference between adjacent pseudo time steps, is the total spatial flux residual of the previous physical time step plus the double time step correction;

[0052] Equation (10) is solved using the scalar pursuit method for tridiagonal strip equations, yielding:

[0053] ;

[0054] ;

[0055] Among them, the subscript Represents three adjacent units, parameters , yes Cell flux residuals, eigenvalues It is a constant in one dimension, so the catch-up method is used to calculate , and so on, and then calculate dimension by dimension, and finally calculate the difference of conservation variables in adjacent pseudo time steps .

[0056] Furthermore, the criterion for the sub-iteration convergence is the absolute value of the maximum amplitude difference between the electric field and the magnetic field in the calculation space, and a maximum number of sub-iteration steps is set in the numerical calculation.

[0057] A time-varying electromagnetic field implicit finite volume calculation device comprises a processor and a memory, wherein a computer program is stored in the memory, and when the computer program is loaded by the processor, the method described in any one of the above items is executed.

[0058] A computer-readable storage medium stores a computer program, wherein the computer program is loaded by a processor and executes any of the above methods.

[0059] A time-varying electromagnetic field implicit finite volume calculation system comprises the time-varying electromagnetic field implicit finite volume calculation device as described above.

[0060] The beneficial effects of the present invention include:

[0061] 1. The present invention provides a time-domain finite volume method with dual (physical and virtual) time iterative advancement and spatial flux implicitness of electromagnetic fields. Unlike explicit methods, the time step is limited and it can solve 2D and 3D Maxwell equations.

[0062] 2. The dual time-step time advancement and spatial flux implicit iteration and coefficient matrix splitting proposed in the present invention are iteratively solved in the local orthogonal coordinate system in alternating directions. The time-domain finite volume method combining dual time-step and spatial flux implicit can obtain the time and space distribution of electromagnetic fields stably, accurately and efficiently, and is suitable for numerical calculation of wave propagation, reflection and other problems.

[0063] 3. The present invention adopts implicit iteration of spatial flux and local orthogonal coordinate system coordinate transformation of Jacobi coefficient matrix and flux residual iterative solution, and replaces sparse matrix inversion with dimension-by-diagonal and tridiagonal pursuit method, which is simple and easy to use in engineering.

[0064] 4. The present invention provides a solution for iteratively calculating conservation variables of time-varying electromagnetic fields using the fully implicit finite volume method in alternating directions: a dual time-stepping formula combining virtual time with physical time for the Maxwell equations is used to provide a dimension-by-dimensional local orthogonal coordinate system iterative matrix implicit spatial flux residual inversion calculation formula and steps, thereby forming a new fully implicit finite volume software for time-varying electromagnetic fields in alternating directions. Compared with the traditional time-domain finite volume method for explicit calculation of time and flux, it can relax the extreme restrictions of the grid and explicit algorithm on the physical time step, improve the time advancement efficiency while ensuring the format accuracy, and achieve the purpose of improving the computing performance. BRIEF DESCRIPTION OF THE DRAWINGS

[0065] In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the drawings required for use in the embodiments or the description of the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying creative labor.

[0066] Figure 1 It is a software flow chart for realizing the method of the present invention;

[0067] Figure 2It is a comparison of the time history of the transverse magnetic field at the same position of the cylindrical (ka=5) TM wave; wherein the horizontal axis Time represents the physical time; the vertical axis Amp represents the amplitude of the scattered electromagnetic field; in the icon, By represents the y component amplitude of the real scattered magnetic induction intensity vector, 4RK represents the 4th order Runge-Kutta method, and DD_ADI_FVTD represents the time domain finite volume method using the local orthogonal system and alternating direction proposed by the present invention;

[0068] Figure 3 is a schematic diagram of the effect of the TM wave physical time step on the cylindrical (ka=5) radar cross section; represents the scattering azimuth, and the deg. in brackets indicates the unit is degree; the ordinate It represents the two-dimensional radar scattering width normalized by wavelength, and dB in brackets represents decibels; in the icon, Exact TM represents the exact solution of 2D TM wave polarization; Computation represents numerical calculation, dt represents the physical time step, and DD_ADI_FVTD represents the time-domain finite volume method of the local orthogonal system and alternating direction proposed by the present invention;

[0069] Figure 4 It is a metal football (f=1.18GHz) computational grid;

[0070] Figure 5 is a metal football (f=1.18GHz) vertically polarized, and the bistatic radar cross section calculated by the method of the present invention is compared with the result of the moment method. represents the scattering azimuth in degrees; the ordinate rcs represents the radar scattering cross section in decibels normalized per square meter; MOM represents the method of moments; DD_ADI_FVTD represents the time-domain finite volume method in the local orthogonal system and alternating direction proposed by the present invention. DETAILED DESCRIPTION

[0071] All features disclosed in all embodiments in this specification, or steps in all methods or processes implicitly disclosed, except for mutually exclusive features and / or steps, can be combined and / or expanded or replaced in any manner.

[0072] The specific implementation process of the present invention is as follows:

[0073] In a preferred embodiment, the present invention provides a fully implicit finite volume calculation method for time-varying electromagnetic field alternating direction, involving a body-fitting conformal space grid, dual time step time advancement, fully implicit electromagnetic flux residual, coordinate transformation based on local orthogonal coordinate system, and alternating direction dimension-by-dimension iterative calculation of electromagnetic numerical methods, especially how to improve the efficiency of solving Maxwell equations in the form of scattered fields, maintain the accuracy of time-varying electromagnetic field calculation, dual time step and main diagonal dominance of spatial flux residual, and fully implicit calculation method of ADI advancement in alternating direction of local orthogonal coordinate system. Specifically, firstly, unlike ADI-FDTD which only uses physical time derivative terms, the dual time method of the present invention introduces a steady virtual time derivative term in the control equation, and linearizes the physical time derivative in each physical time period, so that the physical time advancement step can be selected according to the physical problem without being restricted by stability, and the stability requirement of the calculation is satisfied by the virtual time sub-iteration, and the steady-state calculation of the sub-iteration can use the local time step to accelerate convergence, greatly shortening the calculation time required to obtain a stable electromagnetic field. Secondly, the implicit algorithm usually corresponds to the inverse calculation of a linear algebraic equation group. The present invention adopts the Maxwell equation group in the form of a scattered field, the approximate factorization AF (Approximate Factorization, AF) and linearization of the Jacobi coefficient matrix splitting, and obtains the pursuit method ADI that alternates in different directions. Each dimension is approximately implicitly decomposed, the implicit flux is transformed in the local orthogonal coordinate system and the inverse matrix operation is iteratively calculated in the local orthogonal coordinate system, and then the adjacent step difference of the virtual time step conservative electromagnetic field is iteratively calculated. After the virtual time step iteration converges, the new physical time step conservative electromagnetic field value is obtained.

[0074] Further, refer to Figure 1 The entire software for calculating electromagnetic fields using the finite volume method in the time domain can be divided into three parts according to its structure: preprocessing, electromagnetic field calculation, and postprocessing. Preprocessing mainly includes three modules: grid data input, calculation parameter data input, and control parameter input. It is mainly used to read in grid data, calculation parameter data input, and control parameter files, and on this basis, preprocessing is performed to provide calculation support for electromagnetic field calculation; electromagnetic field calculation includes: spatial electromagnetic field MUSCL format interpolation, unit interface flux calculation, time advancement, and convergence judgment module; postprocessing is mainly used to output the time and space distribution of the electromagnetic field, target surface induced current density, and radar scattering cross section output.

[0075] Next, we combine the two curl equations of the passive Maxwell equations to be numerically simulated. , ,in is the gradient symbol, is the real-numbered scattered field magnetic induction intensity vector, is the real-valued electric displacement vector of the scattered field, is the real-valued electric field intensity vector of the scattered field, is the real-valued scattered field magnetic field intensity vector, is the normalized physical time.

[0076] The dual-time-step time-domain finite volume method is introduced to iteratively calculate the target electromagnetic scattering process. The scattered field form of the time-domain Maxwell equations is expressed in rectangular coordinates as follows:

[0077]

[0078]

[0079]

[0080] in, is the conservation variable of the scattered electromagnetic field, It's virtual time. is the normalized physical time, , , is the electromagnetic flux in rectangular coordinates , , Quantity, is the real scattered magnetic induction intensity vector, is the real-valued electric displacement vector of the scattered field, is the electric field intensity vector of the scattered field, Scattering field magnetic field intensity vector, including subscripts , , The scalars are the rectangular coordinates of the corresponding vectors , , Quantity.

[0081] For objects with complex shapes, the computational space body-fitting multi-block structure grid is used, so there is a coordinate transformation:

[0082]

[0083] in Represents a curvilinear coordinate system Three directions, respectively One of them. The conservation form of Maxwell's equations in the curvilinear coordinate system for numerical simulation is obtained:

[0084]

[0085]

[0086]

[0087]

[0088]

[0089]

[0090]

[0091] Where V is the Jacobian matrix of the coordinate transformation, and the superscript variables represent the values ​​in the curvilinear coordinate system, which are obtained by the coordinate transformation. is the conserved variable in the curvilinear coordinate system, , , That is when Take the direction of the curve coordinate system respectively The electromagnetic flux in the curvilinear coordinate system at one time.

[0092] In order to get rid of the drawback of large computational complexity caused by the limitation of the global minimum physical time step of the traditional explicit time-domain finite volume, a fully implicit dual-time-step calculation method is invented to improve the time advancement efficiency of the target time-varying electromagnetic field. The dual time step is used to ensure the time accuracy while selecting the physical time step according to the physical characteristics. The implicit spatial flux residual is combined to obtain a stable and efficient calculation process, which includes the following steps:

[0093] Step 1: Based on the physical background of the electromagnetic problem being simulated, simulation modeling is performed in combination with boundary condition information.

[0094] Step 2: Use quadrilateral (two-dimensional) or hexahedral structure (three-dimensional) to mesh the simulation model. The mesh is encrypted on the wall and geometric singularity. The mesh gradually becomes sparse away from the scattering wall; numerically calculate the mesh of the corresponding area, output the mesh data file, set and output the boundary condition file. The mesh density ensures 15-20 grid points per wavelength, the wall density is >300 points / per wavelength, and the geometric singularity is encrypted to 50-100 grid points / per wavelength. The two-dimensional grid is advanced by 1 layer according to the right-hand rule in the plane perpendicular to the plane, and is calculated uniformly as a special case of the three-dimensional problem. The mesh data file includes the number of structural grid blocks and the dimensions of each block in the three curvilinear coordinate system.

[0095] Step 3: Preprocessing part, input the target calculation electromagnetic parameters and numerical calculation control parameters. The virtual time sub-iteration is implicit because the CFL number is not subject to the above explicit stability requirements. The physical time step is pre-set. For two-dimensional problems, the physical time step of the incident electromagnetic wave period is dimensionless. The calculation accuracy is good enough, and the three-dimensional problem is selected At the same time, set the sub-iteration convergence standard criterion value (for example, if the maximum amplitude difference of adjacent sub-iteration time steps in the full grid space is <0.001, it is judged to be converged) and set the maximum number of sub-iteration steps (for example: isubmax=25).

[0096] Step 4: Input the mesh data and boundary condition information file to initialize the calculation of the spatial electromagnetic field.

[0097] Step 5: The time-varying electromagnetic field of the dielectric Maxwell equations is iteratively solved using an implicit dual time-stepping method based on the time iterative advancement of the conserved electromagnetic field and the spatial flux residual.

[0098] Step 5-1: The outer physical time step cycle ends when the calculation converges.

[0099] Step 5-2: Inner virtual time step sub-iteration loop, until the sub-iteration converges or the maximum sub-iteration step number condition is met. The time accuracy of the dual-time method unsteady calculation is also restricted by the number of sub-iteration steps in each physical time step. The time accuracy can be guaranteed only when the number of iteration steps is large and the convergence control is strict. If the real time step is small, the sub-iteration can be obtained less. In numerical calculation, in order to prevent the field from changing violently and the number of sub-iteration CFL is too large, the residual is likely to fail to decrease to the convergence, and fall into an infinite loop, and the maximum number of sub-iteration steps is given. On the other hand, the sub-iteration convergence is judged and controlled, and the iteration is terminated as soon as possible under the premise of ensuring a certain time accuracy. The sub-iteration convergence criterion uses the absolute value of the maximum amplitude difference between the electric field and the magnetic field in the calculation space.

[0100] The specific algorithms of step 5-1 and step 5-2 are introduced below: MUSCL (Monotonic Upstream Schemes for Conservation Laws, MUSCL) interpolation is used.

[0101] Add virtual time derivative , the time-domain Maxwell equations to be solved in the form of scattered field are modified to:

[0102]

[0103] Obviously, when When converged, the system of equations is equivalent to the original system of equations, and the constant virtual time The sub-iteration is represented as:

[0104]

[0105]

[0106] in, yes The value of the electromagnetic field conservation variable after virtual time iteration, yes The superscript * represents the value of the conserved variable corresponding to the virtual time sub-iteration. is the flux residual, is the flux residual Added the residual after the physical time derivative term and source term; is the number of physical time steps; is the physical time step, It is The electromagnetic conservation variables at the physical time step, It is Electromagnetic conservation variables at the physical time step; for The corresponding intermediate state flux residual; the physical time derivative of the backward second-order difference is the second-order time accuracy; the constant sub-iteration part adopts the implicit algorithm:

[0107]

[0108] in, is the structural grid curve coordinate system direction 1, is the structural grid curve coordinate system direction 2, is the structural grid curve coordinate system direction 3; They correspond to the curvilinear coordinate system Directional electromagnetic flux; is an implicit control parameter, The fully implicit form of the other parameters corresponds to the explicit and implicit mixed format; subscript is the grid cell number, It is Grid unit Electromagnetic conservation variables for virtual time steps, It is Grid unit Electromagnetic conservation variables for virtual time steps, It is The electromagnetic conservation variables for the physical time step, It is The electromagnetic conservation variables for the physical time step, It is The electromagnetic conservation variables for the physical time step, It is Grid unit The spatial flux residual at the virtual time step, It is Grid unit The spatial flux residual at the virtual time step, It is a virtual time step controlled by stability, calculated by the CFL number and the local grid unit geometric scale and eigenvalue, which is significantly different from the explicit method. Different grid units in the steady-state calculation have different local virtual time sub-iteration steps, thereby accelerating the convergence of the electromagnetic field of the unit.

[0109] Implicit time-marching iteration equations:

[0110]

[0111]

[0112]

[0113]

[0114]

[0115]

[0116] in, is the virtual time step; , , Corresponding to the curvilinear coordinate system Directional difference operators; , , Corresponding to the curvilinear coordinate system Jacobian coefficient matrix of directional electromagnetic flux; is the difference of conservation variables of scattering fields in adjacent physical times; It is Physical time step spatial flux residual, is the total spatial flux residual of the previous physical time step plus the double time step correction.

[0117] The spatial flux is calculated as follows:

[0118] The Steger-Warming splitting method is used to calculate the interface flux of the grid cells.

[0119]

[0120]

[0121]

[0122] Subscript Take the curvilinear coordinate system One of the directions, corresponding That is the corresponding The direction of the electromagnetic flux, Represents a curvilinear coordinate system The electromagnetic flux obtained by splitting the positive eigenvalue in the Steger-Warming splitting in the corresponding direction; Represents a curvilinear coordinate system The electromagnetic flux obtained by splitting the negative eigenvalue in the Steger-Warming splitting in the corresponding direction; is the similarity matrix, are diagonal matrices consisting of positive and negative eigenvalues, respectively. They represent the left and right state variables at the interface respectively, and they can use the MUSCL format to achieve the highest third-order accuracy.

[0123]

[0124]

[0125] in is a limiter, subscript is the grid cell number, Corresponding unit interface, is the control parameter in third-order precision format, and are the post-difference and pre-difference operators respectively; Represented in grid cells The left state electromagnetic field conservation variable at the interface, Represented in grid cells Conservation variables of the right-state electromagnetic field at the interface; It is The conservation variables of the electromagnetic field scattered by each grid cell are: It is The grid cells scatter the electromagnetic field conservation variables.

[0126] The alternating direction implicit iterative solution calculation process is as follows:

[0127]

[0128] in, is the identity matrix, , superscript Number the virtual time step.

[0129] The approximation factor AF is decomposed as follows:

[0130]

[0131] in,

[0132] In order to simplify the characteristic matrix form, a local orthogonal coordinate system XYZ is established in each grid integration unit. The three directions of the local coordinate system are expressed in vector form. , , Represented by, thus transforming the global xyz rectangular coordinate system into a local orthogonal system XYZ:

[0133]

[0134]

[0135]

[0136]

[0137]

[0138] in , , Is a vector The x, y, and z components of , , Is a vector The x, y, and z components of , , Is a vector The x, y, and z components of .

[0139] In the passive case, Maxwell's equations are expressed in this local coordinate system as:

[0140]

[0141] in , , They are The x, y, and z components of .

[0142] The Jacobian coefficient matrix corresponding to the flux in the local coordinate system is:

[0143]

[0144] The local coordinate system corresponds to the similarity matrix:

[0145]

[0146]

[0147] in is the magnetic susceptibility, is the dielectric constant.

[0148] For formula (6), the alternating direction (ADI) calculation is performed, for example, when ,Right now , the volume integral of the differential governing equations for each grid cell is obtained by applying Gauss's theorem:

[0149]

[0150] in, is the local grid cell volume, is the local time step, is the physical time step, Corresponding curvilinear coordinate system Directional difference operator, is the mesh interface area in that direction, is the local electromagnetic flux Jacobian coefficient matrix, is the difference between adjacent pseudo time steps, is the total spatial flux residual of the previous physical time step plus the double time step correction.

[0151] Equation (10) is solved using the scalar pursuit method for tridiagonal strip equations, yielding:

[0152]

[0153]

[0154] The subscript Represents three adjacent units, parameters , yes Element flux residuals, eigenvalues It is a constant in one dimension, so the catch-up method is used to calculate , and so on, and then calculate dimension by dimension, and finally calculate the difference of conservation variables in adjacent pseudo time steps .

[0155] The above is the implicit dual-time-step, alternating-direction calculation of the Maxwell equations to control the iterative process of the time-varying electromagnetic field.

[0156] Step 6: Convergence judgment, post-processing process, output of time and space distribution of electromagnetic field, output of surface induced current and radar cross section spatial distribution data, etc.

[0157] like Figure 2As shown, the calculation conditions are TM polarization, electromagnetic scattering of metal cylinder, ka=5, a is the radius of the cylinder, and the same spatial point is calculated by the explicit 4th-order Runge-Kutta method and the dual-time-step DD_ADI_FVTD method respectively. The time oscillation history of the scattered electromagnetic field quantity is almost the same between the two, indicating that the accuracy of the physical time backward difference and the 2nd-order accurate dual-time DD_ADI_FVTD method meets the requirements of Fourier transform. Moreover, the same convergence standard can be seen in this example. The implicit ADI takes 4 cycles, which converges one cycle earlier than the explicit 4th-order Runge-Kutta method, which saves calculation time and improves calculation efficiency accordingly.

[0158] like Figure 3 As shown in the figure, the influence of different physical time steps of TM waves on the RCS results of a cylinder (ka=5) is compared. k is the wave number, a is the radius of the cylinder, and the calculation grid is 91x46. The object surface grid selects 18 grid points per wavelength, the far field boundary is outside 3 wavelengths, and the radial grid wall is encrypted, with an average of about 15 grid points per wavelength. As the physical time interval increases, the difference between the calculation results and the reference value gradually increases. Under the condition of too large physical time step (dt=0.1), the calculation accuracy of the dual time step method becomes lower. For the two-dimensional problem of the cylinder, the DD_ADI_FVTD software with alternating direction implicit iteration is used. Figure 3 It can be seen that when the physical time step is taken to the order of 0.01, the numerical calculation results are in good agreement with the analytical series solution.

[0159] like Figure 4 As shown in the figure, the computational grid for the electromagnetic scattering calculation of the metal football has a dimension of 40x46x65, and CFL=1.2 is selected. Fully implicit format, the far-field boundary is outside 3 wavelengths, the radial grid wall and the cusp geometric singularity are encrypted, the calculation frequency is f=1.18GHz, and the convergence criterion is that the maximum amplitude difference of the electromagnetic field in adjacent iteration steps in the whole space is <0.008.

[0160] Figure 5 In this example, the dual-station radar scattering cross section DD_ADI_FVTD is compared with the moment method (MOM) in the vertical polarization case at the physical time step dt=0.001. It can be seen from the two polarization calculation results that even when the radar scattering cross section is a very small amount -40dBsm, the two are in good agreement, indicating that for typical three-dimensional problems, the physical time step should be taken to the order of 0.001, and DD_ADI_FVTD can guarantee sufficient calculation accuracy.

[0161] Table 1 compares the time steps required for calculating the above examples using explicit 4th-order Runge-Kutta FVTD and implicit alternating direction DD_ADI_FVTD based on the local orthogonal coordinate system. It can be seen that due to the mesh encryption and the explicit method requiring the entire calculation space to use the minimum step size that meets the stability requirements, the explicit method has a smaller step size, which makes the calculation workload greater.

[0162] Table 1

[0163]

[0164] These numerical examples show that the time-varying electromagnetic field corresponding to the invention, based on the local orthogonal coordinate system and the alternating direction implicit finite volume algorithm DD_ADI_FVTD, can not only ensure the numerical accuracy of the format while relaxing the time step restriction, accurately obtain the time and spatial distribution of the target electromagnetic scattering field, but also improve the time advancement efficiency and save computing costs.

[0165] The units involved in the embodiments of the present invention may be implemented by software or hardware, and the units described may also be arranged in a processor. The names of these units do not, in some cases, limit the units themselves.

[0166] According to one aspect of an embodiment of the present invention, a computer program product or a computer program is provided, the computer program product or the computer program includes a computer instruction, and the computer instruction is stored in a computer-readable storage medium. A processor of a computer device reads the computer instruction from the computer-readable storage medium, and the processor executes the computer instruction, so that the computer device executes the method provided in the above various optional implementations.

[0167] As another aspect, an embodiment of the present invention further provides a computer-readable medium, which may be included in the electronic device described in the above embodiment; or may exist independently without being assembled into the electronic device. The above computer-readable medium carries one or more programs, and when the above one or more programs are executed by an electronic device, the electronic device implements the method described in the above embodiment.

Claims

1. An implicit finite volume calculation method for time-varying electromagnetic fields, characterized in that: include: The time-varying electromagnetic field of the Maxwell equations is iteratively solved by an implicit dual-time-step method based on the time iterative advancement of the conserved electromagnetic field and the spatial flux residual; the outer layer of the simulation model is a physical time step cycle until the calculation converges; the inner layer of the simulation model is a virtual time step sub-iteration cycle until the sub-iteration converges; in each virtual time sub-iteration process, the spatial flux calculation and implicit iterative solution calculation are performed on each grid block and each grid unit in turn, and the conserved electromagnetic field value of the next virtual time sub-iteration step is updated; The method uses an implicit dual time step method based on time iterative advancement of the conserved electromagnetic field and spatial flux residual to iteratively solve the time-varying electromagnetic field of the Maxwell equations; the outer layer of the simulation model is a physical time step loop until the calculation converges; the inner layer of the simulation model is a virtual time step sub-iteration loop until the sub-iteration converges; in each virtual time sub-iteration process, the spatial flux calculation and the implicit iterative solution calculation are performed on each grid block and each grid unit in turn, and the conserved electromagnetic field value of the next virtual time sub-iteration step is updated, which specifically includes the following sub-steps: The process of physical time step loop and virtual time step iteration loop is: Add virtual time derivative , the passive time-domain Maxwell equations to be solved in the form of scattered field are modified to: ；(1) ；(2) in, is the conservation variable of the scattered electromagnetic field, It's virtual time. is the normalized physical time, , , is the electromagnetic flux in rectangular coordinates , , Quantity, is the real scattered magnetic induction intensity vector, is the real-valued electric displacement vector of the scattered field, including the subscript , , The scalars are the rectangular coordinates of the corresponding vectors , , Quantity; when This system of equations approaches the original system of equations, with constant virtual time The subiterations are: ；(3) in, yes The value of the electromagnetic field conservation variable after virtual time iteration, yes The superscript * represents the value of the conserved variable corresponding to the virtual time sub-iteration. is the flux residual, is the number of physical time steps; is the physical time step, It is The electromagnetic conservation variables at the physical time step, It is Electromagnetic conservation variables at the physical time step; for The corresponding intermediate state flux residual; the physical time derivative uses the backward second-order difference; The constant subiteration (3) uses a fully implicit algorithm and performs Virtual time step flux linearization expansion: ；(4) ； ； in, is the virtual time step; , , Corresponding to the curvilinear coordinate system Directional difference operator; , , Corresponding to the curvilinear coordinate system Directional electromagnetic flux Jacobian coefficient matrix; is the difference of conservation variables of scattering fields in adjacent physical times; It is Physical time step spatial flux residual, It is the total spatial flux residual of the previous physical time step plus the double time step correction; It is Grid unit Electromagnetic conservation variables for virtual time steps; It is Grid unit Electromagnetic conservation variables for virtual time steps; , , They are the corresponding curvilinear coordinate systems Directional electromagnetic flux; The calculation process of spatial flux is as follows: Steger-Warming splitting is used to calculate the interface flux of grid cells; The fully implicit time iteration process of the local orthogonal system alternating direction is as follows: ；(5) in, is the identity matrix, is the time proportional coefficient, and , superscript is the virtual time step number, rhs is the spatial flux residual within the virtual time step; The approximation factor AF is decomposed as follows: ；(6) Where h is the virtual time step in the alternating direction, ; A local orthogonal coordinate system XYZ is established in each grid integration unit. The three directions of the local coordinate system are expressed in vector form. , , Represented by, thus transforming the global xyz rectangular coordinate system into a local orthogonal system XYZ: ； ； ； ； ； in, , , Is a vector The x, y, and z components of , , Is a vector The x, y, and z components of , , Is a vector The x, y, and z components of In the passive case, Maxwell's equations are expressed in this local coordinate system as: ；(7) in, , , They are The x, y, and z components of The Jacobian coefficient matrix corresponding to the flux in the local coordinate system is: ； The local coordinate system corresponds to the similarity matrix: ；(8) ；(9) in, is the magnetic susceptibility, is the dielectric constant; For formula (6) to calculate the alternating direction ADI, when ,Right now , the volume integral of the differential governing equations for each grid cell is obtained by applying Gauss's theorem: ；(10) in, is the local grid cell volume, is the local time step, is the physical time step, Corresponding curvilinear coordinate system Directional difference operator, is the mesh interface area in the corresponding direction, is the local electromagnetic flux Jacobian matrix, is the difference between adjacent pseudo time steps, is the total spatial flux residual of the previous physical time step plus the double time step correction; Equation (10) is solved using the scalar pursuit method for tridiagonal strip equations, yielding: ； ； Among them, the subscript Represents three adjacent units, parameters , yes Cell flux residuals, eigenvalues It is a constant in one dimension, so the catch-up method is used to calculate , and so on, and then calculate dimension by dimension, and finally calculate the difference of conservation variables in adjacent pseudo time steps .

2. The implicit finite volume calculation method for time-varying electromagnetic fields according to claim 1 is characterized in that: The criterion for the sub-iteration convergence is the absolute value of the maximum amplitude difference between the electric field and the magnetic field in the calculation space, and a maximum number of sub-iteration steps is set in the numerical calculation.

3. A time-varying electromagnetic field implicit finite volume calculation device, characterized in that: The method comprises a processor and a memory, wherein a computer program is stored in the memory, and when the computer program is loaded by the processor, the method according to any one of claims 1 to 2 is executed.

4. A computer-readable storage medium, characterized in that: A computer program is stored in a readable storage medium, and the computer program is loaded by a processor to execute the method according to any one of claims 1 to 2.

5. A time-varying electromagnetic field implicit finite volume calculation system, characterized in that: The device comprises the time-varying electromagnetic field implicit finite volume calculation device as described in claim 3.

Images