KR101315569B1 - Method for numerical analysis of electromagnetic waves and apparatus for thereof - Google Patents

Abstract

According to the present invention, a surface of a scatterer is divided into a plurality of triangles by using a Rao-Wilton-Glisson (RWG) basis function, and each vertex coordinate of the triangle is obtained. Setting the configured first triangular pair as a power point region and a second triangular pair as a viewpoint point region, and using the vertex coordinates of the power point region and the viewpoint point region, the relative between the power point region and the viewpoint point region Calculating a distance index, comparing the relative distance index with a preset threshold, and using the Galerkin method if the relative distance index is less than the threshold, and interpolating if the relative distance index is less than the threshold. electromagnetic wave numerical analysis comprising calculating each element of the impedance matrix for the scatterer using a point interpolation method ≪ / RTI >
According to the electromagnetic wave numerical analysis method and apparatus according to the present invention, in order to effectively analyze the electromagnetic wave scattering characteristics of any metal scattering body, the moment method applying the electric field integral equation (EFIE) and the Rao-Wilton-Glisson (RWG) basis function ( When calculating the impedance matrix of the Method of Moment, the Galerkin method and the Center-point Interpolation Method can be used to increase the calculation efficiency and maintain the accuracy.

Description

Method for numerical analysis of electromagnetic waves and apparatus for

The present invention relates to an electromagnetic wave numerical analysis method and apparatus, and more particularly to an electromagnetic wave numerical analysis method and apparatus that can effectively analyze the electromagnetic wave scattering characteristics of any metal scattering body.

Various numerical techniques are used to analyze the scattering characteristics of randomly shaped scatterers. In particular, the moment method using integral equations such as electric field integral equation (EFIE), magnetic field integral equation (MFIE), and combined field integral equation (CFIE) for far-zone electromagnetic wave scattering and radiation analysis of moments) is common. The integral equation can be selected according to the material and analysis environment of the scatterer for numerical analysis.

In order to analyze the electromagnetic response characteristics of a scatterer having an arbitrary shape in a three-dimensional space by the moment method, various basis functions for expressing the scatterer surface current can be applied. For example, the pulse basis function has a wide range of use and a simple expression, but there is a limit to an arbitrary structural analysis using a rectangular coordinate system. On the other hand, RWG (Rao-Wilton-Glisson) basis function, which has been widely used and sufficiently verified in recent years, can easily calculate arbitrary surface scatterers and thus calculate optimal surface current. The prior art using the moment method using the pulse function as a basis function is described in the Korean Institute of Communication Sciences, 12 (6), 571-579, 1987, 'Electromagnetic Scattering Analysis of Full Conducting Polygonists Using Extrapolation Approximation'.

The surface current, which is expressed as a finite element for applying the moment method, can be represented by matrix expressions by point matching method and Galerkin method by base function and test function. Here, the Galerkin technique using the same test function as the RWG basis function can be applied for any scatterer analysis.

The Galerkin technique helps to improve the accuracy of calculations between arbitrary source point and observation point areas, but it takes time to calculate complex numerical integrations for all the areas that make up the scatterer surface. This has the growing disadvantage.

The present invention, when calculating the impedance matrix of the method of the moment (Method of Moment) applying the electric field integration Equation (EFIE) and the Rao-Wilton-Glisson (RWW) basis function to effectively analyze the electromagnetic scattering characteristics of any metal scatterer The purpose of this study is to provide a method and apparatus for electromagnetic wave numerical analysis that can improve the calculation efficiency and maintain the accuracy by using the Galerkin method and the center-point interpolation method.

According to the present invention, a surface of a scatterer is divided into a plurality of triangles by using a Rao-Wilton-Glisson (RWG) basis function, and each vertex coordinate of the triangle is obtained. Setting the configured first triangular pair as a power point region and a second triangular pair as a viewpoint point region, and using the vertex coordinates of the power point region and the viewpoint point region, the relative between the power point region and the viewpoint point region Calculating a distance index, comparing the relative distance index with a preset threshold, and using the Galerkin method if the relative distance index is less than the threshold, and interpolating if the relative distance index is less than the threshold. electromagnetic wave numerical analysis comprising calculating each element of the impedance matrix for the scatterer using a point interpolation method Provide a method.

The calculating of the relative distance index may include calculating a relative distance index between the center of gravity of one triangle in the power point region and the center of gravity of one triangle in the observation point region, and the relative distance index R _{adj _ cond} can be calculated by the following equation.

는 상기 전원점 영역 내의 한 개의 삼각형의 무게중심과 상기 관측점 영역 내의 한 개의 삼각형의 무게중심 사이의 거리를 나타낸다.Where L _s is the largest distance from the center of gravity of one triangle in the power point area to each vertex of two triangles in the power point area, L _o is the center of gravity of one triangle in the viewpoint area The largest of the distances to each vertex of the two triangles within the viewpoint area,

Denotes the distance between the center of gravity of one triangle in the power point area and the center of gravity of one triangle in the viewpoint area.

The power point region and the observation point region are adjacent to or spaced apart from each other, and the threshold is a value between 1 and 1.5, and when the relative distance index R _{adj _ cond} is smaller than the threshold, the element is in the impedance matrix. If the relative distance index R _{adj _ cond} is greater than the threshold, the element may be arranged near the upper or lower triangle in the impedance matrix.

In addition, the elements located on the main diagonal in the impedance matrix are impedance values corresponding to the relative distance index of 0, and as the distance between the power point region and the observation point region increases, the calculated element is the impedance matrix. And may be located at a distance from the element on the main diagonal.

The electromagnetic wave numerical analysis method may further include calculating an inverse matrix of the impedance matrix, and calculating a surface current matrix using the inverse matrix and the incident wave matrix of the impedance matrix.

According to the present invention, a splitter for dividing a surface of a scatterer into a plurality of triangles by using a Rao-Wilton-Glisson (RWG) basis function, and obtaining a vertex coordinate of each of the triangles, and an adjacent one of the plurality of triangles. An area setting unit configured to set a first triangle pair composed of two triangles as a power point region, and a second triangle pair as a viewpoint point region, and using the vertex coordinates of the power point region and the viewpoint point region, An exponent calculator for calculating a relative distance index between the observation point regions, a comparator for comparing the relative distance index with a preset threshold, and using the Galerkin method if the relative distance index is less than the threshold, and using the threshold value. If this is the case, Impedon calculates each element of the impedance matrix for the scatterer using the Center-point Interpolation Method. It provides an electromagnetic wave numerical analysis device including a switch operator.

Here, the exponential operator calculates a relative distance index between the center of gravity of one triangle in the power point area and the center of gravity of one triangle in the viewpoint area, and the relative distance index R _{adj _ cond} is Can be calculated through the equation.

Here, L _s is the largest distance of the distance from the center of gravity of one triangle to each vertex in the power point area, L _o is the largest distance of the distance from the center of gravity of one triangle in the viewpoint area to each vertex , R _s ^{c ±} -r _o ^{c ±} │ represents the distance between the center of gravity of one triangle in the power point region and the center of gravity of one triangle in the observation point region.

The electromagnetic wave numerical analysis device may further include an inverse matrix calculator configured to calculate an inverse matrix with respect to the impedance matrix, and a surface current calculator configured to calculate a surface current matrix using the inverse matrix and the incident wave matrix of the impedance matrix. .

According to the electromagnetic wave numerical analysis method and apparatus according to the present invention, in order to effectively analyze the electromagnetic wave scattering characteristics of any metal scattering body, the moment method applying the electric field integral equation (EFIE) and the Rao-Wilton-Glisson (RWG) basis function ( When calculating the impedance matrix of the Method of Moment, the Galerkin method and the Center-point Interpolation Method can be used to increase the calculation efficiency and maintain the accuracy.

1 shows an electromagnetic wave numerical analysis apparatus according to an embodiment of the present invention.
2 is a flowchart illustrating an electromagnetic wave numerical analysis method using FIG. 1.
FIG. 3 is an exemplary diagram in which a spherical scatterer surface is divided into several triangles using the RWG basis function in step S210 of FIG. 2.
FIG. 4 shows an example of the geometry of a triangular pair in step S220 of FIG. 2.
FIG. 5 shows an example of another geometry in step S220 of FIG. 2.
6 is a schematic diagram of an impedance matrix using an improved moment method according to an embodiment of the present invention.
7 shows the numerical integration results for specific triangular regions.
8 shows the calculation speed and the error according to the threshold value change of the relative distance index when the embodiment of the present invention is applied.
9 is a result of calculating the Bistatic RCS of the metal sphere when using the conventional moment method and the improved moment method according to the present invention.
Figure 10 shows the results of the analysis of the monostatic RCS characteristics of the three- / forward wave reflector according to an embodiment of the present invention.

DETAILED DESCRIPTION Hereinafter, exemplary embodiments of the present invention will be described in detail with reference to the accompanying drawings so that those skilled in the art may easily implement the present invention.

BACKGROUND OF THE INVENTION 1. Field of the Invention The present invention relates to a method and apparatus for numerical analysis of electromagnetic waves, and provides an improved moment analysis technique as compared to existing method of moments for electromagnetic numerical analysis in three-dimensional space.

1 shows an electromagnetic wave numerical analysis apparatus according to an embodiment of the present invention. The electromagnetic wave numerical analysis apparatus 100 includes a divider 110, an area setting unit 120, an exponential calculator 130, a comparator 140, an impedance calculator 150, an inverse matrix calculator 160, and a surface current calculator. And 170.

The divider 110 divides the surface of the metal scatterer into a plurality of triangles by using a Rao-Wilton-Glisson (RWG) basis function, and obtains coordinates of each vertex of the divided triangles.

The area setting unit 120 sets a first triangle pair composed of two adjacent triangles among the plurality of triangles as a power point region and a second triangle pair as an observation point region. Since the number of split triangles on the surface of the scatterer is very large, there are several power point regions and observation point regions. The power point area and the viewpoint area are adjacent to or spaced apart from each other.

The exponent calculator 130 calculates a relative distance index between the power point area and the observation point area by using the vertex coordinates of the power point area and the observation point area. The larger the physical separation between the viewpoint area and the vertex area, the larger the relative distance index will be.

The comparison unit 140 compares the relative distance index with a preset threshold. The threshold here has a value between 1 and 1.5, but the present invention is not necessarily limited thereto.

If the relative distance index is less than the threshold, the impedance calculator 150 calculates each element of the impedance matrix for the scatterer using the Galerkin method. Also, if the relative distance index is greater than or equal to the threshold, the impedance calculator 150 calculates each element of the impedance matrix for the scatterer using a center-point interpolation method. That is, the impedance calculation unit 150 selectively uses one of the Galerkin method and the center point interpolation method according to the value of the relative distance index during the impedance calculation.

Basically, the moment method calculates each element of the incident wave matrix [B] and each element of the impedance matrix [Z] according to the scatterer structure by using the matrix equation to obtain the surface current matrix [I], and inverse transformation of the impedance matrix. A series of processes for calculating all the elements of surface current from During this process, most of the time is spent, especially in calculating the impedance matrix.

In this regard, the existing double integration using the Galerkin technique shows a high reliability of the calculation results but takes a considerable amount of calculation time. In addition, the center interpolation method is a method of simplifying the integral calculation, which gives a much better gain in terms of calculation time, but has a disadvantage of inferior calculation accuracy.

In the present invention, when the relative distance index is greater than or equal to the threshold, the center point interpolation method is used instead of the Galerkin method in order to shorten the calculation time of the numerical integration for the arbitrary power point region and the observation point region. The method of the present invention is intended to compensate for the shortcomings of the two techniques, and it is possible to drastically reduce the time required for impedance calculation and to maintain the accuracy of the calculation. The accuracy problem according to the center point interpolation method will be solved by defining a relative distance index between the power point region and the observation point region which will be described later.

As described above, the present invention basically uses the electric field integral equation (EFIE) and the RWG basis function to analyze the electromagnetic wave scattering characteristics, and reduces the analysis time by using the existing Galerkin technique and the center point interpolation method to increase the computational efficiency. While maintaining the accuracy of the calculations.

The inverse matrix calculator 160 calculates an inverse matrix [Z] ⁻¹ for the impedance matrix [Z]. The surface current calculator 170 calculates the surface current matrix [I] using the inverse matrix [Z] ⁻¹ and the incident wave matrix [B] of the impedance matrix. This is due to the formula [B] = [Z] [I].

2 is a flowchart illustrating an electromagnetic wave numerical analysis method using FIG. 1. Hereinafter, the electromagnetic wave numerical analysis method will be described in detail with reference to FIG. 2.

First, the surface of the scatterer is divided into a plurality of triangles by using a Rao-Wilton-Glisson (RWG) basis function, and each vertex coordinate of the triangle is obtained (S210).

In general, the RWG basis function is a method of expressing the surface current flowing through the scattering body divided into triangles by the size and direction of the vertical component passing through the edge of each triangle. FIG. 3 is an exemplary diagram in which a spherical scatterer (metal sphere) surface is divided into several triangles using the RWG basis function in step S210 of FIG. 2. This RWG basis function facilitates the expression of arbitrary scatterers and enables the optimal surface current calculation.

The known simple principle for the moment method using this RWG basis function is as follows.

을 적용해 나타낸 것이다. Equations 1 and 2 below represent the boundary conditions of the surface electric field of any metal scatterer located in three-dimensional space for the EFIE electric field.

This is shown by applying.

는 전계(electric field)이며, 위첨자 i(incident)와 s(scattered)는 입사파와 산란파를 나타낸 것이다. 아래첨자 tan(tangential)은 산란체 표면과 수평인 방향 성분을 나타낸다.

는 magnetic vector potential이고, Φ는 scalar pontential 이다. ω=2πf₀이며 f₀는 중심주파수이다.here,

Is an electric field, and the superscripts i (incident) and s (scattered) represent incident and scattered waves. The subscript tan (tangential) indicates the direction component that is parallel to the scatterer surface.

Is the magnetic vector potential and Φ is the scalar potential. ω = 2πf ₀ and f ₀ is the center frequency.

은 관측점 영역 내 좌표값이며,

은 전원점 영역 내 좌표값이다. S는 적분영역을 나타낸다.ε, μ and σ are the permittivity, permeability and conductivity, respectively. Also,

Is the coordinate in the viewpoint area,

Is the coordinate value in the power point area. S represents an integration region.

Equations 3 and 4 show the definition and divergence characteristics of the RWG basis function, and Equation 5 expresses the current distribution flowing on the scatterer surface as the sum of the microcurrents using the RWG basis function.

,

각 삼각형 내부에 위치하는 벡터 성분이다.Equation 3 represents the definition of the RWG basis function. In the n-th triangle pair structure, Tn ⁺ is a triangle of the (+) region and Tn ⁻ is a triangle of the (−) region. In this case, An ⁺ and An ^- represent the width of each triangle. ln is the length of the intervening edge in the triangular pair structure.

,

Vector component located inside each triangle.

Equation 4 represents the divergence theorem of the RWG basis function as an equation.

는 산란체 표면을 흐르는 표면전류이며, 이때 표면전류는 RWG 기저함수로 분할된 요소들의 합으로서 나타낸다. a_n은 n-번째 RWG 기저함수의 크기를 나타내는 계수이다.here,

Is the surface current flowing through the scatterer surface, where surface current is expressed as the sum of the elements divided by the RWG basis function. a _n is a coefficient representing the size of the n-th RWG basis function.

In this case, the RWG basis function is defined as a vector phase through which the surface current of the scatterer divided into triangles numerically passes through the inner edge of the triangle pair structure.

Using Equation 5 and Equation 1 for the surface current represented by the RWG basis function, the EFIE is summarized as in Equation 6.

은 산란체 표면에 수평한 유닛벡터 성분이다. k는 전파상수(propagation constant)이다.

Is the unit vector component horizontal to the scatterer surface. k is a propagation constant.

이다.here,

to be.

Equation 6 may again be expressed as a determinant such as Equation 7.

Equation 7 has been described above. That is, if the inverse matrix [Z] ^-1 of the impedance matrix and the incident wave matrix [B] are used, the surface current matrix [I] can be calculated inversely.

Here, the element Z _mn in the impedance matrix [Z] can be expressed as Equation 8, and the element in the incident wave matrix [B] can be expressed as Equation 9. m and n are the column and row numbers of the matrix, respectively. Accordingly, Z _mn means an element corresponding to the m th row and the n th column in the impedance matrix [Z].

Zmn is the (m, n) th element of the impedance matrix and m, n is computed from the m-th (observation point area) and n-th (power point area) triangle elements, respectively.

b _m is the m-th element of the [B] matrix representing the electric field of the incident wave.

The configuration of Equations 1 to 9 corresponds to the facts apparent to those skilled in the art in the moment method technology to which the present invention belongs.

In the present invention, the Galerkin technique and the center point interpolation method are mixed with each other, and the center point interpolation method that simplifies the integral equation has been mentioned that there is a problem of inaccuracy. This accuracy increases as the power point area and the observation point area to be calculated get closer. In the present invention, it is necessary to define the distance (relative distance index) between the power point region and the observation point region where accuracy is maintained even using the center point interpolation method.

To this end, first, the first triangle pair consisting of two adjacent triangles among the plurality of triangles divided in step S210 is set as the power point region, and the second triangle pair is set as the observation point region through the region setting unit 120 ( S220).

FIG. 4 shows an example of the geometry of a triangular pair in step S220 of FIG. 2. It can be seen that both the power point region and the observation point region consist of triangle pairs. The power point area and the observation point area may be spaced apart from each other as shown in FIG. 4, or the two areas may be in contact with each other and adjacent to each other (for an embodiment in which the two areas contact each other, see FIG. 5 to be described later).

Next, the exponent calculator 130 calculates a relative distance index between the power point area and the observation point area using the vertex coordinates of the power point area and the observation point area (S230).

This step S230 will be described in detail with reference to FIG. 4 as follows. The exponent calculator 130 calculates a relative distance index between the center of gravity of one triangle in the power point area and the center of gravity of one triangle in the viewpoint area.

For convenience of description, the first triangle (T ^- _n ) and the second triangle (T ⁺ _n ) in the order of the two triangles in the power point area up and down, and the third triangle (in the order of the two triangles in the viewpoint area) T ^- _m ) and the fourth triangle (T ⁺ _m ). Where n and m represent triangle indices.

)과 제3 삼각형(T^- _m)의 무게 중심(

) 사이, 제1 삼각형(T^- _n)의 무게 중심(

)과 제4 삼각형(T⁺ _m)의 무게 중심(

) 사이, 제2 삼각형(T⁺ _n)의 무게 중심(

)과 제3 삼각형(T^- _m)의 무게 중심(

) 사이, 제2 삼각형(T⁺ _n)의 무게 중심(

)과 제4 삼각형(T⁺ _m)의 무게 중심(

) 사이의 상대거리 지수가 각각 연산된다. From this, a total of four relative distance indices are calculated. That is, the center of gravity of the first triangle T ^- _n

) And the center of gravity of the third triangle (T ^- _m )

Between), the center of gravity of the first triangle (T ^- _n )

) And the center of gravity of the fourth triangle (T ⁺ _m )

Between), the center of gravity of the second triangle (T ⁺ _n )

) And the center of gravity of the third triangle (T ^- _m )

Between), the center of gravity of the second triangle (T ⁺ _n )

) And the center of gravity of the fourth triangle (T ⁺ _m )

Relative distance exponents are calculated respectively.

All relative distance indices R _{adj _ cond} are calculated using Equation 10 below.

는 상기 전원점 영역 내의 한 개의 삼각형의 무게중심과 상기 관측점 영역 내의 한 개의 삼각형의 무게중심 사이의 거리를 나타낸다.Where the molecule

Denotes the distance between the center of gravity of one triangle in the power point area and the center of gravity of one triangle in the viewpoint area.

And L _s is the largest distance from the center of gravity of one triangle in the power point area to each vertex of two triangles in the power point area, and L _o is the distance from the center of gravity of one triangle in the viewpoint area. It means the largest distance of each vertex of two triangles in the viewpoint area. With such L _s L _o is simply defined by equation (11).

If max {} is a function of selecting the maximum value, max {a, b, c} = a when (a> b> c).

)과 제4 삼각형(T⁺ _m)의 무게 중심(

) 사이에 대한 상대거리 지수의 연산을 설명하는 예이다. 이러한 경우 L_s는 전원점 영역 내의 제1 삼각형(T^- _n)의 무게 중심(

)으로부터 전원점 영역 내의 두 삼각형의 각 꼭짓점(

,

,

,

) 사이의 거리 중 가장 큰 거리에 해당된다. L_o 또한 L_s와 동일한 원리로 연산되므로 상세한 설명은 생략한다.4 is a center of gravity of the first triangle (T ^- _n ) of the above four for convenience of description

) And the center of gravity of the fourth triangle (T ⁺ _m )

This is an example to explain the calculation of the relative distance index between. In this case, L _s is the center of gravity (T ^- _n ) of the first triangle T ^- _n in the power point region.

From each vertex of the two triangles within

,

,

,

) Is the largest distance between them. Since L _{o is} also calculated on the same principle as L _s , the detailed description is omitted.

FIG. 5 shows an example of another geometry in step S220 of FIG. 2. 5 is a case where the power point region and the observation point region are in contact with each other. In this case, the relative distance index will have a smaller size than in FIG. 4. 5 is an example for explaining the relative distance index between the second triangle and the third triangle.

4 and 5 are exponents of the relative distance between the power point region and the observation point region according to the position (m, n) of the element in the impedance matrix [Z], and each triangle for which integral calculation is performed It can be expressed as the ratio of the center of gravity distance between the elements and the longest distance between the center of gravity and each vertex. For example, when the relative distance index is 1, it means the maximum distance that the two calculation regions will be adjacent to each other.

After calculating the relative distance index, the comparison unit 140 compares the relative distance index with a predetermined threshold (value between 1 and 1.5) (S240).

Here, if the relative distance index is less than the threshold, the impedance calculating unit 150 calculates each element of the impedance matrix for the scatterer using a Galerkin technique (S250). Each element of the impedance matrix for the scatterer is calculated (S260).

Here, the equation for calculating the element value (Z _mn ) of the impedance matrix using the Galerkin technique is shown in Equation 8 above.

According to step S260, the greater the distance between the power point region and the observation point region, the greater the relative distance index and the higher the probability of using the center point interpolation method than the Galerkin method.

6 is a schematic diagram of an impedance matrix using an improved moment method according to an embodiment of the present invention. Elements Z ₁₁ , Z ₂₂ , Z ₃₃ ,... Located on the main diagonal in the impedance matrix of FIG. 6 are impedance values corresponding to the case where the relative distance index is zero.

Here, if the relative distance index R _{adj _ cond} is less than the threshold, the calculated corresponding element is arranged around the main diagonal line (a rectangular dotted line area in the diagonal direction) in the impedance matrix. In addition, when the relative distance index R _{adj _ cond} is larger than the threshold, the calculated element is arranged near the upper triangle or lower triangle in the impedance matrix (triangular dotted area in the upper and lower diagonal directions).

In other words, the impedance value calculated using the Galikin technique is located in the area around the main diagonal of the matrix, and the impedance value calculated by the center point interpolation is located near the upper or lower triangle of the matrix. In other words, the greater the degree of separation between the power point region and the viewpoint region, the more the calculated element is located away from the element on the main diagonal in the impedance matrix. In particular, the uppermost right side and the lowermost right side of the matrix represent the calculated impedance value for the case where the relative distance index is the largest. The greater the distance between the power point area and the viewpoint area, the smaller the impedance value and virtually no numerical integration is needed. In addition, the closer the distance between the power point region and the observation point region is, the larger the numerical integration value becomes, so that the center point interpolation method cannot be used.

As described above, the present invention determines a region requiring double integration of the Galerkin technique and a region to which the center point interpolation is to be applied using the relative distance index when calculating the impedance matrix [Z]. Here, the relative distance index can be set to a value of 1 or more to prevent the case where the area adjacent to the main diagonal is calculated by interpolation.

For the triangular region, the method of calculating impedance by the Galkin technique and the center point interpolation method is general in the art, but hereinafter, the concept thereof will be briefly described.

First, the basic contents of the numerical integration related to the Galerkin technique are explained. Coordinate transformation using the Simplex coordinate system is applied to perform the numerical integration of arbitrary triangular regions. The Simplex coordinate system is a coordinate system that expresses an arbitrary position with three triangle area ratios and three vertices divided around an arbitrary point located in an integrated region. It is possible to convert to a new rectangular coordinate system represented by α and β, and to perform a Gaussian quadrature, which is a numerical integration method, in the transformed coordinate system. Equation 12 converts an integral expression expressed in a rectangular coordinate system into a simplex coordinate system, and A and J (α, β) are Jacobian for triangle area and variable substitution, respectively.

는 임의의 함수 f(r)을 적분하기 위한 임의의 (x,y,z)좌표로 구성된 미소면적이다.(12) converts the coordinates of (x, y, z) into coordinates of (1-α-β, α, β).

Is a small area composed of arbitrary (x, y, z) coordinates for integrating arbitrary function f (r).

Equation 13 is a coordinate transformation equation of an arbitrary position vector using a simplex coordinate system.

,

,

는 삼각형을 이루는 세 꼭짓점 좌표이다.here,

,

,

Are the three vertex coordinates that make up a triangle.

Equation 14 is a Gaussian quadrature numerical integration method using a simplex coordinate system, and N represents the number of points for the integral operation.

은 임의의 적분함수이다. T는 적분영역이고, A는 삼각형의 넓이이다.Equation 14 represents a Gaussian quadrature numerical integration method using a simplex coordinate system.

Is any integral function. T is the integral region and A is the area of the triangle.

Equation 15 below represents a Green function.

은 3차원 자유공간의 Green 함수이다. T'은 전원점 내 적분 영역을 나타낸다. k₀는 자유공간의 전파상수이다.

Is the Green function of three-dimensional free space. T 'represents the integral region in the power point. k ₀ is the propagation constant of free space.

7 shows the numerical integration results for specific triangular regions. FIG. 7 illustrates a numerical integration result of the Green function using Equation 13. As the distance R between the power point rc 'and the observation point ro gets closer, the integration characteristic of the Green function is very unstable. It will have a property.

In addition, when performing integral calculations on pairs of triangular elements to calculate surface currents at specific edges, several analytical integration methods are used to avoid singularities for the same power point and observation point (self-term, R = 0). Should be used. The (1 / R) element of the Green function in the three-dimensional free space has a characteristic that the smaller the R value is, the more divergent it is, and a special technique is required to stabilize it. Among them, the most general method uses the simplified equations (16) and (17) using the Cauchy principal value and L'hopital's theorem.

Where k is a propagation constant and R is the distance between the power point and the observation point.

In the first term on the right side of Equation 16, R = 0 shows the same characteristics as Equation 17, and Gaussian quadrature numerical integration is possible even when R approaches zero. Therefore, the integral for the second term (1 / R) is calculated using the analytical method.

I ₁ and I ₂ separate only the integral terms containing singular solutions, and apply the singular solutions to integrate them.

Equation 18 and Equation 19 are integral expressions separated and separated from the EFIE by only the integral term for (1 / R), and may be expressed as integral terms I ₁ and I ₂ to which the Galerkin technique is applied.

On the other hand, when the power point and the observation point for the integration is located at a relatively far distance (R) relative to the integral region, it is generally possible to simplify the integral calculation by using an interpolation method. For this purpose, a center point interpolation method using a distance between center points and a triangular area is generally used.

In the present invention, the integral calculation of the adjacent area is performed by double integration using the Galerkin technique to maintain the calculation accuracy as much as possible, and by applying the center point interpolation method only for the calculation of the relatively long distance (R), the integral calculation error and the time required can be minimized. Can be. Equation 20 shows an integration method using the center point interpolation method.

는 임의의 적분함수이고, A_m은 m-번째 삼각형의 면적에 해당되고, T_m은 m번째 삼각형을 의미한다.

는 m번째 삼각형의 중심점 좌표를 의미한다.here

Is an arbitrary integration function, A _m corresponds to the area of the m-th triangle, and T _m means the m-th triangle.

Is the coordinate of the center point of the m-th triangle.

Equation 20 corresponds to Equation 21. This is the same as the impedance element value calculation formula shown in Equation (8). Therefore, Equation 20 may correspond to an impedance calculation equation when the center point interpolation method is used.

R _mn is the distance between the m-th and n-th triangles, and the superscript 'c' indicates the center of each triangle. R _mn ^c is the distance between the center points of the m-th and n-th triangles, and c _mn ^c is c _mn calculated from the center point.

When all the elements in the impedance matrix [Z] of FIG. 6 are obtained as described above, the inverse matrix calculator 160 calculates the inverse matrix [Z] ⁻¹ for the impedance matrix [Z] (S270). .

Then, the surface current calculator 170 calculates the surface current matrix [I] using the inverse matrix [Z] ^{− 1} of the impedance matrix and the incident wave matrix [B] (S280). Here, the method of obtaining the incident wave matrix [B] has been described above, which also corresponds to a conventional general operation process.

Using the surface current matrix [I], far field parameters (RCS, etc.) can be analyzed. The monostatic RCS is a radar cross section (RCS) when transmitting and receiving with one antenna, and the bistatic RCS is a radar cross section (when the receiving side is different from the transmitting side) when using two transmitting and receiving antennas.

8 shows the calculation speed and the error according to the threshold value change of the relative distance index when the embodiment of the present invention is applied.

As the threshold of the relative distance index increases, the computational time increases because the area required to use the Galerkin technique increases, while the calculation error shows that there is no significant change above the threshold of 1.5. Therefore, when the threshold value of the relative distance index is 1.5 or more, it can be confirmed that the calculation does not affect the error even if it is calculated using the center point interpolation method. In addition, it can be seen that when the threshold of the relative distance index is 7 or more, the calculation time converges to about 6 minutes. This means that the Galerkin technique is used in all areas. As a result, the threshold of the relative distance index is defined as 1.5, and the calculation time is reduced by three times or more when the improved moment method according to the present invention is used. have.

9 is a result of calculating the Bistatic RCS of the metal sphere with ka = 0.5 when using the existing moment method and the improved moment method according to the present invention.

Where ka (k is the propagation factor and a is the radius of the sphere) is the size of the metal. FIG. 9A illustrates a case where an existing moment method MoM is used, and FIG. 9B illustrates a case where the moment method iMoM of the present invention is used. The center frequency (f ₀ ) for numerical analysis is 10 GHz, with 750 edges of the scatterer surface. In this case, as the metal sphere, an equilateral triangle segmentation method using 12 reference points was used. In addition, an iterative method using the BiCJ-Stab (biconjugate gradient stabilized) method was used to calculate the inverse of the impedance matrix [Z] from the determinant of the moment method.

Equation 22 is used for analysis accuracy analysis using the Mie-series method and the moment method. Equation 22 expresses the rms (root-mean-square) value of 'norm of residuals' which is generally used to analyze the error of two data strings.

Is the data string of the radar cross section area RCS.

The calculation results of FIGS. 9A and 9B and the error analysis results using Equation 22 are 0.24dB (vertical polarization) /0.25dB (horizontal polarization) and 0.19dB / 0.21 in the conventional method and the method of the present invention, respectively. dB error was shown. The results are shown in Table 1.

Method of interpretation
Metal sphere [Z] matrix
size Error (dB) Time Conventional moment method
Galerkin Technique 0.24 / 0.25 7 minutes 55 seconds

750 × 750
Improved Moment Method of the Invention
(Gallkin Technique +
Center point interpolation)
0.19 / 0.21
3 minutes 20 seconds Remarks Reduced calculation time by about 50%

9 is a result of analyzing the object of the concrete structure of FIG. 3 through the embodiment of the present invention to which the improved moment method is applied, and the result of analysis by the conventional numerical method and the improved moment method proposed by the present invention. The results used are very consistent. In addition, it can be seen that the calculation time of the impedance matrix is shortened by less than half compared with the case of using only the existing Galerkin technique.

Hereinafter, the results of analyzing the monostatic RCS characteristics of the three- / forward wave reflector using the moment method of the present invention verified above. The size of the three-sided radio reflector is 4λ ₀ at the right side length ( l ), and the total number of internal corners is 1770. Omnidirectional propagation reflector length (l) is 1λ ₀ of radio wave reflector 20 on three sides and one as an array of information icosahedral geometry of the inner edge 3060 is the total atoms. Physical reference (PO) approximation and commercial software (HFSS) were used as reference data for monostatic RCS characterization of three- / forward-radiation reflectors. Expression

Is the radar cross-sectional area RCS, and θ and φ represent the incident angles of the radio waves. l is the length of the right side of the triangular reflector, and λ ₀ is the wavelength of the incident wave.

Figure 10 shows the results of the analysis of the monostatic RCS characteristics of the three- / forward wave reflector according to an embodiment of the present invention.

10 (a) to 10 (c) show characteristics of the three-sided radio reflector having a size 4λ ₀ relatively close to the PO approximation equation. In this case, PO approximation may be suitable for the analysis of relatively large radio reflectors of 10λ ₀ or more, but the improved moment method according to the present invention may be used to calculate the scattering characteristics having high complexity due to multiple reflection of smaller scatterers. It can be utilized.

10 (d) to (f), the omnidirectional radio reflector showed a consistent result compared to commercial software using a finite element method (FEM) algorithm, and exhibited an error within 0.6 dB at each peak. At this time, the analysis time depends on the performance of the PC (CPU 3.4GHz, RAM 3GB) applied in the experiment.

According to the present invention as described above, it is possible to reduce the calculation time of the scatterer in the three-dimensional space and the burden of hardware through the improved moment method using the Galerkin technique and the center point interpolation method. In addition, the computational speed of a single scatterer allows for the calculation of scatterers of complex structure simultaneously. The accuracy of the numerical analysis of the improved moment method using the present invention can be seen that the accuracy is exactly the same within the error range of about 0.2dB as compared with the Mie-series solution. When the present invention is applied to the research for analyzing the scattering radar backscattering characteristics of the scatterer considering the surface reflection characteristics by combining the image theory and the impedance surface theory, it is not easy to apply the numerical method to the range of remote sensing. Wide application will be possible.

While the present invention has been described with reference to exemplary embodiments, it is to be understood that the invention is not limited to the disclosed embodiments, but, on the contrary, is intended to cover various modifications and equivalent arrangements included within the spirit and scope of the appended claims. Accordingly, the true scope of the present invention should be determined by the technical idea of the appended claims.

100: electromagnetic wave numerical analysis device 110: division
120: region setting unit 130: exponential calculator
140: comparison unit 150: impedance calculation unit
160: inverse matrix calculation unit 170: surface current calculation unit

Claims (10)

Dividing the surface of the scatterer into a plurality of triangles using a Rao-Wilton-Glisson (RWG) basis function, and obtaining coordinates of each vertex of the triangle;
Setting a first triangle pair consisting of two adjacent triangles among the plurality of triangles as a power point region and a second triangle pair as a viewpoint point region;
Calculating a relative distance index between the power point area and the viewpoint area using each vertex coordinate of the power point area and the viewpoint area;
Comparing the relative distance index with a preset threshold; And
Calculating each element of the impedance matrix for the scatterer using the Galerkin method if the relative distance index is less than the threshold, and using the Center-point Interpolation Method if the relative distance index is less than the threshold. Numerical analysis method including electromagnetic waves. The method according to claim 1,
Computing the relative distance index,
Calculate a relative distance index between the center of gravity of one triangle in the power point area and the center of gravity of one triangle in the viewpoint area,
The relative distance index R _{adj_cond} is a numerical method of electromagnetic wave calculation through the following equation:

Where L _s is the largest distance from the center of gravity of one triangle in the power point area to each vertex of two triangles in the power point area, L _o is the center of gravity of one triangle in the viewpoint area The largest of the distances to each vertex of the two triangles within the viewpoint area,

Denotes the distance between the center of gravity of one triangle in the power point area and the center of gravity of one triangle in the viewpoint area. The method according to claim 2,
The power point area and the observation point area are adjacent to or spaced apart from each other,
The threshold is a value between 1 and 1.5,
If the relative distance index R _{adj _ cond} is less than the threshold, the elements are arranged around a major diagonal in the impedance matrix,
If the relative distance index R _{adj _ cond} is greater than the threshold, the corresponding element is arranged near an upper triangle or a lower triangle in the impedance matrix. The method according to claim 3,
The elements located on the main diagonal in the impedance matrix are impedance values corresponding to when the relative distance index is 0,
And the larger the distance between the power point region and the viewpoint region is, the calculated element is located at a position spaced apart from the element on the main diagonal in the impedance matrix. The method according to claim 1,
Calculating an inverse matrix for the impedance matrix; And
Computing a surface current matrix using the inverse matrix and the incident wave matrix of the impedance matrix. A division unit for dividing the surface of the scatterer into a plurality of triangles by using a Rao-Wilton-Glisson (RWG) basis function and obtaining coordinates of each vertex of the triangle;
An area setting unit configured to set a first triangle pair composed of two adjacent triangles among the plurality of triangles as a power point region and a second triangle pair as an observation point region;
An exponent calculator configured to calculate a relative distance index between the power point area and the observation point area using coordinates of each vertex of the power point area and the viewpoint point area;
A comparison unit comparing the relative distance index with a preset threshold; And
Impedance calculator for calculating each element of the impedance matrix for the scatterer using the Galerkin method if the relative distance index is less than the threshold, and using the Center-point Interpolation Method if the relative distance index is less than the threshold Electromagnetic wave numerical analysis device comprising a. The method of claim 6,
The exponential operator,
Calculate a relative distance index between the center of gravity of one triangle in the power point area and the center of gravity of one triangle in the viewpoint area,
The relative distance index R _{adj_cond} is an electromagnetic wave numerical analysis device that calculates through the following equation:

Where L _s is the largest distance from the center of gravity of one triangle in the power point area to each vertex of two triangles in the power point area, L _o is the center of gravity of one triangle in the viewpoint area The largest of the distances to each vertex of the two triangles within the viewpoint area,

Denotes the distance between the center of gravity of one triangle in the power point area and the center of gravity of one triangle in the viewpoint area. The method of claim 7,
The power point area and the observation point area are adjacent to or spaced apart from each other,
The threshold is a value between 1 and 1.5,
If the relative distance index R _{adj _ cond} is less than the threshold, the elements are arranged around a major diagonal in the impedance matrix,
And when the relative distance index R _{adj _ cond} is greater than the threshold, the corresponding element is arranged near an upper triangle or a lower triangle in the impedance matrix. The method according to claim 8,
The elements located on the main diagonal in the impedance matrix are impedance values corresponding to when the relative distance index is 0,
And the greater the distance between the power point region and the viewpoint region is, the calculated element is located at a position spaced apart from the element on the main diagonal in the impedance matrix. The method of claim 6,
An inverse matrix calculator configured to calculate an inverse matrix of the impedance matrix; And
And a surface current calculator configured to calculate a surface current matrix using the inverse matrix and the incident wave matrix of the impedance matrix.

Images