An Asymmetrical Mixed Higher-Order Discontinuous Galerkin Time Domain Method for Electromagnetic Scattering from the Plasma Sheath around a Hypersonic Vehicle

,


Introduction
In recent years, hypersonic weapons near space have played an increasingly important strategic role in warfare.When a hypersonic vehicle reenters the atmosphere, ablative hypersonic viscous fows produces a plasma sheath [1,2], which is a nonuniform plasma fow composed of free electrons, ions, and neutral particles around hypersonic vehicles [3,4].
In early research, an approximate method was used to predict the plasma efect, assuming that the plasma fow had a certain electron density distribution.In [5], Rusch simulated the actual plasma sheath profle of the reentry vehicle using a parabolic electron density profle.In [6], a uniform plasma layer was set on the surface of the spherical model to calculate the interaction between the electromagnetic wave and the target.Later, with the development of computer technology, numerical calculation methods such as fnite-diference time-domain (FDTD), integral equation method (IEM), and physical optics (PO) were widely used in the study of electromagnetic wave in plasma fow.Basically, the simulations of interactions between the plasma and electromagnetic wave are performed by using the particle-in-cell method, which is based on the fnite diference time domain (FDTD) method [7,8] or the FDTD method combined with the efective permittivity model.
In [9], the physical optics method was used to calculate the scattering characteristics of the electromagnetic wave in inhomogeneous plasma sheath from the S-band to the Ku-band frequencies.In [10], a well-conditioned internally combined volume surface integral equation (ICVSIE) for analyzing electromagnetic scattering from perfect electrically conducting (PEC) surfaces coated with negative permittivity plasmas is presented to be applied in various canonical scatterers and a plasma-engulfed reentry vehicle.Scarabosio et al. employ the Eikonal approximation in the large inhomogeneous plasma region [11] and compute radiation and scattering via the equivalence theorem, and the results show that signifcant radio link path losses can be associated with plasma spatial variations (gradients) and collisional losses.In [12], the Drude mode is proven to be adopted to characterize collisional plasma.In the calculation with FDTD, radar cross section (RCS) of nonuniform plasma spheres and conductor spheres coated by the plasma layer are calculated.A piecewise linear JE recursive convolution fnite-diference time-domain (PLJERC-FDTD) was used to study the scattering characteristics of electromagnetic waves in time-varying and inhomogeneous plasma sheath [13].In [14], the high-order FDTD method was used to simulate the RCS (radar cross section) from the plasma fow in a shock tube.To increase the computational efciency, Chen and Wang et al. used the hybrid-implicit-explicit (HIE)-FDTD method [15,16] to simulate the gyrotropic plasma in open regions [17].Chen and Wang proposed to apply the weakly conditionally stable (WCS)-FDTD method [18] to simulate the plasma frequency selective surface [19].Numerical experiments showed that both the HIE-FDTD and WCS-FDTD methods can enhance signifcantly the computational efciency of plasma simulations.However, because of the used Yee cell, the FDTD methods sufer from serious accuracy deterioration when they come to irregular objects.
In [20], electromagnetic wave behavior in plasma with the combined CFD and FD2TD method was investigated.Te distribution of charged particles around the ESA ARD and the complicated behavior of electromagnetic waves, with attenuation and refection, are clarifed in detail.
For DGTD analysis, elements in tetrahedral form can be better ft to irregular objects [21,22].In [23], DGTD is used to solve electromagnetic scattering from hypersonic aircraft with plasma sheath; however, the fight status and the specifc parameters of the plasma sheath were not given in the paper.In [24], an exponential-based DGTD method for the 3-D modeling of hypersonic vehicles with surrounded by a uniform plasma layer is presented, and this method enables the usage of larger time steps compared with the explicit time marching scheme.In the abovementioned DGTD methods, tetrahedral elements are used in the computational domain, which is inefcient when calculating complex geometric objects.More recently, the DGTD with hexahedron have attracted extensive attention for modeling with complex geometrical features [25].In [26], DGTD with the high-order basis function is proposed to solve electromagnetic scattering from hypersonic aircraft with plasma sheath in which dielectric parameter changes with frequency, respectively.A hybrid mesh DGTD algorithm based on virtual elements for tetrahedra and hexahedra is proposed to improve the computational efciency [27].But tetrahedra and hexahedra are applied in diferent subdomains, respectively.
In this study, the Lagrange high-order element is used to read more plasma parameters in the plasma fow region, the efcient Serendicity high-order element is used in the uniform medium area, and the asymmetric high-order elements are constructed as the transition units between two types of elements to calculate the numerical fux.Te highorder points of the Serendicity high-order element close to the shape of aircraft can be mapped to aircraft surface by isoparametric element.Tis article is organized as follows.In Section 2, time-domain Maxwell's equations in the Drude dispersive medium are presented.Section 3 describes the physical and numerical models of asymmetric mixed higherorder mesh.As far as the authors' knowledge, this is the frst time that an asymmetric mixed higher-order mesh DGTD method is developed for the numerical simulation of nonuniform plasma fow target.Finally, numerical experiments of scattering problems are presented in Section 4, and the accuracy of the proposed method is verifed through comparison with the results of Mie theory.It was confrmed that plasma fow under diferent fight altitudes has diferent efects on electromagnetic scattering.

Iterative Equation of DGTD
We recall the propagation of an electromagnetic wave in the dispersive medium.We can use the dielectric coefcient to express the characteristics of the nonmagnetized medium.Te dielectric coefcient related to frequency can be given by where ε ∞ is the relative dielectric constant at infnity frequency, ε 0 is the dielectric coefcient in vacuum, and χ ω is the polarizability function.

Drude Dispersive Medium.
Te nonuniform plasma fow of hypersonic vehicles is usually treated as cold plasma, which is expressed by the Drude dispersive medium.Te thermochemical nonequilibrium model is used in computational fuid dynamics (CFD) to simulate the nonuniform plasma fow.Te concentration of electrons, density of neutral particles, and temperature in the nonuniform plasma fow are important parameters for calculating the dielectric coefcient.In the Drude dispersive medium [28], χ ω is given by where ω p is the plasma collision frequency and ] c is the plasma frequency which are given by where m e is the electron mass, n e is the number density of electrons, n s is the number density of species s, e is the magnitude of the electronic charge equal to 4.80298 × 10 − 10 esu, and k is Boltzmann's constant.Te efective electron-neutral energy exchange cross section is defned by a curve ft of the form.
Te constants for equation ( 4) are presented in Table 1.
2 International Journal of Antennas and Propagation

Formulation in Computational Domain. Te Maxwell equations are
( Substituting equation (1) into equation ( 5), we obtain Space Ω is divided into nonoverlapping and continuous subregions.Each subregion is a separate hexahedral element, and its weight function (also called basis function) is ϕ i (i is the serial number of the unit).Te weighted integral of equation ( 6) in each element is taken and set equal to zero, such that where ϕ i is the weight function.Te numerical fuxes are defned at the boundary between two elements as where  n is the unit outward normal vector in each surface of the element, E * and H * are the numerical fuxes, E + and H + are in the felds of the adjacent element, and E and H are in the felds of the element.κ e , κ h , υ e , and υ h are presented in Table 2.
Substituting ( 9) into (7), we obtain International Journal of Antennas and Propagation

Asymmetrically Mixed Higher-Order Mesh
Te aircraft is a RAM-CII aircraft scaling model in this paper.Te nonuniform fow is simulated by FASTRAN whose numerical algorithm is the fnite volume method based on density [29].
In Figure 1, the plasma fow forms a viscous boundary layer on the aircraft surface when fight altitude is 55 km and velocity altitude is 20 Mach, and the plasma parameters have a large gradient.Te maximum value of plasma frequency ω p appears in the area of the blunt cone head (with an order of magnitude of 2π × 10 11 rad/s ), and the closer to the blunt cone head, the greater the value of plasma frequency ω p (as shown in Figure 1(a).Te plasma frequency in the other regions decreased approximately one order of magnitude, about 2π × 10 10 rad/s.Te maximum value of the collision frequency of υ c appears in the area of the blunt cone head, with an order of magnitude of 2π × 10 10 rad/s (as shown in Figure 1(b)).Te collision frequency υ c in the other regions decreased approximately one order of magnitude around 2π × 10 10 rad/s.
A detailed computational model was reconstructed with mixed high-order hexahedral grids for the DGTD code.In order to read more parameters of the nonuniform plasma fow, the electromagnetic model is divided by mixed highorder hexahedral grids.
A hexahedral element constructs the basis function through isoparametric elements [30,31].Te nodes of the actual mesh and isoparametric cell are in one-to-one correspondence, as shown in Figure 2. Construct three isoparametric unctions (as shown in Figures 2(a)-2(c)) for three actual unctions (as shown in Figures 2(d Te isoparametric element is the standard element in ξηζ coordinates.If there are n basis functions, they can be given by Te basis function ϕ i satisfes the abovementioned conditions, and the mapping formula from local coordinates to actual coordinates is as follows: Te isoparametric element can be mapped to element with arbitrary shape using the abovementioned formula.In the higher-order element, the right-angle edge of the standard element can be mapped to the parabolic boundary of the actual element.
In Figure 3, the 27 nodes high-order Lagrange element is used in the nonuniform plasma fow, which can read more plasma parameters by using more internal nodes.Te effcient 20 nodes Serendipity element is used in the homogeneous medium region, which has no internal nodes, thus ensuring the accuracy of calculation and reducing the amount of computation.In order to calculate the numerical fux between the two types of elements in the DGTD algorithm, an asymmetric high-order element is constructed as a transition unit.

NonUniform Plasma Flow Region Mesh.
Lagrange highorder element is used in the nonuniform plasma fow region.It is assumed that the isoparametric element is divided into r + 1 groups, p + 1 groups, and t + 1 groups of nodes in the x, y, and z directions, respectively, and the nodes are located at the intersection lines of the element.Node i is located in groups m, n, and k in the x, y, and z directions, respectively.Its interpolation functions are Te basis function N i of node i can be obtained by multiplying the Lagrange polynomials in the above: Grid-generation software usually provides only the frst order nodes of an element.In order to efectively ft the higher order element nodes near the boundary of the plasma fow to the contour of the plasma fow, we set more surface points of the plasma fow region as mapping points of the higher-order element by constructing a cylindrical coordinate system.
Firstly, the surface mesh of the fow region is further subdivided separately.Te further subdivided mesh is denser than the working mesh as shown in Figure 4(b).
Ten, the cylindrical coordinate system rθZ in Figure 4 is established in the further mesh.Te actual coordinate points and the higher-order grid nodes close to the fow feld have the same values on the Z and θ axes (as shown in Figure 5).Te transformation relationship between the cylindrical coordinate system and the rectangular coordinate system is as follows: Finally, a suitable cylindrical surface point is found according to the value of Z axis and θ axis corresponding to the higher-order node.

Uniform Medium Region Mesh.
Te uniform medium region is divided by the Serendipity high-order element.Te Serendipity higher-order element is shown in Figure 3(c), whose basis functions are 6 International Journal of Antennas and Propagation

Transition Unit (Asymmetric High Order Element).
Te Lagrange high-order element with 27 nodes and the Serendipity high-order element with 20 nodes are used in this paper.Tere are 9 nodes at each face of the Lagrange high-order element and 8 nodes at each face of the Serendipity high-order element.An asymmetric high-order element with 21 nodes is constructed as a transition element (as shown in Figure 6).
To determine the arbitrary basis function N i as follows, we divide n planes that do not pass through node i but pass through all other nodes of the element.
All basis functions can be deduced in terms of the abovementioned way.

Examples
In the following, two numerical experiments are used to verify the efectiveness of the asymmetrically mixed highorder mesh DGTD method, and the infuence of nonuniform plasma fow on electromagnetic scattering is discussed.
In the case of mixed higher order element, the plasma sphere mesh is mapped to 27 nodes second-order Lagrange element, and the other mesh is mapped to 20 nodes secondorder Serendipity element, an asymmetric high order element with 21 nodes is constructed as a transition unit between two types of elements.
As shown in Figure 7, a very good agreement is achieved among the results obtained using the mixed higher-order element method, second-order Lagrange element method, and the Mie theory.Compared with DGTD with the frst order element, this method has less RCS error, especially when the frequency is 76.2 GHz, 93.5 GHz (as shown in Table 3).

RCS of Inhomogeneous Plasma Flow around Blunt Cone.
In this example, monostatic and bistatic radar cross section (RCS) of the blunted cone surrounded by nonuniform plasma fow is illuminated.We consider a Gaussian pulse incidents from the front of the blunted cone (θ � 90, ϕ � 0).Te blunted cone has a blunted nose of radius 0.01 m, a conical cross section with half-angle of 9 ∘ , and a whole length 0.1 m.
As shown in Figure 8, the 3D geometric model is divided by the hexahedral grids consists of 1067631 hexahedrons.Table 4 shows the memory required and calculation times of 27 nodes Lagrange mesh and mix-order mesh with 8 parallel threads algorithm.Te results demonstrate that mix-order mesh is more efcient, with the required memory being reduced by 26.85%.To truncate the open boundary of the simulation, the uniaxial perfectly matched layer (UPML) with a thickness of 0.0032 m is placed at least 0.006 m away from the plasma sheath.
In Figure 9, the infuence of fight altitude and Mach number on the backward RCS of blunt cones at altitudes of 50 and 55 km and speeds of 20 Mach are analyzed.When the frequency is below 1.8 GHz, the plasma fow enhances the 8 International Journal of Antennas and Propagation    International Journal of Antennas and Propagation backward RCS of the blunt cone greatly.When the collision frequency is greater than the incident wave frequency, reducing the incident wave frequency will increase refection [33].When the frequency is 1.8 GHz∼13 GHz, the plasma fow will reduce the backward RCS of blunt cone.When the frequency is above 13 GHz, the backward RCS of blunt cone with plasma is consistent with that without plasma, indicating that the plasma fow has little infuence on the electromagnetic scattering of the reentry target in the high frequency band.
Aiming at abovementioned characteristics, we discuss the bistatic RCS and the electric feld amplitude of the nonuniform fow feld target, which is at the altitude of 55 km and the velocity of 20 Mach with three frequency bands of 1.0 GHz, 4.0 GHz, and 15.0 GHz.
Te amplitude at 1.0 GHz without and with the plasma fow of XOZ section is given in Figures 10(a) and 10(b).When the frequency is 1.0 GHz, the incident wave cannot penetrate the plasma fow near the nose as show in the partial enlarged view.In the plasma layer near the nose, the blue shadow in Figure 10(b) increases compared to Figure 10(a).Terefore, the backward and lateral RCS of the blunt cone with plasma fow, particularly in 50 km 20 Mach condition, are greater than that of the blunt cone without plasma fow in Figure 10(c).Te interaction between the plasma fow and the incident wave is more obvious, as the altitude decreases.
When the frequency of the incident wave increases to 4.0 GHz, the lateral RCS of the blunt cone with plasma fow are smaller than that of the blunt cone without plasma fow, particularly in 50 km 20 Mach condition, the backward RCS is the reverse of the lateral RCS (as shown Figure 11(c)).Compared with Figure 11(a), the efect of plasma fow on incident waves is relatively strong in Figure 11(b).
When the frequency of the incident wave is 15.0 GHz, the incident wave can penetrate the plasma fow.Te infuence of plasma fow on amplitude decreases in Figure 12(b).Te RCS of the blunt cone with plasma fow is close to that without plasma fow, particularly in 55 km 20 Mach condition (as shown in Figure 12(c)).12 International Journal of Antennas and Propagation

Conclusion
Most current calculation schemes for DGTD analysis are based on a single type of scheme, which is inefcient when calculating no-uniform fow feld.Here, we propose an asymmetrically mixed higher-order mesh DGTD method to address this problem.Tis method took into account the higher-order computational accuracy while reducing computation.Te reliability of the asymmetrically mixed highorder mesh DGTD method is verifed by plasma ball experiments.Compared with the single-type high-order mesh, the proposed method reduces the computational efort and has the advantage over the frst-order mesh to further profle the target.We study the infuence between plasma fow of hypersonic vehicle and electromagnetic wave from the perspective of amplitude and scattering for the frst time.More research ideas are proposed for studying the relationship between fight height and electromagnetic wave frequency.Tis article has considered an asymmetrically mixed higher-order mesh DGTD method for only two types of mesh, the Lagrange high-order mesh with 27 nodes and the Serendipity high-order mesh with 20 nodes, for a type of transition unit.Tis method currently applies only to these two kinds of mesh.Nevertheless, it can also be adapted for other types of Lagrange grids and Serendipity grids and constructs corresponding transition unit.Moreover, when the order of the magnitude gradient of nonuniform medium parameters is relatively large, diferent types of mesh need to be constructed according to the gradient values.Tis will be investigated in future work.

Figure 3 :
Figure 3: A transition unit between two types of elements (a) 27 nodes Lagrange element, (b) the asymmetric high-order element, and (c) 20 nodes serendipity element.

Figure 7 :
Figure 7: Backward scattering of the plasma sphere.

Table 1 :
Constants for curve fts of electron-neutral energy exchange cross section σ es .

Table 3 :
RCS error of mix order element and frst order element.

Table 4 :
Parameter records for both meshes.