Energy characteristics of a nonlinear layer at resonant frequencies of wave scattering and generation

1 Introduction

Dielectric materials with a nonlinear, controllable permittivity are actively studied and practically applied in radiophysics, optics, electronics and instrumentation [1, 2, 3, 4, 5, 6, 7, 8, 9, 10]. For modern electronic equipment it is important to increase the polyfunctionality of devices by taking advantage of the properties of the nonlinear dielectric materials which they are made of. Therefore, the fabrication of modern devices of functional electronics requires a comprehensive familiarization with the features of these materials. The properties of solid, liquid or gaseous nonlinear media are important, allowing them to be used for transmutation of energy or information, detection, modulation, amplification, display, storage and further types of transformation of optical or electrical signals carrying information. In this regard, there is a need for new approaches and methods for creating of effective algorithms to solve problems of the analysis and the control of a physical process. The study of properties of resonant structures possessing an optimal generating ability is of high interest.

The present paper complements our previous works [11, 12, 13, 14, 15, 16]. In the setting of a coupled system approach, a numerical investigation of the energy characteristics of the radiation of the third harmonic by a nonlinear layer is performed on the basis of a block-wise iterative scheme for the computational analysis of the resonance properties of wave radiation by a nonlinear layered object in the vicinity of the eigenvalues of induced eigenvalue problems. The results of the calculations showed that as the strength of the incident field increases at a fixed excitation frequency, a shift in the resonant frequency of the nonlinear layer structure is observed. The shift of the excitation frequency away from the resonant frequencies (corresponding to the induced eigen frequencies) leads to a reduction in the generating abilities of the nonlinear object. When a resonant frequency (determined by the eigen frequency of the induced eigenvalue problem) was selected as the excitation frequency, it was possible to increase the portion of the generated energy significantly.

In this paper, the behaviour of energetic and computational peculiarities is studied with increasing excitation amplitude. Comparative results of computational schemes based on different quadrature rules are demonstrated, where the following aspects are evaluated: the convergence properties of the coupled system approach; the energy balance; the conditionalities of the matrices describing the induced eigen modes, and the processes of resonant wave radiation by the nonlinear structure. This made it possible to verify the degeneration of the matrices of the induced eigenvalue problems, and also to estimate the sensitiveness to the errors of the algorithm describing the material and energetic properties of the wave scattering and generation at the resonant frequencies of the non-linear layered structure.

2 The mathematical model

A nonlinear, nonmagnetic, isotropic, E-polarized

layered, plate-like dielectric object (see Figure 1) with a cubic nonlinear polarization vector

Figure 1

The plate-like structure with nonlinear layers

is investigated. The structure is irradiated by a packet of planar time-harmonic quasi-homogeneous oscillations, where the fulfillment of a so-called condition of phase synchronism is required.

The problem is considered in dimensionless space-time coordinates {x, y, z, t} for a harmonic time dependence exp(−inωt), n ∈ N = {1, 2, 3, . . . }. In these coordinates we have the following parameters: 4πδ denotes the height of the nonlinear layer, ω = κc is the circular frequency (dimensionless), κ is a frequency parameter (dimensionless) such that κ = ω c = 2 π λ (it characterizes the proportion between the real height h of the non-linear layer and the excitation wavelength λ, i.e. h λ = 2 κ δ ) ; c = (ε₀μ₀)^1/2 is the dimensionless speed of light in vacuum. The values of the real quantities {x^′ , y^′ , z^′ , t^′}, ω^′ are obtained by:

Here it is assumed that the incident plane wave packet consists of a highly intense field at a certain basic (or fundamental) frequency κ and of at most two weakly intense fields at the frequencies 2κ and 3κ. As a matter of fact, the highly intense field generates waves at the triple frequency. The weakly intense fields have only a controlling influence on the process of wave radiation. Effects of a static electric field are neglected, i.e. E(r, nκ|n₌₀) = 0. If the excitation has the form

the complex-valued amplitudes of the total scattering and generation fields can be represented in the form

where

denotes the Heaviside function;

are the longitudinal propagation constants;

are the excitation angles of the nonlinear object; “_T” denotes the operation of vector/matrix transposition, “^—” and “^_” indicate the excitation from the upper z > 2πδ and lower z < −2πδ half-spaces, respectively.

The quantities { Φ n κ } n = 1 3 at the multiple frequencies nκ are coupled by the condition of phase synchronism of the waves, see Figure 1 and [12, 13, 14, 15]:

The functions U(nκ; ·) ∈ L₂(−2πδ, 2πδ) that are being sought in (2) satisfy a system of nonlinear Hammerstein integral equations [11, 12, 13, 14, 15, 16]:

Here

are components of the incident fields (1), analogous to { a n κ i n c , b n κ i n c } n = 1 3 , the representation (2) of the total fields, with amplitudes

are the permittivity functions;

are the transverse propagation constants; χ x x ( 1 ) a n d χ x x x x ( 3 ) are the entries of the susceptibility tensors of the nonlinear object;

denotes the cubic susceptibility function of the nonlinear medium, where the Kleinman rule about the equality of the coefficients χ x x x x ( 3 ) at multiple frequencies is taken into account [5, 17]; and “^*” indicates the operation of complex conjugation.

The complex-valued amplitudes of the radiated waves (2) that are being sought have the following form:

The amplitudes { a n κ s c a t / g e n , b n κ s c a t / g e n } n = 1 3 of the fields (5) can be obtained from (3) by:

3 Solution of the coupled nonlinear equation system

The system of the nonlinear Hammerstein integral equations (3) can be solved numerically by means of a quadrature formula approach. As a result we arrive at a system of coupled nonlinear complex algebraic equations of the second kind:

where { { z l } l = 1 N : z 1 = − 2 π δ < … < z l < … < z N = 2 π δ } are the quadrature nodes,

are the amplitude vectors of the solution to be sought, I denotes the identity matrix, U ¯ n κ i n c a n d U _ n κ i n c are the amplitude vectors defined by the incident waves, B_nκ(U_κ , U_2κ , U_3κ) are matrix-valued nonlinear functions, C_κ(U_2κ , U_3κ), C_3κ(U_κ , U_2κ) are the right-hand side vectors determined by the choice of the quadrature formula.

The algorithm for the solution of (6) is realized by means of an iterative scheme, at each step of which the following system of linearized complex algebraic equations of second kind is solved [11, 12, 13, 14, 15, 16]:

where

gives the termination indices of the iterations within the blocks at each of the multiple frequencies nκ, n = 1, 2, 3, for the step number q ∈ N;

is the final step of the iterative process (7) for a prescribed value of the relative computing error ξ > 0.

4 The eigen modes of the induced eigenvalue problems

Applying the principle of analytic continuation, the non-linear problems (6) are continued into the complex space of values of the frequency parameter, thereby allowing an investigation of induced eigenvalue problems [12, 13, 14, 15, 16]. That is, we want to determine the eigenvalues and the nontrivial solutions of the homogeneous non-self-adjoint linear problems with the induced permittivities (4):

where Ω_nκ is the discrete set of eigenvalues that are being sought, H_nκ are two-sheeted Riemannian manifolds,

is the nontrivial vector of the solution that is being sought, corresponding to the eigenvalue κ_n ∈ Ω_nκ,

is the matrix defined by the vectors U_nκ (cf. (6)).

The electric components of the eigen fields can be represented by means of the (nontrivial) solutions of the eigenvalue problems (8) as follows:

Here κ = κ^inc is the basic frequency,

are the eigen field radiation amplitudes,

are the transverse propagation coefficients considered as functions of the complex-valued frequency spectral variables,

are given real values of the longitudinal propagation constants.

The eigen frequencies κ_n ∈ Ω_nκ ⊂ H_nκ of the problems (8) are found from the dispersion equations

by means of Newton’s method. The nontrivial solutions U_κn of the homogeneous systems (8)

are unique up to multiplicative constants and describe the eigen fields (9) of the induced eigenvalue problems.

5 The computational scheme for the resonant processes in the vicinity of eigen frequencies

The study of the behaviour of resonant wave radiation for the variation of the excitation amplitudes at frequencies in the vicinity of the induced eigen frequencies of the nonlinear layered object under consideration is executed by the help of the algorithm [16]:

The termination criterion (10) is the condition

Here are real, positive excitation frequencies, which are determined by the eigen frequencies κ_n ∈ Ω_nκ ⊂ H_nκ,

n = 1, 2, 3, of interest. In the last block of the j-th step of the algorithm (10) it is shown that such a frequency can be selected from one of the following quantities, see [16]:

In the scheme (10), the following steps are successively performed: the solution of the coupled problem (6) (by means of the scheme (7)); the solution of the induced eigenvalue problems (8); the selection of the excitation frequency from the values given in the last block of the algorithm (10).

In the investigations of the energetic processes of resonant wave radiation, due to the closeness of the eigen frequencies of the induced eigenvalue problems, the initial approximation to the excitation frequency κ = κ^inc := κ⁽⁰⁾ in the algorithm (10) can be obtained by solving the corresponding linear eigenvalue problems (8) with a trivial non-linear susceptibility function of the structure α(z) ≡ 0.

6 Numerical investigation of the electrodynamic properties of resonant wave scattering and generation

The energetic properties of the nonlinear layered object can be described using the following notation:

The quantities R n κ + , R n κ − characterize the proportions of the radiated energy to the excitation energy.

The complex-valued Fourier amplitudes of the total radiated fields (5) are denoted by

The fractions of the radiated energy at the frequencies nκ correspond to the quantities

The proportions of the energy of the third harmonic to the energy radiated by the first harmonic and to the total energy are defined as

In the numerical experiments, the Q-factors of the waves to the eigenvalues κ_n of the eigenvalue problems (8) at the frequencies nκ, i.e.

and the relative Q-factor of eigen waves

are of particular interest, see [12, 13, 14, 16, 18,].

Let us consider a nonlinear channelizing (α > 0) layer, see Figure1, with a permittivity given by the constant parameters

δ = 0.5, that is irradiated from above in the direction φ_nκ = 00, n = 1, 2, 3. The incident amplitudes are given by a κ i n c = c o n s t , a 2 κ i n c = a 3 κ i n c = 0 , b κ i n c = b 2 κ i n c = b 3 κ i n c = 0

Figures 2-4 show the results of a computational study of the electrodynamic properties of resonant waves radiated by the nonlinear layer, where the systems in (10) originated from the use of the composite Simpson’s rule and had the dimension 301.

Figure 2

The components of the eigenvalues κ_n = Reκ_n + iImκ_n, n = 1, 3; and the radiation frequencies κ, 3κ, see algorithm (10): (a) for κ = 0.375; (b) for κ = Re_κ1

Figure 3

The energetic properties and the relative Q-factors at κ = 0.375 and κ = Re_κ1: (a) scattering, (b) generation

Figure 4

The absolute amplitudes |U_nκ| and the permittivities Reε_nκ, Imε_nκ, n = 1,3, (a) at the frequency κ = R e κ 1 ( a κ i n c ) and (b) at the radiated frequency 3κ

The computations are based on the algorithm (10), where different variants for the determination of the real excitation frequencies nκ, being close to the spectral branches κ n ( a κ i n c ) ∈ Ω n κ ( a κ i n c ) ⊂ H n κ , n = 1 , 3 , of the induced dielectric structure are examined, too. Figure 2(a) depicts the case κ = R e κ 1 ( a κ i n c ) α ( z ) ≡ 0 ≈ 0.375 which lies in the vicinity of the eigenvalue κ₁ of the eigenvalue problem (8) for a linear layer (α(z) ≡ 0), whereas Figure 2(b) illustrates the case κ = R e κ 1 ( a κ i n c ) , w h e r e κ 1 ( a κ i n c ) is the computed value of the complex eigenvalue of the induced problem (8) for a nonlinear layer (α(z) ≠ 0). A generation effect at the frequency 3κ could be observed in both cases.

The curves R e κ n ( a κ i n c ) , I m κ n ( a κ i n c ) corresponding to the branches of the dispersion characteristics,

as well as the curves of the excitation (scattering) κ ( a κ i n c ) = R e κ 1 ( a κ i n c ) frequencies are shown in Figure 2.

In Figure 2 we see that as the amplitude of the incident field increases, a shift of the branches of eigen frequencies occurs. In the case of a fixed excitation frequency κ = 0.375, as shown in Figure 2(a), the complex branches of the dispersion curves κ n ( a κ i n c ) , n = 1 , 3 , move away from the excitation κ and generation 3κ frequencies. The case of the resonant frequency κ = Re κ 1 ( a κ i n c ) is depicted in Figure 2(b). In this situation, the branches of the complex spectra κ n ( a κ i n c ) , n = 1 , 3 , are located near the resonant frequencies of scattering κ ( a κ i n c ) and generation 3 κ ( a κ i n c ) .

In correspondence to the cases depicted in Figure 2, Figure 3 shows the energetic properties of the structure. From Figure 3(a) we can see that the shift of the eigen frequency R e κ 1 ( a κ i n c ) away from the excitation frequency κ = 0.375 (as elucidated above to Figure 2(a) results in a diminution of the resonance properties of the nonlinear layer, see the behaviour of the values R κ := 0.375 ± ( a κ i n c ) . Inspecting the energy characteristics in Figure 3(a), corresponding to the spectral characteristics in Figure 2(b) in the range of excitation amplitudes from 0.1 to 18.8 at the excitation frequencies κ = R e κ 1 ( a κ i n c ) , we notice that the nonlinear structure essentially maintains its resonance properties of transparency. The corresponding reflection coefficient is close to zero R κ := R e κ 1 + ( a κ i n c ) ≈ 0 , and the corresponding transmission coefficient R κ := R e κ 1 − ( a κ i n c ) slightly decreases monotonously from its maximum possible value ≈ 1 to 0.87381.

We emphasize the following effect. The reduction of transparency of the nonlinear layered object with growing amplitude of the excitation field at the resonant frequency of the induced layered structure (while preserving the non-reflective property) is due to the consumption of energy to generate the higher harmonics.

Thus, the growth of the excitation amplitude a κ i n c of the nonlinear layer (at κ = R e κ 1 ( a κ i n c ) of the incident field) leads to a decrease in the transmission coefficient R κ := R e κ 1 − ( a κ i n c ) , while the reflection coefficient remains close to ₊zero R κ := R e κ 1 + ( a κ i n c ) ≈ 0 , see Figure 3(a). Such a reduction of transparency at κ = R e κ 1 ( a κ i n c ) is caused by the energy exchange from the excitation energy to the generated energy. Among others, the monotonous declining of the relative Q-factor Q { 1 : 3 ; κ := R e κ 1 } ( a κ i n c ) is an indicator of this effect, see Figure 3(a).

The curves in Figure 3(b) demonstrate that the energetic properties of the structure, at the resonant frequency κ = Re κ 1 ( a κ inc ) , are substantially (actually many times) better with growing amplitude a κ i n c in the range of generation than for the fixed frequency κ = 0.375. This is supported by a comparison of the graphs from Figure 3(b) for the quantities

as well as the rates of decrease of the functions

for κ = R e κ 1 ( a κ i n c ) and κ = 0.375, resp.

We emphasize that an indirect analysis of the energy exchange processes occurring during the wave radiation by nonlinear objects is possible, namely by considering the magnitude and rate of decrease of the relative Q-factor Q { 1 : 3 ; κ } ( a κ i n c ) of the oscillations. For example, comparing the graphs of Q { 1 : 3 ; κ } ( a κ i n c ) with the energy characteristics represented by the curves R 3 κ ± ( a κ i n c ) , W { 3 : 1 ; κ } ( a κ i n c ) , W { 3 : 123 ; κ } ( a κ i n c ) , f o r κ = Re κ 1 ( a κ i n c ) a n d κ = 0.375 in Figure 3(b) we see the following. The local decrease in the value of the relative quality factor Q { 1 : 3 ; κ } ( a κ i n c ) leads to a surge of energy radiated in the third harmonic. The higher the rate of decrease of Q { 1 : 3 ; κ } ( a κ i n c ) , the greater the portion of the generated energy.

Figure 4 illustrates the amplitude characteristics |U_nκ| of the scattered/generated fields (5), where

as well as the values of the permittivity (4):

corresponding to the resonant frequencies κ = R e κ 1 ( a κ i n c ) and 3κ. The results are given for the range a κ i n c × z ∈ [0.1, 18.8] × [−2πδ, 2πδ].

The scattered/generated (5) and eigen fields (9) of the dielectric layer in the case of E-polarization under consideration is classified by the type

see [11, 12, 13, 14, 15, 16, 18]. The indices specify the number of local maxima |E_x| (or |U|, since |U| = |E_x|, see (5), (9)) along the coordinate axes in the space filled by the dielectric layer, see Figure 1. Since the oscillations under investigation are

homogeneous along the x-axis and quasi-homogeneous along the y-axis, they are of the type

where the index p indicates the number of local maxima of the function |U| with respect to z on the segment [−2πδ, 2πδ].

In the studied range of amplitudes of the incident field, we observe scattered fields of the type H_0,0,4 (see |U_κ| in Figure 4(a) and generated fields of the type H_0,0,10 (see |U_3κ| in Figure 4(b)).

The nonlinear components ε n κ ( N L ) of the permittivity ε_nκ at each frequency κ = Re κ 1 ( a κ i n c ) and 3κ are determined by the values U(κ; z) and U(3κ, z) of the fields, see (4) and Figure 4. For the nonabsorbing medium under consideration Imε^(L)(z) ≡ 0, taking into account the fact that the susceptibility function α(z) is real-valued, the equality Imε_nκ = I m ε n κ ( N L ) holds, see (4). The increase in the amplitude a κ i n c of the incident field at the frequency κ = Re κ 1 ( a κ i n c ) leads to the generation of a field U(3κ; z). In the considered case, the quantity Imε_κ oscillates and assumes both negative as well as positive values in the layer, see Figure 4(a)). The outlined situation is typical for the expenses (losses) of energy in the nonlinear object (at the excitation frequency κ = R e κ 1 ( a κ i n c ) )that are consumed for the generation of the third harmonic field (at the frequency 3κ) [11, 12, 13, 14, 15, 16]. The generated fields at the frequency 3κ are weak and energetically do not contribute to the generation of new harmonics. Note that Im_ε3κ ≡ 0, see Figure 4(b).

7 Computational specifics of the study

The numerical results for the energetic and material characteristics of resonant wave radiation at multiple frequencies by nonlinear media near the eigen frequencies of scattering κ = R e κ 1 ( a κ i n c ) and generation 3κ are obtained by solving systems of coupled nonlinear algebraic equations on the basis of the algorithm (10). The matrix elements of the problems under investigation are obtained by applying quadrature formulas to systems of nonlinear Hammerstein integral equations in analyzing the processes of wave scattering and generation, as well as to linearized integral equations in the study of the eigen modes of layered objects with induced permittivities.

The calculations are carried out with double precision in the standard arithmetic of complex numbers, the corresponding machine epsilon is equal to

That is, the number of significant decimal digits in the mantissa is not higher than 16.

This section compares computational aspects of some variants of the process (10), using the composite Simpson’s quadrature rule in one case and the higher-accurate composite Weddle’s rule in the other case. In the first case, equation systems of dimension 301 result, whereas in the second case the systems have the dimension 121, see [19] and Figure 5-6.

Figure 5

Examination of the iterative schemes for κ = Re_κ1: (a) The violation W^(Error) of the energy balance; (b) The number Q of iterations in the coupled system, and the curves of the resonant scattering frequencies κ = R e κ 1 ( a κ i n c )

In the range of parameters under investigation, the energy balance ratio

Figure 6

The conditionalities of the matrices at κ = Re_κ1 for: (a) the degeneration of I − B_nκ(κ_n) in the points κ_n ∈ Ω_nκ, n = 1, 3, of the spectrum; (b) the sensitiveness of I − B_{n κ} (U_κ , U2_κ , U3_κ) to errors

for nonabsorbing media I m ε n κ ( L ) ( y , z ) ≡ 0 was computationally checked with a relative iteration error less than 10⁻⁷. The value of the energy balance error

was calculated.

In both of the considered cases of realization of the iteration process (10) with κ = Re_κ1, based on the use of Simpson’s and Weddle’s quadrature rules, the absolute computing error is of acceptable size: |W^(Error)| < 6.1·10⁻⁹, see Figure 5(a). In this, a difference in the number of iterations Q for the step size Δ a κ i n c = 0.1 in the block of the coupled system of the scattering and generation processes is observed only at the final stages of the computations for a κ i n c ∈ [ 15 , 16.9 ] , see Figure 5(b). In addition, the computed values of the resonant scattering frequencies κ = Re κ 1 ( a κ ( i n c ) ) f o r a κ ( i n c ) ∈ [ 0.1 , 16.9 ] (see the curves κ = Re_κ1 [Simpson] and κ = Re_κ1 [Weddle] in Figure 5(b) coinciding up to the graphic accuracy) correspond to the depiction of the equivalent processes shown in Figures 2(b), 3 and 4.

Figure 6 shows estimates of the conditionality of the matrices from (10) obtained by Simpson’s and Weddle’s quadrature rules for linear systems of dimensions 301 and 121, resp. The estimation of the conditionalities of the matrices was carried out using the technique described in [20, 21].

It is possible to characterize the level of degeneration of the matrices I − B_nκ(κ_n) at the computed eigenvalues κ_n ∈ Ω_nκ ⊂ H_nκ, n = 1, 3, see (10) and (8), by

These values reflect large matrix conditionalities, which is an indirect feature of the degeneration (or singularity) of these matrices, see Figure 6(a).

The sensitiveness of the matrices I−B_nκ in the coupled system, see (10) and (6), to errors in the vicinity of the resonant radiation frequencies κ = Re_κ1, 3κ can be described by the help of the values

The logarithm of the matrix conditionality gives some information about how many significant digits are lost in solving the system [21].

After comparison of the experimental results depicted in Figure (6b) by the curves Y_nκ [Simpson] (for systems of dimension 301) and Y_nκ [Weddle] (for systems of dimension 121) for n = 1, 3, we can summarize the following observations. The computational method (10), realized with Weddle’s quadrature rule is more sensitive with respect to computing errors than the same method, but constructed using Simpson’s quadrature rule. In this case, a diminution in the dimension of the system from301 to 121 leads to a deficit of at most one significant digit in the solution, see the difference of the quantities Y_nκ [Weddle] and Y_nκ [Simpson], for n = 1 and n = 3 in Figure 6(b). Since the reduction in the dimension of the equation systems significantly reduces the computing time, this is a significant advantage and opens the way to solving new more complex nonlinear resonant wave scattering and generation problems for both layered as well as layered periodic structures.

The computational algorithm (10) developed in the present work is based on the rigorous solution of systems of coupled nonlinear boundary value problems and the corresponding induced spectral homogeneous boundary value problems of Electrodynamics. The spectral problems, in turn, allow to adjust precisely the resonant frequencies of excitation of the nonlinear objects and, in fact, to control their radiating abilities.

The uniform approach described in [14, 15] in the develeopment of mathematical models for the study of layered and layered-periodic problems of Electrodynamics allows to apply the given computational schemes to the analysis of periodic (grating) problems.

While usually the periodicity of a grating is predetermined by a corresponding periodicity property of the linear component of the permittivity (i.e. by an inherent material property), here the periodicity of the non-linear structure under investigation can additionally be controlled by the incident field [14, 15]. Ultimately, this means that a grating can be created, even if the linear component of the permittivity is not periodic.

In particular, if the nonlinear layer is irradiated by a quasi-homogeneous excitation field, it is transformed into a transversely layered structure. If the same nonlinear layer is excited by a quasi-periodic field, it becomes a layered-periodic (i.e. a transversely layered and longitudinally periodic) grating structure.

8 Conclusion

In the setting of a coupled system approach and taking into special account the behaviour of the complex spectra, a block-iterative scheme for the numerical exploration of the resonance properties of radiated multi-frequency oscillations by a nonlinear object in the vicinity of the eigenvalues of induced eigenvalue problems is discussed.

The results of numerical investigations of the properties of resonant wave radiation at multiple frequencies by a nonlinear object are described. It could be observed that the effect of the transparency diminution of a nonreflecting nonlinear layered structure at the resonant excitation is a consequence of the energy exchange that is expended for the generation of higher harmonics.

It has been demonstrated that it is possible to study the energy exchange processes occurring during multi-frequency wave radiation by nonlinear objects indirectly by means of the spectral properties of those objects induced by the incident field.

A comparison of the computational features, resulting from the use of the composite Simpson’s and Weddle’s quadrature rules (the latter with a higher degree of accuracy) is presented. Results for the energy balance, the convergence rate of the coupled system approach with increasing excitation amplitude are discussed. Estimations of the conditionalities of the matrices describing both the induced eigen modes and the resonant multi-frequency wave radiation processes by a nonlinear layered structure are given. This made it possible to obtain an indirect criterion for testing the degeneration of the matrices of the eigenvalue problems at the points of the spectrum, and also to estimate the sensitiveness to computing errors of the algorithm or the investigation of the processes of wave radiated near resonant frequencies.

It is shown that the use of quadrature formulas of a higher degree of accuracy allows to decrease the dimension of the equation systems to be solved by about 2.5, which saves considerable computing time. This is a significant advantage in solving new nonlinear problems of resonant wave radiation for both layered as well as layered periodic structures. It opens the prospects for the development of new numerical-analytical approaches to the investigation of irregular properties of wave radiation.