CN103235193A - Numerical method of satellite electromagnetic scattering characteristics within millimeter wave band - Google Patents

Abstract

The invention discloses a numerical method of satellite electromagnetic scattering characteristics within a millimeter wave band. Aiming at such a metal and medium mixed structure of a satellite, an electric field integral equation is established only for the metal part of the satellite, an equation does not need to be established for the medium part, and finally, an equation with good conditions is formed, so as to be convenient for iterative solution. Mesh generation is performed only on the metal surface of the satellite, the mesh generation does not need to be performed on the medium coated part, and the strangeness is easily processed. Because of a structural feature of the inclined long satellite, the directionality of a transfer factor component is strong, an acceleration method is adopted for a far field effect on the basis of a multi-level fast multipole algorithm, computational memory requirements are effectively reduced, and the computing time is saved.

Description

The numerical method of satellite electromagnetic scattering characteristic in the millimere-wave band

Technical field

The invention belongs to the numerical method of the interior satellite electromagnetic scattering characteristic of quick computing technique, particularly millimere-wave band of Electromagnetic Scattering of Target characteristic.

Background technology

Satellite is general mainly to be made of metal material, in order to reach anticorrosion, decoration and purpose such as stealthy, can adopt dielectric material to apply (QJ2754-95) usually.For this electromagnetic scattering problems that is applied mixed structure by metal medium of satellite, scholars have proposed a series of electromagnetism numerical computation method, wherein based on method of moment (MoM) (the Ding-Feng Yu of integral equation, Si-Yuan He, Hai-Tao Chen, Guo-Qiang Zhu, " Research on the electromagnetic scattering of3D target coated with anisotropic medium using impedance boundary condition ", Microwave and Optical Technology Letters, 2011,53 (2): 458 – 462.) obtained using the most widely.

At this metal medium mixed structure, the metal part is used as ideal conducting body (PEC) usually and handles, and quilt cover integral Equation Methods (SIE) is come analysis and solution easily, yet for the medium part, no matter be to use the body integral Equation Methods, also be based on the face integral Equation Methods of equivalence principle and their differentiation algorithm, with respect to the electromagnetism computational analysis method of simple metal target, these algorithms all need more internal memory and time resource.And because satellite coating media part is extremely thin, if select for use area to divide equation to analyze the medium part, the singularity between the upper and lower surface is more difficult; If adopt volume to analyze the medium part, can cause the subdivision mesh distortion, finally cause the matrix equation condition poor, be difficult to find the solution.

Satellite in the millimeter wave band belongs to the TV university size objectives, even if adopted just speed-up computation of the quick multistage sub-technology of multilayer, still faces the problem of computational resource consume significant.In the quick multistage subalgorithm of multilayer, along with the increase of the number of plies, the calculated amount of transition matrix can sharply increase.Transition matrix between two far field groups is made up of ball Hankel function, Legnedre series, and transfer factor is defined within on the unit ball of a large amount of angular spectrum components, even adopt interpolation and anti-interpolation method, the calculating of whole transfer process also is very consuming time.

Summary of the invention

The object of the present invention is to provide the numerical method of satellite electromagnetic scattering characteristic in a kind of millimere-wave band, thereby realize obtaining fast the method for electromagnetic scattering characterisitic parameter.

The technical solution that realizes the object of the invention is: the numerical method of satellite electromagnetic scattering characteristic in a kind of millimere-wave band, and step is as follows:

The first step makes uniform plane wave shine on the satellite structure that scribbles film dielectric layer, will produce surface induction surface current J on the surface of satellite _S, surface charge ρ _S, produce body polarization electric current J in the thin-medium _Pol, body polarization magnetic current M _PolWith polarization charge ρ _{S, pol}, according to the electric field boundary condition of perfect conductor, namely the resultant field tangential component of metal surface is 0, obtains scribbling the field integral equation EFIE of the satellite structure target of film dielectric layer, and is as follows

[E ^Inc(r)+E ^Sca(r)] _Tan=0 (1) wherein, and subscript tan represents the tangential component of electric field, E ^IncExpression is radiated at the electromagnetic wave on the target, E ^ScaThe end of a performance that the expression target produces behind electromagnetic wave irradiation, the expression-form of scattered field is:

$E^{\mathrm{sca}}=-{\mathrm{j\ω\μ}}_0\underset{S}{\∫}{J}_S\left(r^{\′}\right)\stackrel{\‾}{G}(r,r^{\′}){\mathrm{dS}}^{\′}-{\mathrm{j\ω\μ}}_0\underset{V}{\∫}{J}_{\mathrm{pol}}\left(r^{\′}\right)\stackrel{\‾}{G}(r,r^{\′}){\mathrm{dV}}^{\′}+\underset{V}{\∫}{M}_V\left(r^{\′}\right)\×\▿G(r,r^{\′}){\mathrm{dV}}^{\′}$

$=-{\mathrm{j\ω\μ}}_0\underset{S}{\∫}{J}_S\left(r^{\′}\right)G(r,r^{\′}){\mathrm{dS}}^{\′}-\frac{1}{{\mathrm{\ϵ}}_0}\underset{S}{\∫}{\mathrm{\ρ}}_s\left(r^{\′}\right)G(r,r^{\′}){\mathrm{dS}}^{\′}$

(2)

$-{\mathrm{j\ω\μ}}_0\underset{V}{\∫}{J}_{\mathrm{pol}}\left(r^{\′}\right)G(r,r^{\′}){\mathrm{dS}}^{\′}-\frac{1}{{\mathrm{\ϵ}}_0}\underset{S}{\∫}{\mathrm{\ρ}}_{s,\mathrm{pol}}\left(r^{\′}\right)G(r,r^{\′}){\mathrm{dS}}^{\′}-\frac{1}{{\mathrm{\ϵ}}_0}\underset{S^{\mathrm{\Δ}}}{\∫}{\mathrm{\ρ}}_{s,\mathrm{pol}}\left(r^{\′}\right)G(r,r^{\′}){\mathrm{dS}}^{\mathrm{\Δ}\′}$

$+\underset{V}{\∫}{M}_V\left(r^{\′}\right)\×\▿G(r,r^{\′}){\mathrm{dV}}^{\′}$

Wherein, V represents the volume element of thin-medium coating, and S represents that the unit, metal surface is the lower surface unit of thin-medium coating, S ^ΔThe upper surface unit of expression thin-medium coating, ω is electromagnetic angular frequency, μ ₀And ε ₀Be respectively magnetic permeability and dielectric parameter in the vacuum, r and r' are respectively the position coordinates in field and source, G (r r') is the Green function of free space, and expression formula is:

$G(r,r^{\′})=\frac{e^{-\mathrm{jk}|r-r^{\′}|}}{4\mathrm{\π}|r-r^{\′}|}---\left(3\right)$

Wherein k is the wave number of free space;

In second step, convert 5 unknown quantitys in the field integral equation of setting up in the first step to the metal surface current density, J _S, 5 unknown quantitys are respectively the metal surface current density, J _S, metal surface electric density ρ _S, dielectric polarized current density J _Pol, dielectric layer upper and lower surface polarization charge densities ρ _{S, pol}With body polarization magnetic current M _Pol, its concrete representation is as follows,

Formula (1) finally becomes:

$[{\mathrm{j\ω\μ}}_0\underset{S}{\∫}{J}_s\left(r^{\′}\right)G(r,r^{\′}){\mathrm{dS}}^{\′}-\frac{1}{{\mathrm{j\ω\ϵ}}_0}\underset{S}{\∫}{\▿}^{\′}\·J_s\left(r^{\′}\right)\▿G(r,r^{\′}){\mathrm{dS}}^{\′}$

$-{\mathrm{j\ω\μ}}_0\mathrm{\κ}\underset{V}{\∫}{\▿}^{\′}\·J_s\left(r^{\′}\right)\hat{n}G(r,r^{\′}){\mathrm{dS}}^{\′}+\frac{\mathrm{\κ}}{{\mathrm{j\ω\ϵ}}_0}\underset{S}{\∫}{\▿}^{\′}\·J_s\left(r^{\′}\right)\▿G(r,r^{\′}){\mathrm{dS}}^{\′}-\frac{\mathrm{\κ}}{{\mathrm{j\ω\ϵ}}_0}\underset{S^{\mathrm{\Δ}}}{\∫}{\▿}^{\′}\·J_s\left(r^{\′}\right)\▿G(r,r^{\′}){\mathrm{dS}}^{\mathrm{\Δ}\′}$

$+{\mathrm{j\ω\μ}}_0({\mathrm{\μ}}_r-1)\underset{V}{\∫}{\hat{n}}^{\′}\×J_S\left(r^{\′}\right)\×\▿G(r,r^{\′}){\mathrm{dV}}^{\′}{]}_{\tan }=E_{\tan }^{\mathrm{inc}}---\left(5\right)$

The 3rd step, equation (5) is carried out method of moment find the solution, obtain matrix equation, adopt multilayer multistage sub-technology acceleration Matrix Solving fast;

The 4th step, when the distance between the field source meets some requirements, only shift the angular spectrum component on the field source line direction, transfer amount is dropped to minimum, distance condition is:

r _mn＞＞3γD (6)

γ is the control parameter, and γ is more little, and more many non-NULL groups can realize that the loss of accordingly result precision is also just more big with quick far-field approximation, and block number far away increases gradually and reduces along with γ;

The 5th step, find the solution matrix equation, obtain current coefficient, calculate the electromagnetic scattering parameter according to reciprocal theorem by current coefficient again.

The method of unknown quantity conversion is as follows in described second step:

Normal direction uncontinuity by polarization current on interface can get

$(J_{\mathrm{pol},1}-J_{\mathrm{pol},2})\·\hat{n}=-{\mathrm{j\ω\ρ}}_S---\left(7\right)$

Have two kinds of interfaces in the object module, the 1st kind: 1 be medium, and 2 be metal, i.e. medium lower surface, i.e. and satellite and coating interface are owing to the interior J of metal _{Pol, 2}=0, therefore on the medium lower surface, can get

$J_{\mathrm{pol},1}\·\hat{n}=-{\mathrm{j\ω\ρ}}_{S,\mathrm{pol}}---\left(8\right)$

The 2nd kind: 1 is air, and 2 is medium, and namely the medium upper surface namely applies and the air interface, same owing to J in the air _{Pol, 1}=0, therefore on the medium upper surface, can get

$J_{\mathrm{pol},2}\·\hat{n}=-{\mathrm{j\ω\ρ}}_{S,\mathrm{pol}}---\left(9\right)$

For polarization current and the lip-deep electric current of satellite are connected, utilize Maxwell equation to get

$J_{\mathrm{pol}}=\mathrm{j\ω}(\mathrm{\ϵ}-{\mathrm{\ϵ}}_0)E---\left(10\right)$

On the satellite surface

Therefore formula (10) can be write as

$J_{\mathrm{pol}}={\mathrm{j\ω\κ\ρ}}_S\hat{n}---\left(11\right)$

Wherein $\mathrm{\κ}=\frac{\mathrm{\ϵ}-{\mathrm{\ϵ}}_0}{\mathrm{\ϵ}};$

With the current continuity equation

${\mathrm{\ρ}}_S=-\frac{1}{\mathrm{j\ω}}{\▿}_S\·J_S---\left(12\right)$

Substitution formula (11),

$J_{\mathrm{pol}}=-\mathrm{\κ}(\▿\·J_S)\hat{n}---\left(13\right)$

Convolution (8), (11) and formula (12) can get at the medium lower surface

${\mathrm{\ρ}}_{S,\mathrm{pol}}=-\frac{\mathrm{\κ}}{\mathrm{j\ω}}\▿\·J_S---\left(14\right)$

In like manner can get at the medium upper surface

${\mathrm{\ρ}}_{S,\mathrm{pol}}=\frac{\mathrm{\κ}}{\mathrm{j\ω}}\▿\·J_S---\left(15\right)$

Ignore normal direction magnetic field, simultaneously tangential magnetic field is a constant, obtains following formula according to the magnetic field boundaries condition of perfect conductor:

$M_{\mathrm{pol}}\left(r^{\′}\right)\≈\mathrm{j\ω}(\mathrm{\μ}-{\mathrm{\μ}}_0)H_t\left(r^{\′}\right)=\mathrm{j\ω}\frac{({\mathrm{\μ}}_r-1)}{{\mathrm{\μ}}_r}{B}_t\left(r^{\′}\right)---\left(16\right)$

μ _rFor the relative permeability of thin-medium, with B _t(r') do distortion:

$B_t\left(r^{\′}\right)=-{\hat{n}}^{\′}\×{\hat{n}}^{\′}\×B_t\left(r^{\′}\right)=-\mathrm{\μ}{\hat{n}}^{\′}\×{\hat{n}}^{\′}\×H_t\left(r^{\′}\right)=-\mathrm{\μ}{\hat{n}}^{\′}\×J_S\left(r^{\′}\right)---\left(17\right)$

(17) formula is updated in (16) obtains M _PolSensitive surface electric current J with the metal surface _sThe expression formula of expression is as follows:

$M_{\mathrm{pol}}\left(r^{\′}\right)=-{\mathrm{j\ω\μ}}_0({\mathrm{\μ}}_r-1){\hat{n}}^{\′}\×J_S\left(r^{\′}\right)---\left(18\right)$

To obtain in (12), (13), (14), (15), (18) substitution (2):

$[{\mathrm{j\ω\μ}}_0\underset{S}{\∫}{J}_s\left(r^{\′}\right)G(r,r^{\′}){\mathrm{dS}}^{\′}-\frac{1}{{\mathrm{j\ω\ϵ}}_0}\underset{S}{\∫}{\▿}^{\′}\·J_s\left(r^{\′}\right)\▿G(r,r^{\′}){\mathrm{dS}}^{\′}$

$-{\mathrm{j\ω\μ}}_0\mathrm{\κ}\underset{V}{\∫}{\▿}^{\′}\·J_s\left(r^{\′}\right)\hat{n}G(r,r^{\′}){\mathrm{dS}}^{\′}+\frac{\mathrm{\κ}}{{\mathrm{j\ω\ϵ}}_0}\underset{S}{\∫}{\▿}^{\′}\·J_s\left(r^{\′}\right)\▿G(r,r^{\′}){\mathrm{dS}}^{\′}-\frac{\mathrm{\κ}}{{\mathrm{j\ω\ϵ}}_0}\underset{S^{\mathrm{\Δ}}}{\∫}{\▿}^{\′}\·J_s\left(r^{\′}\right)\▿G(r,r^{\′}){\mathrm{dS}}^{\mathrm{\Δ}\′}$

$+{\mathrm{j\ω\μ}}_0({\mathrm{\μ}}_r-1)\underset{V}{\∫}{\hat{n}}^{\′}\×J_S\left(r^{\′}\right)\×\▿G(r,r^{\′}){\mathrm{dV}}^{\′}{]}_{\tan }=E_{\tan }^{\mathrm{inc}}---\left(5\right)$

The method that adopts the quick multistage sub-technology of multilayer to accelerate Matrix Solving in described the 3rd step is as follows:

Step 3.1 is set up the tree structure of satellite;

Step 3.2 is supposed in given one deck, and a field point r is arranged among the arbitrary group of m _i, a source point r is arranged among the arbitrary group of n _j, r _m, r _nIt is respectively the central point of a some group and source point group; The space vector of two non-NULL group midfield points and source point is designated as: r _Ij=r _i-r _j=r _Im+ r _Mn+ r _Nj, when field point group and the source point group analyzed do not overlap when non-conterminous yet, satisfy | r _Im+ r _Nj|＜| r _Mn|, the scalar Green function of free space can be write as in the angular spectrum space:

$\frac{e^{-\mathrm{jk}|r_i-r_j|}}{|r_i-r_j|}=\frac{\mathrm{jk}}{4\mathrm{\π}}\underset{S_E}{\∫}{e}^{-\mathrm{jk}\·(r_{\mathrm{im}}+r_{\mathrm{nj}})}{\mathrm{\α}}_{\mathrm{mn}}(\hat{k}\·{\hat{r}}_{\mathrm{mn}})---\left(19\right)$

Whole integration is to be defined in unit ball S _EOn, α _MnBe the transfer factor of two groups, be defined as follows:

${\mathrm{\α}}_{\mathrm{mn}}=\underset{l=0}{\overset{L}{\mathrm{\Σ}}}{(-j)}^l(2l+1)h_l^{\left(2\right)}\left(\mathrm{kr}\right)P_l(\hat{k}\·\hat{r})---\left(20\right)$

Wherein, j _l(kd) be ball shellfish plug function,

Be the second class ball Hankel function, For Legendre function L is the length of blocking of infinite series, L=kD+ β (kD) ^1/3, D is the size of grouping, β is precision parameter, β 〉=2; With (θ, φ) representation unit ball S _ECoordinate, integration are counted and are K _L=2L ², wherein gather L the one dimension Gauss integration of putting altogether in the θ direction, the φ direction is gathered 2L trapezoidal rule integration of counting altogether;

Step 3.3, when the size of two far field effect group switching centre distances and group was determined, the transition matrix that a plurality of transfer factor are constituted added a window function:

$w_l=\left\{\begin{array}{cc}1& l\≤L/2\\ \sin \left(\frac{1}{L}\mathrm{\π}\right)& l>L/2\end{array}\right.---\left(21\right)$

Be multiplied by window function in transfer factor and make near the value of the transition component the two group switching centre line directions with respect to the value of other parts transfer factor precipitous variation be arranged, omit the effect to the less angular spectrum component of transfer process influence.

The method of tree structure of setting up satellite in the described step 3.1 is as follows: satellite is put in the cubic space that can surround it, remember that this volume coordinate is XYZ, this space is designated as the 0th layer of tree structure, then to being divided into two respectively on this space three-dimensional direction, form 8 sub spaces, be designated as ground floor, continue segmentation on the basis of 8 sub spaces, when the size of the thinnest sheaf space is long smaller or equal to the electric wave of setting, then stop segmentation; Whole satellite is transferred on the tree structure, is encoded in the space of each layer and sets up index, preserves the number in each layer non-NULL group space.

The wavelength of setting in the described step 3.1 is between 0.2 wavelength to 0.5 wavelength.

Judge the method for the less angular spectrum component of transfer process influence as follows in the described step 3.3: establish θ _rBe the angle of two effect group line directions and space z axle, φ _rBe the angle of two effect group line direction vectors projection and x axle on the XOY face, θ _eBe the angular spectrum component direction of preset threshold correspondence and the angle of two effect group line directions; In the process of calculating transfer factor, the size of transfer factor mould value increases along with θ and reduces, and extraction and calculating angular spectrum component direction and two effect group switching centre line angular separations are not more than θ _eThe angular spectrum component.

The present invention compared with prior art, its remarkable advantage: (1) equation is set up simple.At this metal medium mixed structure of satellite, only need metal is partly set up the field integral equation, the medium part does not need to set up equation and is described.(2) required computational resource is few.Because satellite is the design feature of length tiltedly, the high directivity of transfer factor component.And in the ray-tracing scheme owing to only need to shift angular spectrum component around the field source line direction, in the quick multistage son of multilayer, all the angular spectrum components on the unit ball are shifted, saved the computing time of transfer factor and required stored memory.Ought satisfy simultaneously under the condition of quick far-field approximation, only need calculate the angular spectrum component of transfer factor on field source group switching centre line direction, further save computational resource under the condition of computational accuracy satisfying.

Description of drawings

Fig. 1 is washing film dielectric layer object delineation of the present invention, and a is the geometric model synoptic diagram, and b is electromagnetic parameter model structure synoptic diagram.

Fig. 2 is the uncontinuity synoptic diagram of the normal component of border of the present invention polarization current.

Fig. 3 is multilayer Fast Multiple Method of the present invention far field effect synoptic diagram.

Fig. 4 is the structure of dummy satellite of the present invention and applies synoptic diagram

Fig. 5 is the two station of the present invention RCS curve synoptic diagram.

Embodiment

The numerical method of satellite electromagnetic scattering characteristic in the millimere-wave band of the present invention, step is as follows:

The first step, the foundation of electromagnetic scattering integral equation, namely based on the electric field boundary condition of perfect conductor, the resultant field tangential component in the metal surface is 0, and resultant field is incident electric field and scattering electric field sum.The incident electric field is known excitation, and uniform plane wave is used as the incident electric field usually, and scattering electric field can be represented with surperficial unknown current to be asked.

In second step, the conversion between the unknown quantity has 5 unknown quantitys in the field integral equation of namely setting up in the first step, be respectively metal surface current density, metal surface electric density, dielectric polarized current density, the upper and lower surface polarization charge densities of dielectric layer.Adopting method of moment to find the solution in the process of this equation, because a test item is only arranged: the metal surface current density obviously is to solve the correct result of this equation, therefore must all use the metal surface current density to represent all unknown sources;

The 3rd step, in conjunction with the characteristics of satellite structure, utilize the ultimate principle of ray tracing, reasonably select the big part of angular spectrum component in the transfer factor, the part that the angular spectrum component is little is given up to fall.

In the 4th step, in conjunction with the characteristics of satellite structure, calculate the far field effect fast in the method for the quick far-field approximation of more high-rise employing.

In the 5th step, matrix equation is found the solution and the electromagnetic scattering CALCULATION OF PARAMETERS.

Below in conjunction with accompanying drawing the present invention is described in further detail.

In conjunction with Fig. 1, the present invention is based on a kind of numerical method of analyzing the numerical method of satellite electromagnetic scattering characteristic in the millimere-wave band, step is as follows:

The first step makes uniform plane wave shine on the satellite structure that scribbles film dielectric layer, will produce surface induction surface current J on the surface of satellite _SWith surface charge ρ _S, produce body polarization electric current J in the thin-medium _Pol, body polarization magnetic current M _Pol, polarization charge ρ _{S, pol}, according to the electric field boundary condition of perfect conductor, namely the resultant field tangential component of metal surface is 0, obtains scribbling the field integral equation (EFIE) of the satellite structure target of film dielectric layer, and is as follows

[E ^Inc(r)+E ^Sca(r)] _Tan=0 (1) wherein, and subscript tan represents the tangential component of electric field, E ^IncExpression is radiated at the electromagnetic wave on the target, E ^ScaThe end of a performance that the expression target produces behind electromagnetic wave irradiation, and scattered field comprises metal and the two-part effect of thin-medium, the expression-form of scattered field is:

$E^{\mathrm{sca}}=-{\mathrm{j\ω\μ}}_0\underset{S}{\∫}{J}_S\left(r^{\′}\right)\stackrel{\‾}{G}(r,r^{\′}){\mathrm{dS}}^{\′}-{\mathrm{j\ω\μ}}_0\underset{V}{\∫}{J}_{\mathrm{pol}}\left(r^{\′}\right)\stackrel{\‾}{G}(r,r^{\′}){\mathrm{dV}}^{\′}+\underset{V}{\∫}{M}_V\left(r^{\′}\right)\×\▿G(r,r^{\′}){\mathrm{dV}}^{\′}$

$=-{\mathrm{j\ω\μ}}_0\underset{S}{\∫}{J}_S\left(r^{\′}\right)G(r,r^{\′}){\mathrm{dS}}^{\′}-\frac{1}{{\mathrm{\ϵ}}_0}\underset{S}{\∫}{\mathrm{\ρ}}_s\left(r^{\′}\right)G(r,r^{\′}){\mathrm{dS}}^{\′}$

(2)

$-{\mathrm{j\ω\μ}}_0\underset{V}{\∫}{J}_{\mathrm{pol}}\left(r^{\′}\right)G(r,r^{\′}){\mathrm{dS}}^{\′}-\frac{1}{{\mathrm{\ϵ}}_0}\underset{S}{\∫}{\mathrm{\ρ}}_{s,\mathrm{pol}}\left(r^{\′}\right)G(r,r^{\′}){\mathrm{dS}}^{\′}-\frac{1}{{\mathrm{\ϵ}}_0}\underset{S^{\mathrm{\Δ}}}{\∫}{\mathrm{\ρ}}_{s,\mathrm{pol}}\left(r^{\′}\right)G(r,r^{\′}){\mathrm{dS}}^{\mathrm{\Δ}\′}$

$+\underset{V}{\∫}{M}_V\left(r^{\′}\right)\×\▿G(r,r^{\′}){\mathrm{dV}}^{\′}$

Wherein V represents the volume element of thin-medium coating, and S represents unit, metal surface (being the lower surface unit of thin-medium coating), S ^ΔThe upper surface unit of expression thin-medium coating, ω is electromagnetic angular frequency, μ ₀And ε ₀Be respectively magnetic permeability and dielectric parameter in the vacuum, r and r' are respectively the position coordinates in field and source, G (r r') is the Green function of free space, and expression formula is:

$G(r,r^{\′})=\frac{e^{-\mathrm{jk}|r-r^{\′}|}}{4\mathrm{\π}|r-r^{\′}|}---\left(3\right)$

Wherein k is the wave number of free space.

Second step had 5 unknown quantitys in the field integral equation of setting up in the first step, be respectively the metal surface current density, J _S, metal surface electric density ρ _S, dielectric polarized current density J _Pol, the upper and lower surface polarization charge densities ρ of dielectric layer _{S, pol}, body polarization magnetic current M _PolIf, adopt method of moment to find the solution this equation, a test item is only arranged: the metal surface current density, J _S, obviously be to solve the correct result of this equation, therefore all unknown sources all must be used for the metal surface current density to represent.

The uncontinuity of the normal component of border polarization current as shown in Figure 2.Normal direction uncontinuity by polarization current on satellite and coating interface, coating and air interface can get

$(J_{\mathrm{pol},1}-J_{\mathrm{pol},2})\·\hat{n}=-{\mathrm{j\ω\ρ}}_S---\left(4\right)$

Have two kinds of interfaces among the present invention in the object module, the 1st kind: 1 be medium, and 2 be metal, i.e. medium lower surface is owing to the interior J of metal _{Pol, 2}=0, therefore on the medium lower surface, can get

$J_{\mathrm{pol},1}\·\hat{n}=-{\mathrm{j\ω\ρ}}_{S,\mathrm{pol}}---\left(5\right)$

The 2nd kind: 1 is air, and 2 is medium, and namely the medium upper surface is same owing to J in the air _{Pol, 1}=0, therefore on the medium upper surface, can get

$J_{\mathrm{pol},2}\·\hat{n}=-{\mathrm{j\ω\ρ}}_{S,\mathrm{pol}}---\left(6\right)$

For polarization current and the lip-deep electric current of satellite are connected, utilize Maxwell equation to get

J _Pol=j ω (ε-ε ₀) E (7) is on the satellite surface

Therefore formula (7) can be write as

$J_{\mathrm{pol}}=-{\mathrm{j\ω\κ\ρ}}_S\hat{n}---\left(8\right)$

Wherein $\mathrm{\κ}=\frac{\mathrm{\ϵ}-{\mathrm{\ϵ}}_0}{\mathrm{\ϵ}}.$

With the current continuity equation

${\mathrm{\ρ}}_S=-\frac{1}{\mathrm{j\ω}}{\▿}_S\·J_S---\left(9\right)$

Substitution formula (8),

$J_{\mathrm{pol}}=-\mathrm{\κ}(\▿\·J_S)\hat{n}---\left(10\right)$

By formula (9) and formula (10) as can be seen the electric density on the metal and the polarized current density in the dielectric form that all has been transformed into the metal covering current density represent;

Convolution (5), (8) and formula (9) can get at the medium lower surface

${\mathrm{\ρ}}_{S,\mathrm{pol}}=-\frac{\mathrm{\κ}}{\mathrm{j\ω}}\▿\·J_S---\left(11\right)$

Be that the satellite coating can get with the air interface at the medium upper surface in like manner

${\mathrm{\ρ}}_{S,\mathrm{pol}}=\frac{\mathrm{\κ}}{\mathrm{j\ω}}\▿\·J_S---\left(12\right)$

Consider that normal direction magnetic field is to ignore, tangential magnetic field can be thought a constant simultaneously, so can obtain following formula according to the magnetic field boundaries condition of perfect conductor:

$M_{\mathrm{pol}}\left(r^{\′}\right)\≈\mathrm{j\ω}(\mathrm{\μ}-{\mathrm{\μ}}_0)H_t\left(r^{\′}\right)=\mathrm{j\ω}\frac{({\mathrm{\μ}}_r-1)}{{\mathrm{\μ}}_r}{B}_t\left(r^{\′}\right)---\left(13\right)$

μ _rRelative permeability for thin-medium.With B _t(r') do simple distortion:

$B_t\left(r^{\′}\right)=-{\hat{n}}^{\′}\×{\hat{n}}^{\′}\×B_t\left(r^{\′}\right)=-\mathrm{\μ}{\hat{n}}^{\′}\×{\hat{n}}^{\′}\×H_t\left(r^{\′}\right)=-\mathrm{\μ}{\hat{n}}^{\′}\×J_S\left(r^{\′}\right)$

(14) formula is updated in (13) can obtains M _PolSensitive surface electric current J with the metal surface _sThe expression formula of expression is as follows:

$M_{\mathrm{pol}}\left(r^{\′}\right)=-{\mathrm{j\ω\μ}}_0({\mathrm{\μ}}_r-1){\hat{n}}^{\′}\×J_S\left(r^{\′}\right)---\left(15\right)$

To can obtain in (9), (10), (111), (12), (15) substitution (2):

$[{\mathrm{j\ω\μ}}_0\underset{S}{\∫}{J}_s\left(r^{\′}\right)G(r,r^{\′}){\mathrm{dS}}^{\′}-\frac{1}{{\mathrm{j\ω\ϵ}}_0}\underset{S}{\∫}{\▿}^{\′}\·J_s\left(r^{\′}\right)\▿G(r,r^{\′}){\mathrm{dS}}^{\′}$

$-{\mathrm{j\ω\μ}}_0\mathrm{\κ}\underset{V}{\∫}{\▿}^{\′}\·J_s\left(r^{\′}\right)\hat{n}G(r,r^{\′}){\mathrm{dS}}^{\′}+\frac{\mathrm{\κ}}{{\mathrm{j\ω\ϵ}}_0}\underset{S}{\∫}{\▿}^{\′}\·J_s\left(r^{\′}\right)\▿G(r,r^{\′}){\mathrm{dS}}^{\′}-\frac{\mathrm{\κ}}{{\mathrm{j\ω\ϵ}}_0}\underset{S^{\mathrm{\Δ}}}{\∫}{\▿}^{\′}\·J_s\left(r^{\′}\right)\▿G(r,r^{\′}){\mathrm{dS}}^{\mathrm{\Δ}\′}$

$+{\mathrm{j\ω\μ}}_0({\mathrm{\μ}}_r-1)\underset{V}{\∫}{\hat{n}}^{\′}\×J_S\left(r^{\′}\right)\×\▿G(r,r^{\′}){\mathrm{dV}}^{\′}{]}_{\tan }=E_{\tan }^{\mathrm{inc}}---\left(16\right)$

The 3rd step, equation (16) is carried out method of moment find the solution, obtain matrix equation.Green function to free space launches according to addition theorem, utilizes the ultimate principle of ray tracing, and the combination area divides the characteristic of equation reasonably to select the big part of angular spectrum component in the transfer factor, and the part that the angular spectrum component is little is given up to fall.

The millimere-wave band satellite belongs to the TV university size objectives, and multilayer multistage sub-technology fast is used to accelerate Matrix Solving.Multilayer fast multistage sub-technology at first needs to set up tree structure: satellite is put in the cubic space that can surround it, remember that this volume coordinate is XYZ, this space is designated as the 0th layer of tree structure, then to being divided into two respectively on this space three-dimensional direction, form 8 sub spaces, be designated as ground floor, think and analogize segmentation, then stop segmentation up to the size of thin sheaf space smaller or equal to certain electric wave long (between 0.2 wavelength to 0.5 wavelength).Whole satellite just is transferred on the tree structure, is encoded in the space of each layer and sets up index, preserves the number of each layer non-NULL group (space), and all calculating is all carried out in the non-NULL group.

Suppose in given one deck among the group m field point r is arranged _iWith organizing among the n source point r is arranged _j, r _m, r _nIt is respectively the central point of a some group and source point group.The space vector of two non-NULL group midfield points and source point is designated as: r _Ij=r _i-r _j=r _Im+ r _Mn+ r _Nj, when field point group and the source point group analyzed do not overlap when non-conterminous yet, satisfy | r _Im+ r _Nj|＜| r _Mn|, the scalar Green function of free space can be write as in the angular spectrum space:

$\frac{e^{-\mathrm{jk}|r_i-r_j|}}{|r_i-r_j|}=\frac{\mathrm{jk}}{4\mathrm{\π}}\underset{S_E}{\∫}{e}^{-\mathrm{jk}\·(r_{\mathrm{im}}+r_{\mathrm{nj}})}{\mathrm{\α}}_{\mathrm{mn}}(\hat{k}\·{\hat{r}}_{\mathrm{mn}})---\left(17\right)$

Whole integration is to be defined in unit ball S _EOn, α _MnBe the transfer factor of two groups, be defined as follows:

${\mathrm{\α}}_{\mathrm{mn}}=\underset{l=0}{\overset{L}{\mathrm{\Σ}}}{(-j)}^l(2l+1)h_l^{\left(2\right)}{\left(\mathrm{kr}\right)P}_l(\hat{k}\·\hat{r})---\left(18\right)$

Wherein, j _l(kd) be ball shellfish plug function,

Be the second class ball Hankel function,

Be Legendre function.L is the length of blocking of infinite series, L=kD+ β (kD) ^1/3, D is the size of grouping, β is empirical parameter, gets 2 among the present invention.With (θ, φ) representation unit spherical coordinates, integration count and be K _L=2L ², wherein gather L the one dimension Gauss integration of putting altogether in the θ direction, the φ direction is gathered 2L trapezoidal rule integration of counting altogether.

When the size of two far field effect group switching centre distances and group was determined, the property class of transfer factor was similar to the antenna pattern characteristic of linear antenna arrays, and transition matrix is added a window function:

$w_l=\left\{\begin{array}{cc}1& l\≤L/2\\ \sin \left(\frac{1}{L}\mathrm{\π}\right)& l>L/2\end{array}\right.---\left(19\right)$

Wherein, L is the number of L discrete point of θ direction under the spherical coordinates, and l is the numbering of each transfer factor component.Concrete way is: be multiplied by the effect of window function in transfer factor, make near the value of the transition component the two group switching centre line directions with respect to the value of other parts transfer factor precipitous variation be arranged, in this case, can omit the effect to the less angular spectrum component of transfer process influence, thereby make the ray multipole obtain higher efficient.Simultaneously, can verify that this processing is very little to the computational accuracy influence.Then need to set the scope of the angular spectrum component that needs transfer.As shown in Figure 3, θ _rBe the angle of two effect group line directions and z axle, φ _rBe the angle of two effect group line direction vectors projection and x axle on the XOY face, θ _eBe the angular spectrum component direction of preset threshold correspondence and the angle of two effect group line directions.In calculating the process of transfer factor, the size of transfer factor angular spectrum component is along with reducing with the increase of the angle of far field effect group switching centre line direction.Can judge the size of conical area the transfer factor setting threshold like this.In concrete operations, only need extraction and calculating angular spectrum component direction and two effect group switching centre line angular separations to be not more than θ _eThe angular spectrum component.When proof of algorithm, can obtain the effect that this part angular spectrum component can contact two far field effect groups well.

In the 4th step, the combination area divides the characteristic of equation, calculates the far field effect fast in the method for the quick far-field approximation of more high-rise employing.In the far-field approximation method, when the distance between the field source met some requirements, we can only calculate the angular spectrum component on the field source line direction, like this transfer amount were dropped to minimum fast.Adopt the criterion condition of quick far-field approximation to be between employing group and the group among the present invention:

r _mn＞3γD (20)

Wherein D is the electric size of this layer group, and γ is taken as 1.5 for the empirical value that the combination area divides the equation characteristics to arrange among the present invention.γ is more little, and more many non-NULL groups can realize that the loss of accordingly result precision also will be more big with quick far-field approximation, and block number far away increases gradually and reduces along with γ.The quick multistage submethod of quick far-field approximation, ray tracing and multilayer can be combined thus, can adopt three kinds of fast methods to handle the effect between each layer group and the group based on the distance between the group of far field: 1 Fast Multiple Method be accelerated the far field effect; 2 ray multipoles accelerate the far field effect; 3 quick far-field approximations are accelerated the far field effect.

The 5th step, find the solution matrix equation, obtain current coefficient, calculate the electromagnetic scattering parameter according to reciprocal theorem by current coefficient again.

For correctness and the validity of verifying the inventive method, provided the example of calculating Aerospace Satellite RCS in the millimere-wave band below, and result of calculation compares with traditional algorithm all, finely identical.

The unitized construction of Aerospace Satellite is made up of two parts, and it is 2m cube celestial body that there is a length of side centre, and at xoy two length being arranged is the solar panel of 10.87m.Four sides, the front and back up and down coated media of celestial body is (3.5 ,-0.008).The coated media of solar panel is (11.9,0.0).The incident wave frequency is 35GHz, incident wave line of propagation θ=0 ° → 180 °, φ=0 °.Subdivision is of a size of 0.14 wavelength, and the unknown quantity number is 125160048, and the thinnest layer packet size: 0.2 wavelength, the number of plies are 12.As can be seen from Figure 4, add quick far-field approximation and ray multipole after, compare with classic method, RCS result is identical substantially, satisfies accuracy requirement.And in table 1 as can be seen, added quick far-field approximation and ray multipole after, be to use the number of plies of quick far-field approximation at LC=8(LC) time, saved 49.6% than original program the computing time of program; Internal memory has saved 52.8%.Therefore far-field approximation and ray multipole can be saved computing time and internal memory fast.

Table 1 calculates internal memory and time statistical form.

Claims (6)

1. the numerical method of satellite electromagnetic scattering characteristic in the millimere-wave band is characterized in that step is as follows:

The first step makes uniform plane wave shine on the satellite structure that scribbles film dielectric layer, will produce surface induction surface current J on the surface of satellite _s, surface charge ρ _S, produce body polarization electric current J in the thin-medium _Pol, body polarization magnetic current M _PolWith polarization charge p _{S, pol}, according to the electric field boundary condition of perfect conductor, namely the resultant field tangential component of metal surface is 0, obtains scribbling the field integral equation EFIE of the satellite structure target of film dielectric layer, and is as follows

[E ^inc(r)+E ^sca(r)] _tan＝0

(1)

Wherein, subscript tan represents the tangential component of electric field, E ^IncExpression is radiated at the electromagnetic wave on the target, E ^ScaThe end of a performance that the expression target produces behind electromagnetic wave irradiation, the expression-form of scattered field is:

$E^{\mathrm{sca}}=-{\mathrm{j\ω\μ}}_0\underset{S}{\∫}{J}_S\left(r^{\′}\right)\stackrel{\‾}{G}(r,r^{\′}){\mathrm{dS}}^{\′}-{\mathrm{j\ω\μ}}_0\underset{V}{\∫}{J}_{\mathrm{pol}}\left(r^{\′}\right)\stackrel{\‾}{G}(r,r^{\′}){\mathrm{dV}}^{\′}+\underset{V}{\∫}{M}_V\left(r^{\′}\right)\×\▿G(r,r^{\′}){\mathrm{dV}}^{\′}$

$=-{\mathrm{j\ω\μ}}_0\underset{S}{\∫}{J}_S\left(r^{\′}\right)G(r,r^{\′}){\mathrm{dS}}^{\′}-\frac{1}{{\mathrm{\ϵ}}_0}\underset{S}{\∫}{\mathrm{\ρ}}_s\left(r^{\′}\right)G(r,r^{\′}){\mathrm{dS}}^{\′}$

$-{\mathrm{j\ω\μ}}_0\underset{V}{\∫}{J}_{\mathrm{pol}}\left(r^{\′}\right)G(r,r^{\′}){\mathrm{dS}}^{\′}-\frac{1}{{\mathrm{\ϵ}}_0}\underset{S}{\∫}{\mathrm{\ρ}}_{s,\mathrm{pol}}\left(r^{\′}\right)G(r,r^{\′}){\mathrm{dS}}^{\′}-\frac{1}{{\mathrm{\ϵ}}_0}\underset{S^{\mathrm{\Δ}}}{\∫}{\mathrm{\ρ}}_{s,\mathrm{pol}}\left(r^{\′}\right)G(r,r^{\′}){\mathrm{dS}}^{\mathrm{\Δ}\′}$

$+\underset{V}{\∫}{M}_V\left(r^{\′}\right)\×\▿G(r,r^{\′}){\mathrm{dV}}^{\′}---\left(2\right)$

Wherein, V represents the volume element of thin-medium coating, and S represents that the unit, metal surface is the lower surface unit of thin-medium coating, S ^ΔThe upper surface unit of expression thin-medium coating, ω is electromagnetic angular frequency, μ ₀And ε ₀Be respectively magnetic permeability and dielectric parameter in the vacuum, r and r' are respectively the position coordinates in field and source, G (r r') is the Green function of free space, and expression formula is:

$G(r,r^{\′})=\frac{e^{-\mathrm{jk}|r-r^{\′}|}}{4\mathrm{\π}|r-r^{\′}|}---\left(3\right)$

Wherein k is the wave number of free space;

In second step, convert 5 unknown quantitys in the field integral equation of setting up in the first step to the metal surface current density, J _S, 5 unknown quantitys are respectively the metal surface current density, J _S, metal surface electric density ρ _S, dielectric polarized current density J _Pol, dielectric layer upper and lower surface polarization charge densities ρ _{S, pol}With body polarization magnetic current M _Pol, its concrete representation is as follows,

Formula (1) finally becomes:

$[{\mathrm{j\ω\μ}}_0\underset{S}{\∫}{J}_s\left(r^{\′}\right)G(r,r^{\′}){\mathrm{dS}}^{\′}-\frac{1}{{\mathrm{j\ω\ϵ}}_0}\underset{S}{\∫}{\▿}^{\′}\·J_s\left(r^{\′}\right)\▿G(r,r^{\′}){\mathrm{dS}}^{\′}$

$-{\mathrm{j\ω\μ}}_0\mathrm{\κ}\underset{V}{\∫}{\▿}^{\′}\·J_s\left(r^{\′}\right)\hat{n}G(r,r^{\′}){\mathrm{dS}}^{\′}+\frac{\mathrm{\κ}}{{\mathrm{j\ω\ϵ}}_0}\underset{S}{\∫}{\▿}^{\′}\·J_s\left(r^{\′}\right)\▿G(r,r^{\′}){\mathrm{dS}}^{\′}-\frac{\mathrm{\κ}}{{\mathrm{j\ω\ϵ}}_0}\underset{S^{\mathrm{\Δ}}}{\∫}{\▿}^{\′}\·J_s\left(r^{\′}\right)\▿G(r,r^{\′}){\mathrm{dS}}^{\mathrm{\Δ}\′}$

$+{\mathrm{j\ω\μ}}_0({\mathrm{\μ}}_r-1)\underset{V}{\∫}{\hat{n}}^{\′}\×J_S\left(r^{\′}\right)\×\▿G(r,r^{\′}){\mathrm{dV}}^{\′}{]}_{\tan }=E_{\tan }^{\mathrm{inc}}---\left(5\right)$

The 3rd step, equation (5) is carried out method of moment find the solution, obtain matrix equation, adopt multilayer multistage sub-technology acceleration Matrix Solving fast;

The 4th step, when the distance between the field source meets some requirements, only shift the angular spectrum component on the field source line direction, transfer amount is dropped to minimum, distance condition is:

r _mn＞＞3γD

(6)

γ is the control parameter, and γ is more little, and more many non-NULL groups can realize that the loss of accordingly result precision is also just more big with quick far-field approximation, and block number far away increases gradually and reduces along with γ;

The 5th step, find the solution matrix equation, obtain current coefficient, calculate the electromagnetic scattering parameter according to reciprocal theorem by current coefficient again.

2. the numerical method of satellite electromagnetic scattering characteristic in the millimere-wave band according to claim 1 is characterized in that, the method for unknown quantity conversion is as follows in described second step:

Normal direction uncontinuity by polarization current on interface can get

$(J_{\mathrm{pol},1}-J_{\mathrm{pol},2})\·\hat{n}=-{\mathrm{j\ω\ρ}}_S---\left(7\right)$

Have two kinds of interfaces in the object module, the 1st kind: 1 be medium, and 2 be metal, i.e. medium lower surface, i.e. and satellite and coating interface are owing to the interior J of metal _{Pol, 2}=0, therefore on the medium lower surface, can get

$J_{\mathrm{pol},1}\·\hat{n}=-{\mathrm{j\ω\ρ}}_{S,\mathrm{pol}}---\left(8\right)$

The 2nd kind: 1 is air, and 2 is medium, and namely the medium upper surface namely applies and the air interface, same owing to J in the air _{Pol, 1}=0, therefore on the medium upper surface, can get

$J_{\mathrm{pol},2}\·\hat{n}=-{\mathrm{j\ω\ρ}}_{S,\mathrm{pol}}---\left(9\right)$

For polarization current and the lip-deep electric current of satellite are connected, utilize Maxwell equation to get

J _pol=jω(ε-ε ₀)E (10)

On the satellite surface

Therefore formula (10) can be write as

$J_{\mathrm{pol}}={\mathrm{j\ω\κ\ρ}}_S\hat{n}---\left(11\right)$

Wherein $\mathrm{\κ}=\frac{\mathrm{\ϵ}-{\mathrm{\ϵ}}_0}{\mathrm{\ϵ}};$

With the current continuity equation

${\mathrm{\ρ}}_S=-\frac{1}{\mathrm{j\ω}}{\▿}_S\·J_S---\left(12\right)$

Substitution formula (11),

$J_{\mathrm{pol}}=-\mathrm{\κ}(\▿\·J_S)\hat{n}---\left(13\right)$

Convolution (8), (11) and formula (12) can get at the medium lower surface

${\mathrm{\ρ}}_{S,\mathrm{pol}}=-\frac{\mathrm{\κ}}{\mathrm{j\ω}}\▿\·J_S---\left(14\right)$

In like manner can get at the medium upper surface

${\mathrm{\ρ}}_{S,\mathrm{pol}}=\frac{\mathrm{\κ}}{\mathrm{j\ω}}\▿\·J_S---\left(15\right)$

Ignore normal direction magnetic field, simultaneously tangential magnetic field is a constant, obtains following formula according to the magnetic field boundaries condition of perfect conductor:

$M_{\mathrm{pol}}\left(r^{\′}\right)\≈\mathrm{j\ω}(\mathrm{\μ}-{\mathrm{\μ}}_0)H_t\left(r^{\′}\right)=\mathrm{j\ω}\frac{({\mathrm{\μ}}_r-1)}{{\mathrm{\μ}}_r}{B}_t\left(r^{\′}\right)---\left(16\right)$

μ _rFor the relative permeability of thin-medium, with B _t(r') do distortion:

$B_t\left(r^{\′}\right)=-{\hat{n}}^{\′}\×{\hat{n}}^{\′}\×B_t\left(r^{\′}\right)=-\mathrm{\μ}{\hat{n}}^{\′}\×{\hat{n}}^{\′}\×H_t\left(r^{\′}\right)---\left(17\right)$

(17) formula is updated in (16) obtains M _PolSensitive surface electric current J with the metal surface _sThe expression formula of expression

As follows:

$M_{\mathrm{pol}}\left(r^{\′}\right)=-{\mathrm{j\ω\μ}}_0({\mathrm{\μ}}_r-1){\hat{n}}^{\′}\×J_S\left(r^{\′}\right)---\left(18\right)$

To obtain in (12), (13), (14), (15), (18) substitution (2):

$[{\mathrm{j\ω\μ}}_0\underset{S}{\∫}{J}_s\left(r^{\′}\right)G(r,r^{\′}){\mathrm{dS}}^{\′}-\frac{1}{{\mathrm{j\ω\ϵ}}_0}\underset{S}{\∫}{\▿}^{\′}\·J_s\left(r^{\′}\right)\▿G(r,r^{\′}){\mathrm{dS}}^{\′}$

$-{\mathrm{j\ω\μ}}_0\mathrm{\κ}\underset{V}{\∫}{\▿}^{\′}\·J_s\left(r^{\′}\right)\hat{n}G(r,r^{\′}){\mathrm{dS}}^{\′}+\frac{\mathrm{\κ}}{{\mathrm{j\ω\ϵ}}_0}\underset{S}{\∫}{\▿}^{\′}\·J_s\left(r^{\′}\right)\▿G(r,r^{\′}){\mathrm{dS}}^{\′}-\frac{\mathrm{\κ}}{{\mathrm{j\ω\ϵ}}_0}\underset{S^{\mathrm{\Δ}}}{\∫}{\▿}^{\′}\·J_s\left(r^{\′}\right)\▿G(r,r^{\′}){\mathrm{dS}}^{\mathrm{\Δ}\′}$

$+{\mathrm{j\ω\μ}}_0({\mathrm{\μ}}_r-1)\underset{V}{\∫}{\hat{n}}^{\′}\×J_S\left(r^{\′}\right)\×\▿G(r,r^{\′}){\mathrm{dV}}^{\′}{]}_{\tan }=E_{\tan }^{\mathrm{inc}}---\left(5\right)$

3. the numerical method of satellite electromagnetic scattering characteristic in the millimere-wave band according to claim 1 is characterized in that, the method that adopts the quick multistage sub-technology of multilayer to accelerate Matrix Solving in described the 3rd step is as follows:

Step 3.1 is set up the tree structure of satellite;

Step 3.2 is supposed in given one deck, and a field point r is arranged among the arbitrary group of m _i, a source point r is arranged among the arbitrary group of n _j, r _m, r _nIt is respectively the central point of a some group and source point group; The space vector of two non-NULL group midfield points and source point is designated as: r _Ij=r _i-r _j=r _Im+ r _Mn+ r _Nj, when field point group and the source point group analyzed do not overlap when non-conterminous yet, satisfy | r _Im+ r _Nj|＜| r _Mn|, the scalar Green function of free space can be write as in the angular spectrum space:

$\frac{e^{-\mathrm{jk}|r_i-r_j|}}{|r_i-r_j|}=\frac{\mathrm{jk}}{4\mathrm{\π}}\underset{S_E}{\∫}{e}^{-\mathrm{jk}\·(r_{\mathrm{im}}+r_{\mathrm{nj}})}{\mathrm{\α}}_{\mathrm{mn}}(\hat{k}\·{\hat{r}}_{\mathrm{mn}})---\left(19\right)$

Whole integration is to be defined in unit ball S _EOn, α _MnBe the transfer factor of two groups, be defined as follows:

${\mathrm{\α}}_{\mathrm{mn}}=\underset{l=0}{\overset{L}{\mathrm{\Σ}}}{(-j)}^l(2l+1)h_l^{\left(2\right)}\left(\mathrm{kr}\right)P_l(\hat{k}\·\hat{r})---\left(20\right)$

Wherein, j _l(kd) be ball shellfish plug function,

Be the second class ball Hankel function,

Be Legendre function; L is the length of blocking of infinite series, L=kD+ β (kD) ^1/3, D is the size of grouping, β is precision parameter, β 〉=2; With (θ, φ) representation unit ball S _ECoordinate, integration are counted and are K _L=2L ², wherein gather L the one dimension Gauss integration of putting altogether in the θ direction, the φ direction is gathered 2L trapezoidal rule integration of counting altogether;

Step 3.3, when the size of two far field effect group switching centre distances and group was determined, the transition matrix that a plurality of transfer factor are constituted added a window function:

$w_l=\left\{\begin{array}{cc}1& l\≤L/2\\ \sin \left(\frac{1}{L}\mathrm{\π}\right)& l>L/2\end{array}\right.---\left(21\right)$

Be multiplied by window function in transfer factor and make near the value of the transition component the two group switching centre line directions with respect to the value of other parts transfer factor precipitous variation be arranged, omit the effect to the less angular spectrum component of transfer process influence.

4. according to the numerical method of satellite electromagnetic scattering characteristic in claim 1 or the 3 described millimere-wave bands, it is characterized in that, the method of tree structure of setting up satellite in the described step 3.1 is as follows: satellite is put in the cubic space that can surround it, remember that this volume coordinate is XYZ, this space is designated as the 0th layer of tree structure, then to being divided into two respectively on this space three-dimensional direction, form 8 sub spaces, be designated as ground floor, continue segmentation on the basis of 8 sub spaces, when the size of the thinnest sheaf space is long smaller or equal to the electric wave of setting, then stop segmentation; Whole satellite is transferred on the tree structure, is encoded in the space of each layer and sets up index, preserves the number in each layer non-NULL group space.

5. according to the numerical method of satellite electromagnetic scattering characteristic in claim 1 or the 3 described millimere-wave bands, it is characterized in that: the wavelength of setting in the described step 3.1 is between 0.2 wavelength to 0.5 wavelength.

6. according to the numerical method of satellite electromagnetic scattering characteristic in claim 1 or the 3 described millimere-wave bands, it is characterized in that, judge the method for the less angular spectrum component of transfer process influence as follows in the described step 3.3: establish θ _rBe the angle of two effect group line directions and space z axle, φ _rBe the angle of two effect group line direction vectors projection and x axle on the XOY face, θ _eBe the angular spectrum component direction of preset threshold correspondence and the angle of two effect group line directions; In the process of calculating transfer factor, the size of transfer factor mould value increases along with θ and reduces, and extraction and calculating angular spectrum component direction and two effect group switching centre line angular separations are not more than θ _eThe angular spectrum component.

Images