Computing the Electric and Magnetic Matrix Green’s Functions in a Rectangular Parallelepiped with a Perfect Conducting Boundary

Abstract

A method for the approximate computation of frequency-dependent magnetic and electric matrix Green’s functions in a rectangular parallelepiped with a perfect conducting boundary is suggested in the paper. This method is based on approximation (regularization) of the Dirac delta function and its derivatives, which appear in the differential equations for magnetic and electric Green’s functions, and the Fourier series expansion meta-approach for solving the elliptic boundary value problems. The elements of approximate Green’s functions are found explicitly in the form of the Fourier series with a finite number of terms. The convergence analysis for finding the number of the terms is given. The computational experiments have confirmed the robustness of the method.

1. Introduction

In recent years a lot of attention has been devoted to derivation of the electric and magnetic fields inside bounded domains with perfect conducting boundaries [1–4]. The study of the distribution of electromagnetic fields inside a real indoor environment (cabinets, desks, etc.) has the great practical interest (see, e.g., [5]). The theory (existence, uniqueness, and stability estimate theorems) of time-dependent Maxwell’s equations in bounded domains with perfect conducting boundary conditions and smooth current and charge densities has been described by Dautray and Lions in [6]. These theoretical results form the background for many computational methods such as finite element method, finite-difference method, finite-difference time domain method, modal expansion method, and Nitsche type method, which have been handled for computation of the electric and magnetic fields in bounded domains with a perfect conducting boundary [1–3, 6–8].

We need to note that the solution (both analytical and numerical) of the initial value and initial boundary value problems in electrodynamics is often facilitated by introduction of electric and magnetic Green’s functions [9–17], but the derivations of electric and magnetic Green’s functions have been made in a free space or unbounded domain. The computation of time-dependent electric and magnetic Green’s functions for the bounded domain with perfect conducting conditions has not been achieved so far.

In our paper we suggest a new method for the approximate computation of frequency-dependent magnetic and electric matrix Green’s functions in a rectangular parallelepiped with perfect conducting boundary. This method consists of the following. The equations for magnetic Green’s function are written in the special form which does not contain elements of electric Green’s function. These equations are elliptic partial differential equations and contain the components of the generalized vector functions (distributions) curl_𝑥[𝛿(𝑥−𝑥⁰)𝐞^𝑠], 𝑠=1,2,3 in their right hand sides. We approximate components of these generalized vector functions by the Fourier series with a finite number of terms and then applied standard Fourier series expansion meta-approach for solving the boundary value problem for the elliptic partial differential equations in the bounded region. After calculation of all columns of magnetic matrix Green’s function the approximate computation of the columns of electric Green’s function is given by explicit calculation of the differential operator curl_𝑥 from each column of magnetic matrix Green’s function with respect to 3D space variable 𝑥=(𝑥₁,𝑥₂,𝑥₃).

2. Frequency-Dependent Magnetic and Electric Green’s Functions

Let 𝑏₁,𝑏₂,𝑏₃ be given positive numbers and let 𝑉=𝑥∈𝐑³∶0<𝑥₁<𝑏₁,0<𝑥₂<𝑏₂,0<𝑥₃<𝑏₃ (1) be a rectangular parallelepiped; Γ is the boundary of 𝑉; 𝜔∈𝐑 (frequency, 𝜔>0), let 𝑥⁰=(𝑥⁰₁,𝑥⁰₂,𝑥⁰₃)∈𝑉 be 1D and 3D parameters, respectively. The 3×3 matrices 𝐆^𝐻𝑥,𝜔,𝑥⁰=𝐻¹₁𝑥,𝜔,𝑥⁰𝐻²₁𝑥,𝜔,𝑥⁰𝐻³₁𝑥,𝜔,𝑥⁰𝐻¹₂𝑥,𝜔,𝑥⁰𝐻²₂𝑥,𝜔,𝑥⁰𝐻³₂𝑥,𝜔,𝑥⁰𝐻¹₃𝑥,𝜔,𝑥⁰𝐻²₃𝑥,𝜔,𝑥⁰𝐻³₃𝑥,𝜔,𝑥⁰,𝐆^𝐸𝑥,𝜔,𝑥⁰=𝐸¹₁𝑥,𝜔,𝑥⁰𝐸²₁𝑥,𝜔,𝑥⁰𝐸³₁𝑥,𝜔,𝑥⁰𝐸¹₂𝑥,𝜔,𝑥⁰𝐸²₂𝑥,𝜔,𝑥⁰𝐸³₂𝑥,𝜔,𝑥⁰𝐸¹₃𝑥,𝜔,𝑥⁰𝐸²₃𝑥,𝜔,𝑥⁰𝐸³₃𝑥,𝜔,𝑥⁰, (2) whose columns 𝐇^𝑠(𝑥,𝜔,𝑥⁰), 𝐄^𝑠(𝑥,𝜔,𝑥⁰), and 𝑠=1,2,3 satisfy curl_𝑥𝐇^𝑠𝑥,𝜔,𝑥⁰=−𝑖𝜔𝜖𝐄^𝑠𝑥,𝜔,𝑥⁰+𝛿𝑥−𝑥⁰𝐞^𝑠,𝑥∈𝑉,curl_𝑥𝐄^𝑠𝑥,𝜔,𝑥⁰=𝑖𝜇𝜔𝐇^𝑠𝑥,𝜔,𝑥⁰,𝐄^𝑠×𝐧^→_Γ=0,𝐇^𝑠⋅𝐧^→_Γ=0,𝑠=1,2,3, (3)(4)(5) are called frequency-dependent magnetic and electric Green’s functions in 𝑉, respectively. Here 𝑖 is the imaginary unit (𝑖²=−1); 𝜖>0 and 𝜇>0 are given constants characterizing the permittivity and permeability of the electromagnetic medium, which is located inside 𝑉; 𝐞^𝑠 are the basis vectors from Euclidean space 𝐑³ (𝐞¹=(1,0,0)^𝑇,𝐞²=(0,1,0)^𝑇,𝐞³=(0,0,1)^𝑇); 𝑥⁰=(𝑥⁰₁,𝑥⁰₂,𝑥⁰₃)∈𝑉 is the 3D parameter; 𝛿(𝑥−𝑥⁰)=𝛿(𝑥₁−𝑥⁰₁)𝛿(𝑥₂−𝑥⁰₂)𝛿(𝑥₃−𝑥⁰₃), and 𝛿(𝑥_𝑗−𝑥⁰_𝑗) are the Dirac delta distribution concentrated at 𝑥_𝑗=𝑥⁰_𝑗 and 𝑗=1,2,3; 𝐧^→ is an outward unit normal vector on Γ. The relations in (5) are called the perfect conducting boundary conditions.

Remark 1. We note that div_𝑥𝐇^𝑠(𝑥,𝜔,𝑥⁰)=0 for all 𝑥∈𝑉 is a consequence of (4).

Applying curl_𝑥 to (3) and using (4), we find 𝜔²𝜖𝜇𝐇^𝑠−curl_𝑥curl_𝑥𝐇^𝑠=−curl_𝑥𝛿𝑥−𝑥⁰𝐞^𝑠,𝑥∈𝑉. (6) Taking the cross product of both sides of (3) and 𝐧^→, and using (5) we get curl_𝑥𝐇^𝑠×𝐧^→_Γ=0,(𝐇^𝑠⋅𝐧^→)_Γ=0. (7)

The following equation follows from (3) 𝐄𝑥,𝜔,𝑥⁰=𝑖𝜔𝜖curl_𝑥𝐇^𝑠𝑥,𝜔,𝑥⁰−𝛿𝑥−𝑥⁰𝐞^𝑠. (8)

Remark 2. If 𝐇^𝑠 is a vector function satisfying (6) and (7) and 𝐄^𝑠 is a vector function satisfying (8) for a fixed value of 𝑠=1,2,3, then 𝐇^𝑠 and 𝐄^𝑠 satisfy (3)–(5); that is, 𝐇^𝑠 and 𝐄^𝑠 are columns of frequency-dependent magnetic and electric Green’s functions. Really, from (8), we find curl_𝑥𝐇^𝑠(𝑥,𝜔,𝑥⁰)=−𝑖𝜔𝜖𝐄^𝑠(𝑥,𝜔,𝑥⁰)+𝛿(𝑥−𝑥⁰)𝐞^𝑠. Substituting the expression for curl_𝑥𝐇^𝑠 into (6) we have curl_𝑥𝐄^𝑠(𝑥,𝜔,𝑥⁰)=𝑖𝜇𝜔𝐇^𝑠(𝑥,𝜔,𝑥⁰). Moreover, taking the cross product of (8) and 𝐧^→, letting 𝑥∈Γ and, using (7) we find that 𝐄^𝑠 satisfies the first condition of (5). Therefore, we obtain that if 𝐇^𝑠 and 𝐄^𝑠 satisfy (6), (7), and (8) then 𝐇^𝑠 and 𝐄^𝑠 satisfy (3)–(5) also; that is, 𝐇^𝑠 and 𝐄^𝑠 are 𝑠-columns of frequency-dependent electric and magnetic Green’s functions.

Remark 3. Let 𝐇^𝑠 and 𝐄^𝑠 satisfy (6) and (8). Then div_𝑥𝐇^𝑠=0, 𝑥∈𝑉.

3. The Eigenvalue-Eigenfunction Problem for the Operator Δ_𝑥I₃ in a Rectangular Parallelepiped

Let us consider the problem of finding all values of 𝜆 (eigenvalues) for which there exists a nonzero vector-valued function 𝐔(𝑥)=(𝑈₁(𝑥),𝑈₂(𝑥),𝑈₃(𝑥))^𝑇 (vector-valued eigenfunction) satisfying Δ_𝑥𝐈₃𝐔(𝑥)+𝜆𝐔(𝑥)=0,𝑥∈𝑉,curl_𝑥𝐔(𝑥)×𝐧^→_Γ=0,𝐔⋅𝐧^→_Γ=0, (9)(10) where 𝐈₃ is the 3×3 identity matrix.

In this section we show that the eigenvalue-eigenfunction problem (9) and (10) can be decomposed into three eigenvalue-eigenfunction problems for the scalar-valued Laplace operator in the rectangular parallelepiped with the mixed boundary conditions. These subproblems can be solved by a standard technique [18].

Lemma 4. Let the boundary Γ of the parallelepiped 𝑉 be presented as Γ=Γ₁∪Γ₂∪Γ₃, where Γ₁, Γ₂, and Γ₃ are the sets defined by Γ₁=𝑥=𝑥₁,𝑥₂,𝑥₃∶𝑥₁=0∨𝑥₁=𝑏₁,0<𝑥₂<𝑏₂,0<𝑥₃<𝑏₃,Γ₂=𝑥=𝑥₁,𝑥₂,𝑥₃∶𝑥₂=0∨𝑥₂=𝑏₂,0<𝑥₁<𝑏₁,0<𝑥₃<𝑏₃,Γ₃=𝑥=𝑥₁,𝑥₂,𝑥₃∶𝑥₃=0∨𝑥₃=𝑏₃,0<𝑥₁<𝑏₁,0<𝑥₂<𝑏₂. (11) The boundary condition in (10) can be written as 𝛼_𝑗(𝑥)𝜕𝑈_𝑗(𝑥)𝜕𝑛+𝛽_𝑗(𝑥)𝑈_𝑗(𝑥)_Γ=0, (12) where the functions 𝛼_𝑗(𝑥) and 𝛽_𝑗(𝑥), 𝑗=1,2,3 are defined by 𝛼_𝑗(𝑥)=0,𝑥∈Γ_𝑗,1,𝑥∈Γ⧵Γ_𝑗;𝛽_𝑗(𝑥)=1,𝑥∈Γ_𝑗,0,𝑥∈Γ⧵Γ_𝑗. (13)

Proof. Let us prove this for 𝑗=1 and let us consider the part of the boundary where 𝑥₁=𝑏₁. Thus, the outward normal unit vector is 𝐞¹=(1,0,0)^𝑇. The second condition of (10) reads 𝑈₁𝑏₁,𝑥₂,𝑥₃=0,for0<𝑥₂<𝑏₂,0<𝑥₃<𝑏₃. (14) The first condition of (10) takes the form 𝜕𝑈₂𝜕𝑥₁𝑏₁,𝑥₂,𝑥₃−𝜕𝜕𝑥₂𝑈₁𝑏₁,𝑥₂,𝑥₃=0,𝜕𝑈₃𝜕𝑥₁𝑏₁,𝑥₂,𝑥₃−𝜕𝜕𝑥₃𝑈₁𝑏₁,𝑥₂,𝑥₃=0. (15) Combining the last three equalities we can write the boundary condition (10) in the form 𝛼_𝑗(𝑥)𝜕𝑈_𝑗(𝑥)𝜕𝑛+𝛽_𝑗(𝑥)𝑈_𝑗(𝑥)_𝑥1=𝑏1=0,𝑗=1,2,3, (16) where 𝛼₁=0, 𝛽₁=1; 𝛼₂=1, and 𝛽₂=0; 𝛼₃=1 and 𝛽₃=0. Similar arguments hold for the lower part of Γ₁ and for the parts of the boundaries Γ₂ and Γ₃ leading to (12).

Let 𝜆_𝑘𝑚𝑛 denote the eigenvalues and 𝐔^𝑘𝑚𝑛(𝑥) the corresponding vector-valued eigenfunctions of (9) and (10). The eigenfunctions 𝐔^𝑘𝑚𝑛(𝑥) can be written in the form 𝐔^𝑘𝑚𝑛(𝑥)=𝑈₁^𝑘𝑚𝑛(𝑥)𝐞¹+𝑈₂^𝑘𝑚𝑛(𝑥)𝐞²+𝑈₃^𝑘𝑚𝑛(𝑥)𝐞³, (17) and 𝜆_𝑘𝑚𝑛, 𝑈_𝑗^𝑘𝑚𝑛(𝑥), and 𝑗=1,2,3 satisfy Δ_𝑥𝑈_𝑗^𝑘𝑚𝑛(𝑥)+𝜆_𝑘𝑚𝑛𝑈_𝑗^𝑘𝑚𝑛(𝑥)=0,𝑥∈𝑉;𝛼_𝑗(𝑥)𝜕𝑈_𝑗^𝑘𝑚𝑛(𝑥)𝜕𝑛+𝛽_𝑗(𝑥)𝑈_𝑗^𝑘𝑚𝑛(𝑥)_𝑥1=𝑏1=0. (18)(19) Applying the standard technique of the spectral theory (see, e.g., [18]) we find all eigenvalues and corresponding eigenfunctions of the problem (18) and (19) for each 𝑗 as follows.

For 𝑗=1, 𝜆_𝑘𝑚𝑛=𝑘𝜋𝑏₁²+𝑚𝜋𝑏₂²+𝑛𝜋𝑏₃²,𝑈₁^𝑘𝑚𝑛(𝑥)=𝐴_𝑘𝑚𝑛sin𝑘𝜋𝑏₁𝑥₁cos𝑚𝜋𝑏₂𝑥₂cos𝑛𝜋𝑏₃𝑥₃, (20)𝑘=1,2,…;  𝑛,𝑚=0,1,…  are eigenvalues and corresponding eigenfunctions.

For 𝑗=2, 𝜆_𝑘𝑚𝑛=𝑘𝜋𝑏₁²+𝑚𝜋𝑏₂²+𝑛𝜋𝑏₃²,𝑈₂^𝑘𝑚𝑛(𝑥)=𝐵_𝑘𝑚𝑛cos𝑘𝜋𝑏₁𝑥₁sin𝑚𝜋𝑏₂𝑥₂cos𝑛𝜋𝑏₃𝑥₃, (21)𝑚=1,2,…; 𝑘,𝑛=0,1,… are eigenvalues and corresponding eigenfunctions.

For 𝑗=3, 𝜆_𝑘𝑚𝑛=𝑘𝜋𝑏₁²+𝑚𝜋𝑏₂²+𝑛𝜋𝑏₃²,𝑈₃^𝑘𝑚𝑛(𝑥)=𝐶_𝑘𝑚𝑛cos𝑘𝜋𝑏₁𝑥₁cos𝑚𝜋𝑏₂𝑥₂sin𝑛𝜋𝑏₃𝑥₃, (22)𝑛=1,2,…; 𝑘,𝑚=0,1,… are eigenvalues and corresponding eigenfunctions.

Here 𝐴_𝑘𝑚𝑛, 𝐵_𝑘𝑚𝑛, and 𝐶_𝑘𝑚𝑛 are constants, such that ^𝑏3₀^𝑏2₀^𝑏1₀𝑈_𝑗^𝑘𝑚𝑛(𝑥)⋅𝑈^{𝑘′𝑚′𝑛′}_𝑗(𝑥)=()≠𝑘𝑑𝑥0,if𝑘,𝑚,𝑛^′,𝑚^′,𝑛^′1,if𝑘=𝑘^′,𝑚=𝑚^′,𝑛=𝑛^′. (23)𝑘,𝑚,𝑛=0,1,…; 𝑗=1,2,3. That is, 𝐴_𝑘𝑚𝑛=𝐵_𝑘𝑚𝑛=𝐶_𝑘𝑚𝑛=8/𝑏₁𝑏₂𝑏₃, 𝐴_𝑘0𝑛=𝐴_𝑘𝑚0=𝐵_0𝑚𝑛=𝐵_𝑘𝑚0=𝐶_0𝑚𝑛=𝐶_𝑘0𝑛=4/𝑏₁𝑏₂𝑏₃ 𝐴_𝑘00=𝐵_0𝑚0=𝐶_00𝑛=2/𝑏₁𝑏₂𝑏₃, and 𝑘=1,2,…; 𝑚=1,2,…; 𝑛=1,2,….

Remark 5. Let 𝑈_𝑗^𝑘𝑚𝑛(𝑥), 𝑗=1,2,3, be eigenfunctions of (18) and (19) defined above. By the direct computation we can show that 𝐔^𝑘𝑚𝑛(𝑥)=(𝑈₁^𝑘𝑚𝑛(𝑥),𝑈₂^𝑘𝑚𝑛(𝑥),𝑈₃^𝑘𝑚𝑛(𝑥))^𝑇 satisfies the boundary conditions (10) for the rectangular parallelepiped 𝑉.

Remark 6. For each 𝑗=1,2,3 the system of eigenfunctions {𝑈_𝑗^𝑘𝑚𝑛(𝑥)} is orthogonal and complete in the space of the square integrable functions over 𝑉. That is, any scalar-valued function ℎ(𝑥)∈𝐶^∞₀(𝑉) can be written in the form of uniformly convergent series over 𝑉 (see, e.g., [18]) ℎ(𝑥)=^∞∞_𝑘=0^∞_𝑚=0𝑛=0ℎ_𝑗^𝑘𝑚𝑛𝑈_𝑗^𝑘𝑚𝑛(𝑥), (24) where ℎ_𝑗^𝑘𝑚𝑛=_𝑉ℎ(𝑥)𝑈_𝑗^𝑘𝑚𝑛(𝑥)𝑑𝑥. (25) The constants ℎ_𝑗^𝑘𝑚𝑛 defined in (25) are called the Fourier coefficients of ℎ(𝑥) relative to the basis functions {𝑈_𝑗^𝑘𝑚𝑛(𝑥)}.

Lemma 7 (see [18]). Let ℎ(𝑥)∈𝐶^∞₀(𝑉). Then _𝑉|ℎ(𝑥)|²𝑑𝑥=^∞∞_𝑘=0^∞_𝑚=0𝑛=0ℎ_𝑗^{𝑘𝑚𝑛2}<∞ (26) holds. Here ℎ_𝑗^𝑘𝑚𝑛 are the Fourier coefficients (25). The equality (26) is known as Parseval’s equality. The proof of this lemma can be found, for example, in [18].

As a result, any vector function Ψ(𝑥)=(𝜓₁(𝑥),𝜓₂(𝑥),𝜓₃(𝑥))^𝑇 with components from 𝐶^∞₀(𝑉) can be presented in the form Ψ(𝑥)=^∞∞_𝑘=0^∞_𝑚=0𝑛=0𝜓₁^𝑘𝑚𝑛𝑈₁^𝑘𝑚𝑛(𝑥)𝐞¹+𝜓₂^𝑘𝑚𝑛𝑈₂^𝑘𝑚𝑛(𝑥)𝐞²+𝜓₃^𝑘𝑚𝑛𝑈₃^𝑘𝑚𝑛(𝑥)𝐞³,∀𝑥∈𝑉, (27) where 𝜓_𝑗^𝑘𝑚𝑛=_𝑉𝜓_𝑗(𝑥)𝑈_𝑗^𝑘𝑚𝑛(𝑥)𝑑𝑥,𝑗=1,2,3. (28) We note that if 𝜓_𝑗(𝑥)∈𝐶^∞₀(𝑉) then ∑^∞_𝑘=0∑^∞_𝑚=0∑^∞_𝑛=0𝜓_𝑗^𝑘𝑚𝑛𝑈_𝑗^𝑘𝑚𝑛(𝑥); 𝑗=1,2,3, is uniformly convergent series over 𝑉.

4. Computation of Frequency-Dependent Magnetic Matrix Green’s Function

4.1. The Regularization of curl_𝑥[𝛿(𝑥−𝑥⁰)𝐞^→𝑠]

Let 𝑠=1, 𝑥=(𝑥₁,𝑥₂,𝑥₃)∈𝐑³, 𝑥⁰=(𝑥⁰₁,𝑥⁰₂,𝑥⁰₃)∈𝐑³, 𝛿(𝑥−𝑥⁰)=𝛿(𝑥₁−𝑥⁰₁)𝛿(𝑥₂−𝑥⁰₂)𝛿(𝑥₃−𝑥⁰₃) be the Dirac delta distribution. The vector distribution −curl_𝑥[𝛿(𝑥−𝑥⁰)𝐞¹]=curl_𝑥0[𝛿(𝑥−𝑥⁰)𝐞¹] can be written as follows: −curl_𝑥𝛿𝑥−𝑥⁰𝐞¹=curl_𝑥0𝛿𝑥−𝑥⁰𝐞¹=0⋅𝐞¹+𝜕𝛿𝑥−𝑥⁰𝜕𝑥⁰₃𝐞²−𝜕𝛿𝑥−𝑥⁰𝜕𝑥⁰₂𝐞³. (29) Using formulas (26) and (28) we find the Fourier series expansion of curl_𝑥0[𝛿(𝑥−𝑥⁰)𝐞¹] in the form of (27), where 𝜓₁^𝑘𝑚𝑛=0,𝜓₂^𝑘𝑚𝑛=𝜕𝑈₂^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₃,𝜓₃^𝑘𝑚𝑛=−𝜕𝑈₃^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₂. (30) (See, e.g., [19–22].)

Lemma 8. Let 𝑥⁰∈𝑉 be fixed, and ℎ(𝑥)∈𝐶^∞₀(𝑉). Then the series ^∞∞_𝑘=0^∞_𝑚=1𝑛=0𝜕𝑈₂^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₃ℎ₂^𝑘𝑚𝑛,_{∞∞𝑘=0∞𝑚=0𝑛=1}𝜕𝑈₃^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₂ℎ₃^𝑘𝑚𝑛 (31)(32) are absolutely convergent. Here ℎ₂^𝑘𝑚𝑛 and ℎ₃^𝑘𝑚𝑛 are defined by (25) for the given function ℎ(𝑥).

The proof of Lemma 8 can be done by standard arguments of calculus (see Appendix A).

Using Lemma 8 the distributions 𝜕𝛿(𝑥−𝑥⁰)/𝜕𝑥⁰₃ and 𝜕𝛿(𝑥−𝑥⁰)/𝜕𝑥⁰₂ can be defined by 𝜕𝛿𝑥−𝑥⁰𝜕𝑥⁰₃(𝑥=𝑥,ℎ)𝜕ℎ⁰𝜕𝑥⁰₃=^∞∞_𝑘=0^∞_𝑚=1𝑛=0𝜕𝑈₂^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₃ℎ₂^𝑘𝑚𝑛,𝜕𝛿𝑥−𝑥⁰𝜕𝑥⁰₂(𝑥=𝑥,ℎ)𝜕ℎ⁰𝜕𝑥⁰₂=^∞∞_𝑘=0^∞_𝑚=0𝑛=1𝜕𝑈₃^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₂ℎ₃^𝑘𝑚𝑛, (33) for any function ℎ(𝑥)∈𝐶^∞₀(𝑉).

Let us consider the infinitely differentiable functions 𝜓₂^1,𝑁(𝑥,𝑥⁰), 𝜓₃^1,𝑁(𝑥,𝑥⁰) defined by 𝜓₂^1,𝑁𝑥,𝑥⁰=^𝑁𝑁_𝑘=0^𝑁_𝑚=1𝑛=0𝜕𝑈₂^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₃𝑈₂^𝑘𝑚𝑛(𝑥),𝜓₃^1,𝑁𝑥,𝑥⁰=−^𝑁𝑁_𝑘=0^𝑁_𝑚=0𝑛=1𝜕𝑈₃^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₂𝑈₃^𝑘𝑚𝑛(𝑥). (34)(35)

Lemma 9. Let 𝑥⁰∈𝑉 be fixed, and let functions 𝜓₂^1,𝑁(𝑥,𝑥⁰) and 𝜓₃^1,𝑁(𝑥,𝑥⁰) be defined by (34) and (35). Then for any ℎ(𝑥)∈𝐶^∞₀(𝑉)lim_{𝑁→∞𝑉}𝜓₂^1,𝑁𝑥,𝑥⁰ℎ(𝑥)𝑑𝑥=𝜕𝛿𝑥−𝑥⁰𝜕𝑥⁰₃,ℎ(𝑥),lim_{𝑁→∞𝑉}𝜓₃^1,𝑁𝑥,𝑥⁰ℎ(𝑥)𝑑𝑥=−𝜕𝛿𝑥−𝑥⁰𝜕𝑥⁰₂,ℎ(𝑥). (36)(37)

The proof of Lemma 7 follows from Lemma 4 and standard arguments of the theory of distributions (see, e.g., [19, 20]).

Remark 10. Equalities (36) and (37) mean that lim_𝑁→∞𝜓₂^1,𝑁𝑥,𝑥⁰=𝜕𝛿𝑥−𝑥⁰𝜕𝑥⁰₃,lim_𝑁→∞𝜓₃^1,𝑁𝑥,𝑥⁰=−𝜕𝛿𝑥−𝑥⁰𝜕𝑥⁰₂ (38) in the space of distributions 𝒟^′(𝑉) (see, e.g., [19, 20]).

Remark 11. The functions 𝜓₂^1,𝑁(𝑥,𝑥⁰) and 𝜓₃^1,𝑁(𝑥,𝑥⁰) can be taken as regularizations (approximations) of the distributions 𝜕𝛿(𝑥−𝑥⁰)/𝜕𝑥⁰₃ and −𝜕𝛿(𝑥−𝑥⁰)/𝜕𝑥⁰₂, respectively. Here 𝑁 is a parameter of the regularization (see, e.g., [19, 20]).

Remark 12. Let Ψ^1,𝑁𝑥,𝑥⁰=0,𝜓₂^1,𝑁,𝜓₃^1,𝑁𝑇=^𝑁𝑁_𝑘=0^𝑁_𝑚=0𝑛=0𝜕𝑈₂^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₃𝑈₂^𝑘𝑚𝑛(𝑥)𝐞²−𝜕𝑈₃^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₂𝑈₃^𝑘𝑚𝑛(𝑥)𝐞³. (39) Using Lemma 9 and Remarks 10 and 11 we find that lim_𝑁→∞Ψ^1,𝑁𝑥,𝑥⁰=curl_𝑥0𝛿𝑥−𝑥⁰𝐞¹ (40) in the space of distributions 𝒟^′(𝑉) and therefore the vector function Ψ^1,𝑁(𝑥,𝑥⁰), with infinitely differentiable components, can be taken as a regularization of the generalized vector function curl_𝑥0[𝛿(𝑥−𝑥⁰)𝐞¹].
Similarly, we find that the vector functions Ψ^2,𝑁(𝑥,𝑥⁰) and Ψ^3,𝑁(𝑥,𝑥⁰), defined by Ψ^2,𝑁𝑥,𝑥⁰=^𝑁𝑁_𝑘=0^𝑁_𝑚=0𝑛=0−𝜕𝑈₁^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₃𝑈₁^𝑘𝑚𝑛(𝑥)𝐞¹+𝜕𝑈₃^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₁𝑈₃^𝑘𝑚𝑛(𝑥)𝐞³,Ψ^3,𝑁𝑥,𝑥⁰=^𝑁𝑁_𝑘=0^𝑁_𝑚=0𝑛=0𝜕𝑈₁^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₂𝑈₁^𝑘𝑚𝑛(𝑥)𝐞¹−𝜕𝑈₂^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₁𝑈₂^𝑘𝑚𝑛(𝑥)𝐞², (41)(42) can be taken as the regularizations of curl_𝑥0[𝛿(𝑥−𝑥⁰)𝐞²] and curl_𝑥0[𝛿(𝑥−𝑥⁰)𝐞³].

4.2. Solving (6) and (7) with the Regularized Inhomogeneous Term

Let us consider (6) and (7) for frequency-dependent magnetic Green’s function. We replace curl_𝑥[𝛿(𝑥−𝑥⁰)𝐞^𝑠] by their regularizations defined by (39), (41), and (42) for 𝑠=1,2,3, respectively. Moreover, (6) and (7) are replaced by the following equations: Δ_𝑥𝐈_𝟑𝐇^𝑠,𝑁+𝜔²𝜖𝜇𝐇^𝑠,𝑁=Ψ^𝑠,𝑁𝑥,𝑥⁰,curl_𝑥𝐇^𝑠,𝑁×𝐧^→_Γ=0,𝐇^𝑠,𝑁⋅𝐧^→_Γ=0. (43)(44) We will find solutions of the problem (43) and (44) in the form 𝐇^𝑠,𝑁𝑥,𝜔,𝑥⁰=^𝑁𝑁_𝑘=0^𝑁_𝑚=0𝑛=0𝑌₁^{𝑠,𝑘𝑚𝑛}𝜔,𝑥⁰𝑈₁^𝑘𝑚𝑛(𝑥)𝐞¹+𝑌₂^{𝑠,𝑘𝑚𝑛}𝜔,𝑥⁰𝑈₂^𝑘𝑚𝑛(𝑥)𝐞²+𝑌₃^{𝑠,𝑘𝑚𝑛}𝜔,𝑥⁰𝑈₃^𝑘𝑚𝑛(𝑥)𝐞³, (45) where 𝑥∈𝑉; 𝑈_𝑗^𝑘𝑚𝑛(𝑥) are eigenfunctions of the eigenvalue-eigenfunction problem (18) and (19); 𝑌_𝑗^{𝑠,𝑘𝑚𝑛}(𝜔,𝑥⁰) and 𝑗=1,2,3, are the unknown Fourier coefficients.

Substituting (45) into (43) and using orthogonality of 𝑈_𝑗^𝑘𝑚𝑛(𝑥), we find 𝜔²𝜖𝜇−𝜆_𝑘𝑚𝑛𝑌_𝑗^{𝑠,𝑘𝑚𝑛}𝜔,𝑥⁰=𝜓_𝑗^{𝑠,𝑘𝑚𝑛}𝑥⁰, (46) for 𝑗=1,2,3; 𝑠=1,2,3; 𝑘=0,1,…, 𝑚=0,1,…, and 𝑛=0,1,….

Let 𝜔²≠𝜆_𝑘𝑚𝑛/𝜖𝜇; then 𝑌_𝑗^{𝑠,𝑘𝑚𝑛}(𝜔,𝑥⁰) is uniquely determined by 𝑌_𝑗^{𝑠,𝑘𝑚𝑛}𝜔,𝑥⁰=𝜓_𝑗^{𝑠,𝑘𝑚𝑛}𝑥⁰𝜔²𝜖𝜇−𝜆_𝑘𝑚𝑛, (47) where 𝜓_𝑗^{𝑠,𝑘𝑚𝑛}𝑥⁰=_𝑉Ψ^𝑠_𝑗𝑥,𝑥⁰𝑈_𝑗^𝑘𝑚𝑛(𝑥)𝑑𝑥, (48)𝑘,𝑚,𝑛=0,1,…. As a result we obtain the following explicit formula: 𝐇^1,𝑁𝑥,𝜔,𝑥⁰=^𝑁𝑁_𝑘=0^𝑁_𝑚=0𝑛=01𝜔²𝜖𝜇−𝜆_𝑘𝑚𝑛×𝜕𝑈₂^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₃𝑈₂^𝑘𝑚𝑛(𝑥)𝐞²−𝜕𝑈₃^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₂𝑈₃^𝑘𝑚𝑛(𝑥)𝐞³, (49) for all 𝑥∈𝑉, 𝑥⁰∈𝑉, and all parameter 𝜔 satisfying 𝜔²𝜖𝜇≠𝜆_𝑘𝑚𝑛   (𝑘=0,1,…; 𝑚=0,1,…; 𝑛=0,1,…).

Similarly, the solutions of (43) and (44) for 𝑠=2 and 𝑠=3 are obtained in the form 𝐇^2,𝑁𝑥,𝜔,𝑥⁰=^𝑁𝑁_𝑘=0^𝑁_𝑚=0𝑛=01𝜔²𝜖𝜇−𝜆_𝑘𝑚𝑛×𝜕𝑈₁^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₃𝑈₁^𝑘𝑚𝑛(𝑥)𝐞¹−𝜕𝑈₃^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₁𝑈₃^𝑘𝑚𝑛(𝑥)𝐞³,𝐇^3,𝑁𝑥,𝜔,𝑥⁰=^𝑁𝑁_𝑘=0^𝑁_𝑚=0𝑛=01𝜔²𝜖𝜇−𝜆_𝑘𝑚𝑛×𝜕𝑈₁^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₂𝑈₁^𝑘𝑚𝑛(𝑥)𝐞¹−𝜕𝑈₂^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₁𝑈₂^𝑘𝑚𝑛(𝑥)𝐞². (50)(51) Using the explicit formulas for 𝐇^𝑠,𝑁(𝑥,𝜔,𝑥⁰) and 𝑈_𝑗^𝑘𝑚𝑛(𝑥) we can check directly that 𝐇^𝑠,𝑁(𝑥,𝜔,𝑥⁰) satisfy the boundary condition (44).

Remark 13. By direct computation we can show that the vector functions 𝐇^𝑠,𝑁(𝑥,𝜔,𝑥⁰), defined by, (49)–(51), satisfy div_𝑥𝐇^𝑠,𝑁𝑥,𝜔,𝑥⁰=0,𝑠=1,2,3. (52)

4.3. Frequency-Dependent Magnetic Green’s Function and Its Regularization

In this section we show that lim_𝑁→∞𝐻_𝑗^𝑠,𝑁(𝑥,𝜔,𝑥⁰); (𝑠=1,2,3;𝑗=1,2,3) exists in the space of distributions 𝒟^′(𝑉), and the distributions 𝐻^𝑠_𝑗=lim_𝑁→∞𝐻_𝑗^𝑠,𝑁 are the entries of matrix-valued Green’s function 𝐆^𝐻(𝑥,𝜔,𝑥⁰).

Let us consider for example, 𝑠=1, 𝑗=2 and an arbitrary function ℎ(𝑥)∈𝐶^∞₀(𝑉) of the form (24) and (25). The integrable over 𝑉 function 𝐻_𝑗^𝑠,𝑁(𝑥,𝑡,𝑥⁰) defines a regular functional (distribution) by 𝐻₂^1,𝑁𝑥,𝜔,𝑥⁰(𝑥=,ℎ)_{𝑁𝑁𝑘=0𝑁𝑚=1𝑛=0}1𝜔²𝜖𝜇−𝜆_𝑘𝑚𝑛𝜕𝑈₂^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₃ℎ₂^𝑘𝑚𝑛. (53) Using Lemma 8 we find that ∑^∞_𝑘=0∑^∞_𝑚=1∑^∞_𝑛=0(1/(𝜔²𝜖𝜇−𝜆_𝑘𝑚𝑛))(𝜕𝑈₂^𝑘𝑚𝑛(𝑥⁰)/𝜕𝑥⁰₃)ℎ₂^𝑘𝑚𝑛 is absolutely convergent for any ℎ(𝑥)∈𝐶^∞₀(𝑉). This means that lim_𝑁→∞⟨𝐻₂^1,𝑁(𝑥,𝜔,𝑥⁰),ℎ(𝑥)⟩ exists for any ℎ(𝑥)∈𝐶^∞₀(𝑉) and defines the linear continuous functional 𝐻₂^1,𝑁(𝑥,𝜔,𝑥⁰) over 𝐶^∞₀(𝑉).

Similarly, we can find 𝐻^𝑠_𝑗(𝑥,𝜔,𝑥⁰) as lim_𝑁→∞𝐻_𝑗^𝑠,𝑁(𝑥,𝜔,𝑥⁰) in the complete space 𝒟^′(𝑉) (see, e.g., [19, 20]) for 𝑠=1,2,3 and 𝑗=1,2,3.

Moreover, since 𝐇^𝑠,𝑁(𝑥,𝜔,𝑥⁰)=(𝐻₁^𝑠,𝑁(𝑥,𝜔,𝑥⁰),𝐻₂^𝑠,𝑁(𝑥,𝜔,𝑥⁰),𝐻₃^𝑠,𝑁(𝑥,𝑡𝜔,𝑥⁰))^𝑇 satisfies equalities (43) and (44) then 𝐇^𝑠(𝑥,𝜔,𝑥⁰)=(𝐻^𝑠₁(𝑥,𝜔,𝑥⁰),𝐻^𝑠₂(𝑥,𝜔,𝑥⁰),𝐻^𝑠₃(𝑥,𝜔,𝑥⁰))^𝑇 satisfies the following equalities: 𝜔²𝜖𝜇𝐇^𝑠−curl_𝑥curl_𝑥𝐇^𝑠=−curl_𝑥𝛿𝑥−𝑥⁰𝐞^𝑠,curl_𝑥𝐇^𝑠×𝐧^→_Γ=0,𝐇^𝑠⋅𝐧^→_Γ=0, (54) that is, equalities (6) and (7). This means that 𝐇^𝑠(𝑥,𝜔,𝑥⁰) is the 𝑠-column of frequency-dependent magnetic matrix Green’s function.

5. Approximate Solution of Electric Green’s Function

Let columns of magnetic Green’s function be derived approximately; that is, let solutions 𝐇^𝑠,𝑁, 𝑠=1,2,3 of the problem (43) and (44) be computed by formulae (49), (50), and (51). For the approximate derivation of 𝐄^𝑠(𝑥,𝜔,𝑥⁰), satisfying (8), we replace curl_𝑥𝐇^𝑠(𝑥,𝜔,𝑥⁰) by curl_𝑥𝐇^𝑠,𝑁(𝑥,𝜔,𝑥⁰) and 𝛿(𝑥−𝑥⁰)𝐞^𝑠 by its regularization ∑^𝑁_𝑘=0∑^𝑁_𝑚=0∑^𝑁_𝑛=0𝑈_𝑠^𝑘𝑚𝑛(𝑥⁰)𝑈_𝑠^𝑘𝑚𝑛(𝑥)𝐞^𝑠 in (8) for 𝑠=1,2,3, respectively, and find an explicit formula for 𝐄^𝑠,𝑁(𝑥,𝜔,𝑥⁰) as follows: 𝐄^𝑠,𝑁𝑥,𝜔,𝑥⁰=𝑖𝜔𝜖curl_𝑥𝐇^𝑠,𝑁𝑥,𝜔,𝑥⁰−^𝑁𝑁_𝑘=0^𝑁_𝑚=0𝑛=0𝑈_𝑠^𝑘𝑚𝑛𝑥⁰𝑈_𝑠^𝑘𝑚𝑛(𝑥)𝐞^𝑠, (55) for all 𝑥∈𝑉, 𝑥⁰∈𝑉, and all parameter 𝜔 satisfying 𝜔²𝜖𝜇≠𝜆_𝑘𝑚𝑛.

5.1. The Regularization of 𝐄^𝑠(𝑥,𝜔,𝑥⁰)

Let us consider the differentiable function 𝐸₂^2,𝑁(𝑥,𝜔,𝑥⁰), given by (55). This function defines a linear continuous regular functional over 𝐶^∞₀(𝑉) by the following relation (see, e.g., [19, 20]): 𝐸₂^2,𝑁(𝑥𝑖,ℎ)=−𝜔𝜖^𝑁𝑁_𝑘=0^𝑁_𝑚=0𝑛=01𝜔²𝜖𝜇−𝜆_𝑘𝑚𝑛×𝜕𝑈₁^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₃×_𝑉𝑈₁^𝑘𝑚𝑛(𝑥)𝜕ℎ(𝑥)𝜕𝑥₃+𝑑𝑥𝜕𝑈₃^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₁×_𝑉𝑈₃^𝑘𝑚𝑛(𝑥)𝜕ℎ(𝑥)𝜕𝑥₁+𝑖𝑑𝑥𝜔𝜖^𝑁𝑁_𝑘=0^𝑁_𝑚=0𝑛=0𝑈₂^𝑘𝑚𝑛𝑥⁰ℎ₂^𝑘𝑚𝑛, (56) for an arbitrary ℎ(𝑥)∈𝐶^∞₀(𝑉). The series in the right-hand side of the last equality is absolutely convergent for any ℎ(𝑥)∈𝐶^∞₀(𝑉) (see, Appendix B). This means that lim_𝑁→∞⟨𝐸₂^2,𝑁(𝑥,𝜔,𝑥⁰),ℎ(𝑥)⟩ exists for any ℎ(𝑥)∈𝐶^∞₀(𝑉) and defines a linear continuous regular functional over 𝐶^∞₀(𝑉). We denote this functional as 𝐸²₂(𝑥,𝜔,𝑥⁰); that is, 𝐸²₂(𝑥,𝜔,𝑥⁰)=lim_𝑁→∞𝐸₂^2,𝑁(𝑥,𝜔,𝑥⁰) in the complete space 𝒟^′(𝑉). Similarly, we can find 𝐸^𝑠_𝑗(𝑥,𝜔,𝑥⁰) as lim_𝑁→∞𝐸_𝑗^𝑠,𝑁(𝑥,𝜔,𝑥⁰) in the space 𝒟^′(𝑉) for 𝑠=1,2,3 and 𝑗=1,2,3.

The 𝑠-column of electric Green’s function is 𝐄^𝑠𝑥,𝜔,𝑥⁰=𝑖𝜔𝜖curl_𝑥𝐇^𝑠𝑥,𝜔,𝑥⁰−𝛿𝑥−𝑥⁰𝐞^𝑠. (57) Let now 𝐄^𝑠(𝑥,𝜔,𝑥⁰) = (𝐸^𝑠₁(𝑥,𝜔,𝑥⁰), 𝐸^𝑠₂(𝑥,𝜔,𝑥⁰), 𝐸^𝑠₃(𝑥,𝜔,𝑥⁰))^𝑇 be the generalized vector function, whose components are distributions defined as 𝐸^𝑠_𝑗=lim_𝑁→∞𝐸_𝑗^𝑠,𝑁, 𝑗=1,2,3. Since 𝐄^𝑠,𝑁(𝑥,𝜔,𝑥⁰)=(𝐸₁^𝑠,𝑁(𝑥,𝜔,𝑥⁰),𝐸₂^𝑠,𝑁(𝑥,𝜔,𝑥⁰),𝐸₃^𝑠,𝑁(𝑥,𝜔,𝑥⁰))^𝑇 satisfies (55), then the generalized vector function 𝐄^𝑠(𝑥,𝜔,𝑥⁰) satisfies (4).

Remark 14. We note that the generalized vector functions 𝐇^𝑠(𝑥,𝜔,𝑥⁰), 𝐄^𝑠(𝑥,𝜔,𝑥⁰), defined in Sections 4.3 and 5.1, satisfy (6), (7), and (8). It follows that 𝐇^𝑠(𝑥,𝜔,𝑥⁰) and 𝐄^𝑠(𝑥,𝜔,𝑥⁰) satisfy (3)–(5), and therefore, 𝐇^𝑠(𝑥,𝜔,𝑥⁰) and 𝐄^𝑠(𝑥,𝜔,𝑥⁰) are 𝑠-columns of frequency-dependent electric and magnetic Green’s functions 𝐆^𝐸(𝑥,𝜔,𝑥⁰) and 𝐆^𝐻(𝑥,𝜔,𝑥⁰), respectively.

6. Computational Experiments

In this section, we present computational experiments to show applicability of the approximate computation of frequency-dependent magnetic and electric matrix Green’s functions in the rectangular parallelepiped. Here we consider the parallelepiped 𝑉=𝑥∈𝐑³∶0<𝑥₁<3,0<𝑥₂<4,0<𝑥₃<2.5 (58) which contains a dielectric material characterized by the electric permittivity 𝜖 and the magnetic permeability 𝜇. We choose 𝜖 and 𝜇 as 𝑐=1/√𝜖𝜇   (𝑐≈3⋅10⁸ m/s is the speed of light) and 𝜔 is equal to the value of 𝑐. The boundary of 𝑉, denoted as Γ, is a perfect conducting boundary on which the tangential components of the electric field component and normal component of the magnetic field component vanish. There is an electric current 𝐉(𝑥,𝜔) inside 𝑉 which produces the electric and magnetic fields. In our experiments we take 𝐉^𝑠(𝑥,𝜔), 𝑠=1,2,3, as 𝐉^𝑠(𝑥,𝜔)=𝛿𝑥−𝑥⁰𝐞^𝑠, (59)

where 𝐞^𝑠 are the basis vectors from 𝐑³; 𝑥⁰=(𝑥⁰₁,𝑥⁰₂,𝑥⁰₃) is the 3D parameter from 𝑉 and chosen as 𝑥⁰₁=2, 𝑥⁰₂=3, 𝑥⁰₃=1.5; 𝛿(𝑥−𝑥⁰)=𝛿(𝑥₁−𝑥⁰₁)𝛿(𝑥₂−𝑥⁰₂)𝛿(𝑥₃−𝑥⁰₃), 𝛿(𝑥_𝑗−𝑥⁰_𝑗) is the Dirac delta distribution at 𝑥_𝑗=𝑥⁰_𝑗, 𝑗=1,2,3.

Two types of computational experiments will be presented. The first one is related to convergence and error analysis of our method; the second one is the computation of the elements of matrix magnetic and electric Green’s functions in the rectangular parallelepiped 𝑉 and the visualization of the computed magnetic and electric fields components by MATLAB graphic tools.

6.1. Convergence Analysis

Green’s function of the magnetic field is 3×3 matrix whose elements are distributions. Generally speaking, distributions do not have values at fixed points [19, 20]. The elements of regularized (approximate) matrix Green’s function, which we have computed by our method, are classical differentiable and integrable functions over 𝑉 having the values at fixed points. These classical functions define the regular distributions by the standard rule (see, [19, 20]). The comparison of the distributions can be made by the comparison of their integral characteristics. Some integral characteristics of the elements of matrix Green’s function of the magnetic field can be found explicitly. For example, let 𝐻¹₂(𝑥,𝜔,𝑥⁰) be a solution of (6) and let (7) for 𝑠=1 and 𝑗=2 and let ℎ¹₂(𝑥₃,𝜔,𝑥⁰) be defined by ℎ¹₂𝑥₃,𝜔,𝑥⁰=^𝑏1₀^𝑏2₀𝐻¹₂𝑥,𝜔,𝑥⁰cos𝜋𝑥₁𝑏₁sin𝜋𝑥₂𝑏₂𝑑𝑥₂𝑑𝑥₁. (60) We can show that ℎ¹₂=(𝑥₃,𝜔,𝑥⁰) satisfies 𝜕²ℎ¹₂𝜕𝑥²₃+𝑘²ℎ¹₂𝜕=−𝜕𝑥₃𝛿𝑥₃−𝑥⁰₃cos𝜋𝑥⁰₁𝑏₁sin𝜋𝑥⁰₂𝑏₂,0<𝑥₃<𝑏₃,𝑑ℎ¹₂𝑑𝑥_3𝑥3=0=0,𝑑ℎ¹₂𝑑𝑥_3𝑥3=𝑏3=0, (61)(62) where 𝑘²=𝜔²𝜖𝜇−(𝜋/𝑏₁)²−(𝜋/𝑏₂)². A solution of (61) and (62) can be written in the form ℎ¹₂𝑥₃,𝜔,𝑥⁰=sin𝑘𝑏₃−𝑥⁰₃sin𝑘𝑏₃cos𝑘𝑥₃×cos𝜋𝑥⁰₁𝑏₁sin𝜋𝑥⁰₂𝑏₂,0<𝑥₃<𝑥⁰₃,−sin𝑘𝑥⁰₃sin𝑘𝑏₃cos𝑘𝑏₃−𝑥₃×cos𝜋𝑥⁰₁𝑏₁sin𝜋𝑥⁰₂𝑏₂,𝑥⁰₃<𝑥₃<𝑏₃. (63) We can compare the values of ℎ¹₂(𝑥₃,𝜔,𝑥⁰), as the values of the integral characteristics of 𝐻¹₂, with the values of the integral characteristic ℎ₂^1,𝑁𝑥₃,𝜔,𝑥⁰=^𝑏1₀^𝑏2₀𝐻₂^1,𝑁𝑥,𝜔,𝑥⁰cos𝜋𝑥₁𝑏₁sin𝜋𝑥₂𝑏₂𝑑𝑥₂𝑑𝑥₁=^𝑁_𝑛=1−1𝑘²−𝑛𝜋/𝑏₃²2𝑛𝜋𝑏₃²×cos𝜋𝑥⁰₁𝑏₁sin𝜋𝑥⁰₂𝑏₂sin𝑛𝜋𝑥⁰₃𝑏₃cos𝑛𝜋𝑥₃𝑏₃ (64) of 𝐻₂^1,𝑁. Here 𝐻₂^1,𝑁 is the element of approximate matrix Green’s function computed by our method and 𝑁 is the parameter of the approximation corresponding to the number of the terms in each series of the right hand side of (64).

In Figures 1, 2, 3, and 4 and Table 1 the values of ℎ₂^1,𝑁 are compared with the values of ℎ¹₂ for fixed 𝑥⁰₁=2, 𝑥⁰₂=3, 𝑥⁰₃=1.5, and 𝜖 and 𝜇 such that 𝑐=1/√𝜖𝜇 and 𝜔=𝑐 and different numbers 𝑁. The maximum of relative error 𝑟_𝑁(𝑥₃)=|ℎ¹₂−ℎ₂^1,𝑁|/|ℎ¹₂| at the point 𝑥₃=1.3 of the interval [0,2.5] for 𝑁=15,50,100,150 is 𝑟₁₅=0.1806, 𝑟₅₀=0.0625, 𝑟₁₀₀=0.0316, and 𝑟₁₅₀=0.0211, respectively. Making the error analysis we have concluded that the biggest 𝑁 corresponds to the smallest error between values of ℎ₂^1,𝑁 and ℎ¹₂. We have taken 𝑁=100 as a reasonable value between accuracy and the time for our computation in MATLAB.

Table 1 

The accuracy of the calculation of the magnetic field for 𝑥⁰₁=2; 𝑥⁰₂=3; 𝑥⁰₃=1.5; 𝑟_𝑁=|ℎ¹₂−ℎ₂^1,𝑁|/|ℎ¹₂|; 𝑁=15, 50, 100, and 150.

Figure 1 

The graphs of ℎ12 and ℎ21,𝑁 for 𝑁=15, 𝑥01=2, 𝑥02=3, and 𝑥03=1.5.

Figure 2 

The graphs of ℎ12 and ℎ21,𝑁 for 𝑁=50, 𝑥01=2, 𝑥02=3, and 𝑥03=1.5.

Figure 3 

The graphs of ℎ12 and ℎ21,𝑁 for 𝑁=100, 𝑥01=2, 𝑥02=3, and 𝑥03=1.5.

Figure 4 

The graphs of ℎ12 and ℎ21,𝑁 for 𝑁=150, 𝑥01=2, 𝑥02=3, and 𝑥03=1.5.

6.2. Approximate Computation and Visualization of the Elements of Matrix Magnetic and Electric Green’s Functions

The method described in Sections 4 and 5 has been applied for the approximate computation of frequency-dependent magnetic and electric Green’s matrices in the rectangular parallelepiped 𝑉={𝑥∈𝐑³∶0<𝑥₁<𝑏₁,0<𝑥₂<𝑏₂,0<𝑥₃<𝑏₃}, respectively. We note that each column of magnetic and electric frequency-dependent Green’s matrices is the magnetic and electric fields arising from the pulse dipole polarized in the direction of the basis vector of 𝐑³. For computation we take 𝑏₁=3, 𝑏₂=4, and 𝑏₃=2.5 and 𝜀 and 𝜇 such that 𝑐=1/√𝜀𝜇, 𝜔=𝑐 (𝑐 is the value of the light velocity), 𝑥⁰₁=2, 𝑥⁰₂=3, and 𝑥⁰₃=1.5. Computation of 𝐻_𝑗^𝑠,𝑁 and 𝐸_𝑗^𝑠,𝑁 has been implemented in MATLAB and their graphs are presented. Figures 5(a), 5(b), 6(a), and 6(b) present the second and third components of the first column of the regularized magnetic matrix, respectively, 𝐻₂^1,100(𝑥,𝜔,𝑥⁰) and 𝐻₃^1,100(𝑥,𝜔,𝑥⁰). Figures 7(a)-7(b), 8(a)-8(b), and 9(a)-9(b) present each component of the first column of the regularized electric matrix, respectively, 𝐸₁^1,100(𝑥,𝜔,𝑥⁰), 𝐸₂^1,100(𝑥,𝜔,𝑥⁰), and 𝐸₃^1,100(𝑥,𝜔,𝑥⁰). The views of the surfaces 𝑧=𝐻₂^1,𝑁, 𝑧=𝐻₃^1,𝑁, 𝑧=𝐸₁^1,𝑁, 𝑧=𝐸₂^1,𝑁, and 𝑧=𝐸₃^1,𝑁 (𝑁=100) from the top of 𝑧-axis are presented in Figures 5(a), 6(a), 7(a), 8(a), and 9(a), where 𝑥₁-axis is vertical and 𝑥₂-axis is a horizontal one. Figures 5(b), 6(b), 7(b), 8(b), and 9(b) are 3D images of the surfaces 𝐻₂^1,100(𝑥,𝜔,𝑥⁰), 𝐻₃^1,100(𝑥,𝜔,𝑥⁰), 𝐸₁^1,100(𝑥,𝜔,𝑥⁰), 𝐸₂^1,100(𝑥,𝜔,𝑥⁰), and 𝐸₃^1,100(𝑥,𝜔,𝑥⁰) for 𝑥₃=1.5, respectively. Here, 𝑥₁ and 𝑥₂ are horizontal axes and 𝑧-axis is the magnitude.

(a) 𝐻21,𝑁(𝑥1,𝑥2,1.5,𝜔,𝑥0), 𝑁=100

(b) 𝐻21,𝑁(𝑥1,𝑥2,1.5,𝜔,𝑥0), 𝑁=100

Figure 5 

The second component of the first column of regularized magnetic Green’s function.

(a) 𝐻31,𝑁(𝑥1,𝑥2,1.5,𝜔,𝑥0), 𝑁=100

(b) 𝐻31,𝑁(𝑥1,𝑥2,1.5,𝜔,𝑥0), 𝑁=100

Figure 6 

The third component of the first column of regularized magnetic Green’s function.

(a) 𝐸11,𝑁(𝑥1,𝑥2,1.5,𝜔,𝑥0), 𝑁=100

(b) 𝐸11,𝑁(𝑥1,𝑥2,1.5,𝜔,𝑥0), 𝑁=100

Figure 7 

The first component of the first column of regularized electric Green’s function.

(a) 𝐸21,𝑁(𝑥1,𝑥2,1.5,𝜔,𝑥0), 𝑁=100

(b) 𝐸21,𝑁(𝑥1,𝑥2,1.5,𝜔,𝑥0), 𝑁=100

Figure 8 

The second component of the first column of regularized electric Green’s function.

(a) 𝐸31,𝑁(𝑥1,𝑥2,1.5,𝜔,𝑥0), 𝑁=100

(b) 𝐸31,𝑁(𝑥1,𝑥2,1.5,𝜔,𝑥0), 𝑁=100

Figure 9 

The third component of the first column of regularized electric Green’s function.

7. Conclusion

The new method for the approximate computation of frequency-dependent magnetic and electric Green’s matrices in the rectangular parallelepiped with a perfect conducting boundary has been developed. This method has the following steps.

The first step is determination of the eigenvalues and corresponding eigenfunctions of the vector Laplace differential operator of the form Δ_𝑥𝐈₃ (where Δ is the Laplace operator and 𝐈₃ is the 3×3 identity matrix) in a parallelepiped subject to some boundary conditions. The components of these vector eigenfunctions form an orthogonal set which is complete in the space of the square integrable functions over the considered parallelepiped.

The second step is an approximation (regularization) of the generalized vector functions (distributions) curl_𝑥[𝛿(𝑥−𝑥⁰)𝐞^𝑠]; 𝑠=1,2,3, appear in the differential equations for the columns of magnetic Green’s matrix as the free terms. Here 𝐞^𝑠 are the basis vectors of 𝐑³ (𝐞¹=(1,0,0), 𝐞²=(0,1,0), 𝐞³=(0,0,1)); 𝑥⁰=(𝑥⁰₁,𝑥⁰₂,𝑥⁰₃) is a 3D parameter from the considered parallelepiped, 𝛿(𝑥−𝑥⁰)=𝛿(𝑥₁−𝑥⁰₁)𝛿(𝑥₂−𝑥⁰₂)𝛿(𝑥₃−𝑥⁰₃), and 𝛿(𝑥_𝑗−𝑥⁰_𝑗) is the Dirac delta function considered at 𝑥_𝑗=𝑥⁰_𝑗 and 𝑗=1,2,3. We find the Fourier coefficients of these free terms relative to the obtained set of eigenfunctions and approximate them by the Fourier series with a finite number of terms. We note that the Dirac delta function is very often used for modeling the point source in physics and engineering (see, [19, 20, 23]). The Dirac delta function does not have point-wise values and can not be drawn. For this reason, a regularization (approximation) of the Dirac delta function in the form of the classical function which has point-wise values is usually used for drawing and computation. For the Dirac delta function, as a singular generalized function, there are different types of approximation (see, e.g., [19, 20, 23]). Unfortunately, the known approximations of the Dirac delta function are not convenient for our approach. For this reason we have suggested a new approximation which is very useful for the computation of considered Green’s function.

The third step is an approximate computation of elements of Green’s matrix in the form of the Fourier series with a finite number of terms relative to the considered set of eigenfunctions. All terms of this series contain known functions which are given explicitly. The simple implementation of our method for computing magnetic Green’s function in a rectangular parallelepiped is based on the obtained presentation and does not contain any type of discretization.

Using explicit formulas for the approximate elements of magnetic Green’s matrix we derive explicitly the approximate elements of electric Green’s matrix. Finally, we have obtained the columns of magnetic and electric matrix Green’s function in the space of generalized functions (distributions) and their regularizations (approximations). The method has been implemented in MATLAB. The computational experiments have confirmed its robustness. We can note also that application, for example, finite elements or finite difference approaches can be done if we approximate the generalized vector functions curl_𝑥[𝛿(𝑥−𝑥⁰)𝐞^𝑠] by vector functions with components which are classical functions. We note here that numerical computations of the values of magnetic Green’s function by equations, where curl_𝑥[𝛿(𝑥−𝑥⁰)𝐞^𝑠] are replaced by their approximation, will contain the error depending on the grid of discretization. For this reason our method of the computation of the magnetic Green function has the accuracy which is significantly better than other methods based on any type of discretization.

Appendices

A. Proof of Lemma 8

Let 𝐷^𝛼=𝜕^{𝛼1+𝛼2+𝛼3}/𝜕^𝛼1_𝑥1𝜕^𝛼2_𝑥2𝜕^𝛼3_𝑥3 be the differential operator, where 𝛼=(𝛼₁,𝛼₂,𝛼₃) and 𝛼₁,𝛼₂,𝛼₃ are nonnegative integers. We have _𝑉𝐷^𝛼ℎ(𝑥)²𝑑𝑥<∞, (A.1) for an arbitrary ℎ(𝑥)∈𝐶^∞₀(𝑉).

Using Parseval’s equality (26) we find _𝑉𝐷^𝛼ℎ(𝑥)²𝑑𝑥=^∞∞_𝑘=0^∞_𝑚=0𝑛=0ℎ^𝑘𝑚𝑛_𝑗,𝛼²<∞. (A.2) Here ℎ^𝑘𝑚𝑛_𝑗,𝛼 denotes the Fourier coefficients of the function 𝐷^𝛼ℎ(𝑥) relative to basis functions {𝑈_𝑗^𝑘𝑚𝑛(𝑥)}, 𝑘,𝑚,𝑛=0,1,…; 𝑗=1,2,3.

For even positive integers 𝛼₁,𝛼₂,𝛼₃ we have ℎ^𝑘𝑚𝑛_𝑗,𝛼=(−1)^|𝛼|/2𝑘𝜋𝑏₁^𝛼1𝑚𝜋𝑏₂^𝛼2𝑛𝜋𝑏₃^𝛼3ℎ_𝑗^𝑘𝑚𝑛, (A.3) where |𝛼|=𝛼₁+𝛼₂+𝛼₃.

Using (21) and (A.3) we have 𝜕𝑈₂^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₃ℎ₂^𝑘𝑚𝑛(𝑥)≤𝑛𝜋𝑏₃𝐵_𝑘𝑚𝑛ℎ₂^𝑘𝑚𝑛≤2𝑐₁1𝑘^𝛼1𝑚^𝛼2𝑛^𝛼3−1ℎ^𝑘𝑚𝑛_2,𝛼≤𝑐₁1𝑘^2𝛼1𝑚^2𝛼2𝑛^2𝛼3−2+ℎ^𝑘𝑚𝑛_2,𝛼², (A.4) where 𝑐₁=𝜋^1−|𝛼|2/(𝑏₁𝑏₂𝑏₃)𝑏^𝛼1₁𝑏^𝛼2₂𝑏^𝛼3₃⁻¹.

The convergence of ∑^∞_𝑘=1∑^∞_𝑚=1∑^∞_𝑛=1|ℎ^𝑘𝑚𝑛_2,𝛼|² follows from (A.2). The series ^∞∞_𝑘=1^∞_𝑚=1𝑛=11𝑘^2𝛼1𝑚^2𝛼2𝑛^2𝛼3−2 (A.5) converges for 2𝛼₁>1, 2𝛼₂>1, and 2(𝛼₃−1)>1. Therefore, taking even positive integers 𝛼₁,𝛼₂,𝛼₃ and using equality (A.4), we find that the series (31) is convergent for any fixed 𝑥⁰∈𝑉 and ℎ(𝑥)∈𝐶^∞₀(𝑉). The convergence of (32) can be established similarly.

Lemma 8 is proved.

B. Convergence of Series (56)

Let us consider the convergence of the following series for any fixed 𝑥⁰∈𝑉 and any function ℎ(𝑥)∈𝐶^∞₀(𝑉): −𝑖𝜔𝜖^𝑁𝑁_𝑘=0^𝑁_𝑚=0𝑛=01𝜔²𝜖𝜇−𝜆_𝑘𝑚𝑛𝜕𝑈₁^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₃×_𝑉𝑈₁^𝑘𝑚𝑛(𝑥)𝜕ℎ(𝑥)𝜕𝑥₃𝑑𝑥,−𝑖𝜔𝜖^𝑁𝑁_𝑘=0^𝑁_𝑚=0𝑛=01𝜔²𝜖𝜇−𝜆_𝑘𝑚𝑛𝜕𝑈₃^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₁×_𝑉𝑉𝑈₃^𝑘𝑚𝑛(𝑥)𝜕ℎ(𝑥)𝜕𝑥₁𝑑𝑥,𝑖𝜔𝜖^𝑁𝑁_𝑘=0^𝑁_𝑚=0𝑛=01𝜔²𝜖𝜇−𝜆_𝑘𝑚𝑛𝑈₂^𝑘𝑚𝑛𝑥⁰ℎ₂^𝑘𝑚𝑛. (B.1)(B.2)(B.3) Let ℎ(𝑥) be an arbitrary function from 𝐶^∞₀(𝑉) and 𝜙(𝑥)=𝜕ℎ(𝑥)/𝜕𝑥₃; 𝛼=(𝛼₁,𝛼₂,𝛼₃), where 𝛼_𝑗, 𝑗=1,2,3 are natural numbers, such that 2𝛼₁>1, 2𝛼₂>1, and 2(𝛼₃−1)>1. Let 𝜙₁^𝑘𝑚𝑛and𝜙^𝑘𝑚𝑛_1,𝛼 be the Fourier coefficients of 𝜙(𝑥) and 𝐷^𝛼𝜙(𝑥) relative to the basis {𝑈₁^𝑘𝑚𝑛(𝑥)}, 𝑘=1,2,…; 𝑚,𝑛=0,1,2,3,…. Then we have 1𝜔𝜖𝜔²𝜖𝜇−𝜆_𝑘𝑚𝑛𝜕𝑈₁^𝑘𝑚𝑛𝑥⁰𝜕𝑥⁰₃×_𝑉𝑈₁^𝑘𝑚𝑛(𝑥)𝜙(𝑥)≤𝐴𝑑𝑥_𝑘𝑚𝑛𝜔𝜖𝜆_𝑘𝑚𝑛−𝜔²𝜖𝜇𝑛𝜋𝑏₃𝜙₁^𝑘𝑚𝑛≤𝐴_𝑘𝑚𝑛𝜔𝜖𝜋𝑛𝑏₃𝜙₁^𝑘𝑚𝑛≤𝑐₂1𝑘^2𝛼1𝑚^2𝛼2𝑛^2(𝛼3−1)+𝜙^𝑘𝑚𝑛_1,𝛼², (B.4) where 𝑐₂=(𝜋^1−|𝛼|/𝜔𝜖)2/(𝑏₁𝑏₂𝑏₃)𝑏^𝛼1₁𝑏^𝛼2₂𝑏^(𝛼3₃⁻¹⁾. Using the technique similar to Appendix A we find that ∑^∞_𝑘=1∑^∞_𝑚=1∑^∞_𝑛=1𝑐₂[1/𝑘^2𝛼1𝑚^2𝛼2𝑛^2(𝛼3−1)+|𝜙^𝑘𝑚𝑛_1,𝛼|²] converges and, as a result of it, series (B.1) converges absolutely for any fixed 𝑥⁰∈𝑉 and any function ℎ(𝑥)∈𝐶^∞₀(𝑉). The convergence of series (B.2) and (B.3) for any fixed 𝑥⁰∈𝑉 and an arbitrary function ℎ(𝑥)∈𝐶^∞₀(𝑉) can be established in a similar way. Convergence of series (B.1), (B.2), and (B.3) implies the absolute convergence of the series in (56) for any fixed 𝑥⁰∈𝑉 and arbitrary ℎ(𝑥)∈𝐶^∞₀(𝑉).

Conflict of Interests

The authors declare that there is no conflict of interests regarding the publication of this paper.

Acknowledgment

This work was supported by the Scientific and Technological Research Council of Turkey (TÜBITAK) under the Research Grant no. 2211.