Domain Decomposition Methods for Finite Periodic Structures in Electromagnetics

(1)

Master Thesis

Domain Decomposition Methods for Finite Periodic Structures in

Electromagnetics

Gustav Casselbrant

Stockholm, Sweden 2012 XR-EE-ETK 2012:006

(2)

Domain Decomposition Methods for Finite Periodic Structures in Electromagnetics

A Degree Project in Electromagnetic Theory

GUSTAV CASSELBRANT

Master’s Thesis at Department of Electromagnetic Engineering Supervisors: Anders Ellgardt & Henrik Holter

Examiner: Martin Norgren

XR-EE-ETK 2012:006

(3)

(4)

iii

Abstract

Long computational times persist on being a limiting factor when designing complex electromagnetic structures. Even though computers get more and more powerful, the most efficient improvement for reducing the necessary time when solving a numerical problem is to change or optimise the used numerical method. One numerical method introduced lately is the domain decomposition method (DDM).

This thesis explores the available algorithms for computational electromagnetics which are based on DDM. The study is concentrated on algorithms available for computational electromagnetics. The aim of the study is to find and further study one appropriate method for analysing finite periodic electromagnetic structures like frequency selective surfaces and array antennas.

After an introducing study of DDM, the FACTOPO method, is selected and studied in detail. The study includes an implementation of the method, which is performed in three versions. The first deals with a 1-Dimensional waveguide problem which follows the FACTOPO method for determining the scattering parameters of a subdomain using admittance parameters of the corresponding domain. By limiting the subdomain with perfect electric conductor (PEC) the computational problem is efficiently reduced. The second implementation continues to explore the 1-Dimensional problem but introduces the ability of importing scattering parameters of a specific domain from an electromagnetic simulation software. The third, and last, implementation is an extension of the second version which introduces a 2-D structure based on a 4-port unit cell.

The theoretical study and the results from the implementation show that that domain decomposition is a promising numerical technique which, when used properly, can improve the numerical simulation software for electromagnetic structures. It is also concluded that domain decomposition is especially promising for periodic electromagnetic structures.

(5)

Introduction

Even though we today possess powerful supercomputers, complex electromagnetic problems still continue to demand long computational time (days or even weeks).

This is especially true when considering that during design processes for electromagnetic structures, similar problems are often solved repetitively, in order to optimise the design. Over the years a large number of numerical algorithms for reducing the necessary computational time have been introduced. Most of these algorithms have been applied to one out of four major electromagnetic computational methods used today or to a combination of them.

During the last couple of decades numerical algorithms based on domain decomposition have become increasingly more popular. Originally a domain decomposition algorithm (DDM) was a method for solving partial differential equations (PDE) introduced by H. Schwartz in 1870 [36]. Even though the method was not intended as a numerical method, it has recently become a popular choice due to the fact that it is highly appropriate for computational parallelization.

1.1 Purpose

The purpose of this thesis is to summarise the current research and to evaluate the value of using methods based on domain decomposition in order to analyse and design finite periodic electromagnetic structures like frequency selective surfaces and array antennas. The purpose has been divided into the following objectives:

• Perform a theoretical study which summarises the current research status for domain decomposition.

• Select one promising domain decomposition method and implement one simple and one more advanced example of the method using electromagnetic simulation software and matlab.

• Present the project advances and its final result both orally and in writing.

1

(9)

2 CHAPTER 1. INTRODUCTION

1.2 Methodology

The thesis project has been divided into four main phases. Section 1.2.1 - 1.2.4 describe these phases and how the work has been performed within them.

1.2.1 Planning

An important part of the planning was the division of the project purpose into project objectives. The objectives were then further divided into requirements for each phase. Based on the requirements, tasks could be identified and used for producing a time plan. The project plan was the result of the planning phase.

1.2.2 Theoretical Study

The first part of the theoretical study was to study two preparative areas. The first preparative area covered three of the major numerical methods for conducting electromagnetic simulations. These are essential for the thesis since any implementation of domain decomposition in a electromagnetic perspective will be based on one or a combination of these methods. The second preparative area was to study domain decomposition in a general numerical analysis perspective. These two preparative areas were studied mainly using literature.

The second part of the theoretical study concentrated on the different imple- mentations of domain decomposition methods for solving electromagnetic problems.

This part was mainly performed using published articles which were found mainly using Google Scholar. Most articles were available through the IEEE and Science direct databases. During the study it turned out that many methods were mentioned in only a few scientific articles. Therefore the first task was to identify and list all methods mentioned in the articles. The second task was to determine which methods could be categorised as established methods. In order to do so a selection criteria had to be chosen. The criteria taken was that if a method had been mentioned in at least four articles it was to be considered established. Before neglecting any methods, a search for more articles of every method listed was performed. The last step in the study for established methods was to summarise the identified methods.

An important milestone during the literature study was the selection of a domain decomposition method to be used during the implementation phase of the project.

The decision was mainly motivated by the intended application for the final result as specified by the project supervisor Henrik Holter. The motivation of the selection is listed in section 3.2.

The third part of the main study was to comprehensibly describe the selected domain decomposition method. An important part when describing the method was to develop a simple example which could illustrate the different steps in the method.

(10)

1.2. METHODOLOGY 3

1.2.3 Implementation

The implementation of the selected method was performed using Matlab. The development process of the algorithm was more of evolving nature where small changes between the versions lead to an increasing complexity. Therefore a number of different problems, with different complexity have been treated. The report will deal with one analytical example, one 1-D example and one 2-D example. The analytical example exemplifies how the scattering parameters are determined by the FACTOPO method. The two later examples makes use of data exported from CST Microwave Studio.

1.2.4 Reporting

Even though the official reporting phase did not start until the end of the project the reality is that the documentation as well as the reporting have been performed in parallel with the project. Besides the thesis report two shorter reports (included in the thesis report as Chapter 4 & 5) have been produced together with the project plan. The project plan as well as the theoretical study were also presented at the Department of Electromagnetic Engineering at KTH.

During the official report phase the last parts of the thesis report will be finished.

Together with the presentation of the thesis they will mark the end of this thesis project.

(11)

(12)

Chapter 2

Computational Electromagnetics

Traditionally electromagnetic structures were designed using analytical electrody- namic theory [16][26][4]. Analytical solutions methods are still important for un- derstanding the dynamics of simple structures. They are also commonly used for solving problems with non complicated structures. However when the problem structure is more complicated, or when it is necessary to determine the character- istics without simplifications, the analytical formulas get rather complex and time consuming. Therefore it is natural applying numerical methods in order to simulate the design performance.

In the electromagnetic community the numerical methods are commonly referred to as Computational Electromagnetics. In similarity with most other numerical ex- pertise areas the four dominating methods are: the Method of Moments (MoM), the Finite Element Method (FEM), the Finite Difference Time Domain (FDTD) method and the Asymptotic method. These methods are today widely used both within industrial and within academical contexts and constitute the absolute foun- dation of most available commercial simulation softwares. This chapter describes the first three methods and how they are implemented in an electromagnetic perspective. The fourth method, the Asymptotic Method, is a method used for solving very large problem. It is not commonly used together with domain decomposition and will therefore not be covered by this thesis.

2.1 The Method of Moments

Perhaps the Boundary Element Method (BEM) is the traditionally most commonly used numerical method within the electromagnetic community. However within the community BEM is more known as the Method of Moments (MoM). The method was first introduced by R. F. Harington in 1967 [25]. In the article the author presents the method of moments as a concept for reducing the electromagnetic field equations into matrix equations. In order to demonstrate the method the article also includes examples for computation of the electromagnetic fields caused by surfaces and wires.

5

(13)

6 CHAPTER 2. COMPUTATIONAL ELECTROMAGNETICS

In contrast to most other numerical methods the MoM uses integral equations instead of differential equations in order to describe the dynamics of the problem.

The method further utilises the fact that the electromagnetic fields are directly related to the current distribution of an object. Integral equations describing this relationship is achieved by the use of Green’s functions which relates the electric field at an arbitrary observation point to the current distribution at an arbitrary source point [33]. By the use of orthogonal expansion the MoM expresses the integral equations as a system of linear equations. The method further uses test functions in order to invoke the electric field boundary conditions. Using the matrix notation the expansion coefficients associated with the basis functions are determined and thereby the current distribution can be found.

The most well known implementation of MoM is the Numerical Electromagnetic Code (NEC). The software was written in FORTRAN during the 1970:s. The development was led by G. J. Burke and A. J. Poggio at the Lawrence Livermore National Laboratory and was sponsored by Naval Oceans Systems Center and the Air Force Weapons Laboratory. The theory behind as well as the code itself is described in [11] and [12] while the usage of the code is described in [13].

2.2 The Finite Element Method

The development of the Finite Element Method (FEM) is credited to R. L. Courant who added an appendix to one of his papers [19] where he described how variational methods could be used in potential theory. In the description Courant introduced the use of piecewise-linear approximants using sets of triangles. These triangles were refereed to as elements by Courant and have, since then, been characteristic for the method.

The use of finite element method in the computational electromagnetic community began in the late sixties. One early paper, written by P. L. Bahrani [3] in 1968, was written on the subject of treating Helmholtz’s equation using FEM. A concurrent paper was released the year after by P. Silvester [35] on the subject of solving homogeneous waveguide problems using FEM.

Today the use of FEM is well spread in the computational electromagnetic community. A prominent writer in the subject is J. Jin who formulates the basic steps of the algorithm as follows [27]:

1. Make a subdivision of the computational domain. For a 2D problem this is usually done using triangular or rectangular shaped elements.

2. Select the interpolation function. Usually a polynomial function of first, second or higher order is chosen.

3. Formulate the problem through system equations. The two most popular solvers that are used for FEM is the Galerkin’s Method and the Ritz method.

4. Solve the problem by solving the system equations.

(14)

2.3. THE FINITE-DIFFERENCE TIME-DOMAIN METHOD 7

2.3 The Finite-Difference Time-Domain Method

The Finite-Difference Time-Domain (FDTD) method was first introduced by K. S. Yee [56] who presented the algorithm now known as the Yee algorithm. The method was given its current name, Finite Differential Time Domain, as well as its acronym (FD-TD at first but later FDTD) by Allen Taflove in his first article on the regarding the method [40]. According to [34] the Yee algorithm can be summarised as follows:

1. Rewrite Ampere’s and Faraday’s laws to finite difference form.

2. Discretize space and time.

3. Express the future electric and magnetic field using the difference equations based on passed (known) values.

4. Evaluate the magnetic fields one time-step into the future.

5. Evaluate the electric fields one time-step into the future.

6. Repeat the last two steps until the fields are known for the desired time period.

Differential Time Domain (FDTD) it today one of the primary methods in Compu- tational Electromagnetics and is thoroughly described in [41]. Some of the reasons for its popularity are the flexibility together with its easiness to implement. The fact that it operates in time domain means that it provides a solution for all fre- quencies only limited by the length of the sapling period as specified by the Nyquist theorem. This makes it capable of handling impulsive and non linear behaviour in a natural way. On the drawback side is the fact that FDTD has a time step limit which for three dimensions can be expressed as ∆t < h/c √

3, where ∆t is the time step, h is the cell size and c the speed of light. This limitation mainly applies for problems with time scales much longer than the wave propagation time for the region of interest [10].

(15)

(16)

Chapter 3

Domain Decomposition

Domain decomposition methods (DDM) were originally developed for solving partial differential equations (PDE). The methods solve a problem by decomposing the initial computational domain into smaller, less complicated, subdomains which are solved independently. The solution of the full problem is then achieved by assem- bling the solutions of the subdomains. The reason for the recent popularity of the method is that the subdomains are easily solved in parallel. This advantage makes the method suitable of todays computers which in general utilise parallel processors.

Three strong arguments, presented in [36], for using DDM are:

• simpler to parallelize with higher performance than the original problem

• provides simpler and less complicated problems and geometries

• good convergence properties

The earliest known method utilising domain decomposition, the Schwartz alternat- ing method, was introduced by the German mathematician H. Schwartz in 1870.

Though the method was not intended as a numerical method it may be used for solving elliptic boundary conditions. A numerical implementation of the method could be summarised as follows:

1. Start by determining an approximate solution of the problem for the full initial domain.

2. Decompose the initial domain into subdomains.

3. Solve the the problem of each subdomain using the information in the overlapping regions as boundary values.

4. Assemble the solutions of the subdomains.

5. Compare the improvement of the solution compared with the last solution.

9

(17)

10 CHAPTER 3. DOMAIN DECOMPOSITION

6. If the problem has not reach the desired convergence condition, iterate from step 2.

A characterisation between different DDMs can be done by using two fundamental properties. The first property is used if the method is designed as an iterative of direct solvers. Iterative solvers has been given much attention lately since they demand less computational memory and often are more efficient than the direct version. Theses iterative solvers are usually based on the Schwartz methods. The direct methods, also known as substructuring or tearing methods, on the other hand do not demand the full problem to be solved for a change in one subdomain.

An example of a substructuring method is the Schur complement method. It uses the Schur complement, which is a linear algebraic method, in order to solve the subdomain problems. A description of the major steps of the algorithm, based on Algorithm 4.1.1 in [36], follows below:

1. Partition of the structure into nonoverlapping subdomains 2. Compute the Shur complement of each subdomain

3. Merge the Schur complements 4. Factor the Schur complement

The second property of a DDM is the decomposing strategy of the initial domain which can be performed with either non-overlapping or overlapping subdomains.

Just as the naming of the two categories implies, each division is based on how the decomposition of the original domain is performed. Figure 3.1 displays an overlapping decomposition and a non overlapping decomposition of the respective initial domain Ω. The popularity of methods based on the two decomposition strategies has changed back and forth over time. Generally it can be said that the overlapping methods have better convergence properties since neighbouring domain share more information. On the other hand, since subdomains often are solved in parallel and the communication in between the processors often is the limiting factor, a faster convergence does not necessarily mean a faster solution.

3.1 Domain Decomposition and Electromagnetics

During the years a number of techniques involving Domain Decomposition for solving numerical electromagnetic problems have been presented. An early attempt was made by A. Taflove and K. Umashankar in 1989 [45] where they introduced spatial decomposition as an approach to reduce the dimensions of the matrices involved when solving MoM systems. They also emphasised the improved condition of the resulting matrices which also contributed to the reduction of the computational bur- den. This early attempt was further followed by additional articles [44, 46] which

(18)

3.1. DOMAIN DECOMPOSITION AND ELECTROMAGNETICS 11

Ω₁ Ω₂

(a) A non overlapping decomposition, Ω = Ω1∪ Ω2, Ω1∩ Ω2= 0

Ω₁ Ω₂

(b) A overlapping decomposition, Ω = Ω1∪Ω2, Ω1∩Ω26= 0

Figure 3.1: Examples of a overlapping respectively a non overlapping decomposition of an initial domain.

continued to explore the possibilities of spatial decomposition together with the method of moments.

During the same time frame B. Despres also started to explore the possibilities of using Domain Decomposition in order to solve numerical electromagnetic problems.

His first article on the subject [20] was released in 1991. It was later continued in a second article [21] released in 1993. In these articles he deals with solving the Helmholtz problem using a DDM. The articles also introduce the use of Robin type transmission conditions, a method which later would be followed by many others [38, 18, 28, 32, 39, 31] Between the two articles B. Despres, together with a couple of co-writers, released an article [22] in 1992 where they dealt with the harmonic Maxwell equations. The method introduced here is a direct development of the one B. Despres earlier introduced in his article.

Today there are a number of domain decomposition methods used for solving numerical electromagnetic problem. Below follows a summary of the more common ones.

3.1.1 Complex Source Expansion

The Complex Source Expansion (CSE) was introduced in an article in 2008 [15]

by G. Carli et. al. The method is a development of the Complex Source Point Expansion (CSP) which in turn is based on the concept of CSP beams [24]. The CSP beam concept is established as a way of constituting an exact solution of the Maxwells equations. The CSP method was first described by K. Tap et. al. in 2007 [43] and then shortly after in a second article regarding the subject [42]. The CSP method divides the object surface in subsectional basis function elements which in turn are grouped together. The interactions between close neighbour groups as well as elements within the groups are then determined using conventional e.g. FEM or MoM. The interaction between well separated groups are described using complex source point beams [24]. This means that the fields of each basis functions enclosed

(19)

group are expanded. By then using the surface equivalent theorem an equivalent spherical surface can be used to analytically describe the group as a complex point value.

The developments of the CSP method introduced by G. Carli et. al. in the CSE method provide a technique which is more adaptable to electromagnetic structures. The writers have continued to develop the method further in [14, 30]. The method decomposes a complex domain into simpler interacting subdomains. The interactions are described by using a network formalism with CSP beams as prop- agators. The interactions inside each subdomain are described by using scattering matrices. This further provides that the interactions between all subdomains and the source field can be described by correctly combining the scattering matrices and the source vector. One advantage with the method is that the subdomains can be characterised independently once and for all. This means that different interactions between the subdomains can be described by correctly combining the scattering matrices. Another advantage is that the interactions within the subdomains can be determined by using different electromagnetic computational methods.

3.1.2 Factorisation Domain Decomposition Method

The factorisation domain decomposition method (FACTOPO) was first released by A. Barka et. al. in 2000[9]. The development of the method had primarily been made at the Electromagnetism and Radar Department of the Office National d’Etudes et de recherches Aerospatiales (ONERA), Toulouse. The basic was to develop a method for collaborative studies involving independent manufacturers.

In the article A. Barka et. al. describe the method when used for determining the radar cross section for a fighter jet engine intake. This usage of the method is later also covered in [37].

The method works by dividing the initial domain into multiple subdomains. The electromagnetic fields at each interface can then be expanded using basis functions in the form of waveguide modes. By characterised each subdomain, using scattering matrices, the interactions in between subdomains and the exterior region can be determined by combining the scattering matrices into a Global Scattering Matrix (GSM).

In [5], released 2005, two of the original writers describe how the method can be generalised and thus used as a solver for large electromagnetic problems. The article describes how the method can be used to optimise the positioning of a monopole antenna on a Fokker 100. This is achieved by utilising the advantage of the fact that most of the subdomains remain unchanged and therefore do not need to be recalcu- lated when the antenna position is changed. The generalisations of the method are achieved by the use of deltagaps as basis functions for expanding the fields on the interfaces. The method is further presented in a more comprehensive manner in [6].

In [7] the method is used in order to evaluate the performance of a pod antenna installer on a motor glider.

An interesting development of the method is presented in [8]. In this article

(20)

3.1. DOMAIN DECOMPOSITION AND ELECTROMAGNETICS 13

A. Barka and F. X. Roux describe how the FETI-DP algorithm can be used to determine the scattering parameters of the subdomains.

As can be noted the FACTOPO method is very similar to CSE. The most significant difference is the choice of basis functions describing the relationship in between the domains. In the CSE case the choice of CSP beams makes the method appropriate for problems with well separated subdomains. FACTOPO on the other hand can make use of many different basis functions, which makes it appropriate for problems where the subdomains are more tightly assembled.

3.1.3 Domain Decomposition - Finite Difference Time Domain

The method of Domain Decomposition - Finite Difference Time Domain (DD- FDTD) was presented by F. Xu and W. Hong in 2004 [52]. In the first article the writers introduce the method for analysing complicated horn antennas. The method works by decomposing the structure into overlapping subdomains. The computation of each subdomain is carried out independently by the use of a local mesh and time step. The exchange of data in between the regions is performed in the overlapping regions. Since each region has its local mesh and coordinate system, the meshes do not coincide the data exchange that is performed by interpolating it from one region to another. An advantage with the method is that it provides the possibility of having non uniform meshes. In this way it is possible to catch the electromagnetic fields of fine structural details without the necessity of using the same fine mesh over the full structure.

In a second article [51], also released in 2004, F. Xu and W. Hong develop their method further by presenting an additional application. Here they make use of a concept for directly computing the transient far field response of a two dimensional structure, earlier described in [53], in order to create an effective method of analysing scattering problems for two dimensional structures. This method is later continued in[54] where a version for three dimensional structures is presented.

The next application of the method is presented in 2006 [55]. Here F. Xu and W. Hong together with K. Wu present a hybrid method, based on DD-FDTD, for analysing substrate integrated circuits.

3.1.4 Finite Element Tearing and Interconnecting

The Finite Element Tearing and Interconnecting method (FETI) was first introduced in 1991 [23] by C. Farhat and F. X. Roux. The two writers were both involved in the development of a FEM DDM within the mechanical engineering community.

The general idea was to create a parallel method where each sub domain could be assigned a separate processor. This meant that each domain could be solved quickly.

However they saw that the bottleneck with these methods was the slow information exchange in between the processors during the computation. In the article they suggested a FEM DDM [1] called FETI, that uses hybrid variational principle

(21)

based on Lagrange Multipliers. By this they achieved finite element domains that are disconnected and also avoided the problem with ”floating” domains.

The explorations of the possibilities using FETI algorithm were continued by the mechanical community and were later adopted by the electrophysics community.

An early attempt was made by C. T. Wolfe et.al in 2000 [50]. They introduced the method as a way of solving the vector wave equation and used perfectly matched layer (PML) in order to limit the problem domain.

The algorithm proposed by Wolfe can be summarised in the following steps:

1. Divide the computational domain into nonoverlapping subdomains 2. Discretisize the electric field testing functions in each subdomain

3. Solve each subdomain using Lagrange multipliers by enforcing the continuity of the tangential field across each boundary

4. Compute the electric field of each subdomain

5. Gather the subsolutions to determine the resulting fields

The FETI algorithm proved to be effective but could cause floating subdomains [1]. To deal with these early problems the FETI algorithm has been developed further. One of the more common version today in the electromagnetic community is the DP-FETI (which also may be referred to as FETI-DP and FETI- DPEM) [49, 47, 59, 29, 58]. The algorithm was presented in 2004 [49] and introduces the Robin condition across the interfaces as well as the mortar method for non-confirming triangulations.

Besides the DP-FETI algorithm, a very promising version of the FETI algorithm is the method presented in [48]. In similarity to DP-FETI this method also makes use of ROBIN transmission conditions at the interfaces. However it also introduces a ”set of ”cement” unknowns on the inter-domain interfaces” [48] These unknowns allow the grid to be non-conforming/nonmatching. The method is further developed in [57] which also makes use of the results presented in [28].

3.2 Selecting a Domain Decomposition Method

The method chosen for implementation is the FACTOPO algorithm. The motivation why can be summarised as follows.

The FACTOPO method

• provides a solution of the system dynamics which is independent of the excitation since the excitation is not introduced the resulting resulting fields pattern for a specific excitation is desired.

• does not demand the solving of the full initial domain if a subdomain is changed.

(22)

3.2. SELECTING A DOMAIN DECOMPOSITION METHOD 15

• can make use of different numerical solvers for different subdomains.

• is highly suitable for small near field problems.

The first two arguments are important since they make it possible to make changes in a design, in order to find an optimum, without having to solve the full problem. The second is a practical one because it makes the method independent of the solver. This is an important feature since different solvers are suitable for different parts in a design. The last motivation differs FACTOPO from the CSE method due to the fact that CSE is more suitable for far field dependencies. This is a direct consequence of CSE using source point expansion which makes use of traditional methods for determining the close range dependencies.

(23)

(24)

Chapter 4

FACTOPO - The Selected Method

4.1 Introduction to the FACTOPO Method

As earlier described in Section 3.1.2 the fundamental part of the FACTOPO method is the possibility to characterise each subdomain using scattering matrices. This further provides the possibility of describing the subvolume’s interactions with each other as well as with the surroundings by the use of a global scattering matrix (GSM). The main steps in the algorithm are as follows:

1. Decompose the initial computational volume into sub volumes 2. Determine the scattering matrices for each of the sub volumes 3. Assemble the scattering matrices in to a GSM

4. Solve the global problem by inverting the GSM

An important feature of the method is how the scattering matrices of the subdomains are determined. This is done by first determining the admittance matrix of the corresponding domain while limiting the domain using PEC. This limits the computational volume for each subdomain effectively. The admittance parameters are then converted to scattering parameters. The remainder of this chapter deals with the described main steps of the algorithm in a more comprehensive manner.

4.2 Notations and Indexes

Notation Explanation

Ω Initial computational domain

P Number of computational subdomains Continued on the next page

17

(25)

18 CHAPTER 4. FACTOPO - THE SELECTED METHOD

Notation Explanation

V_i Subdomain i

M Number of interfaces

Γ The total interface surrounding the current domain

Γ_j A part of the total Interface Γ corresponding to the interface towards a specific neighboring subdomain

Γ₀ The perfectly conducting part of the total interface Γ N Number of interfaces for a certain domain

E~^j The tangential electric field of interface Γ_j H~^j The tangential magnetic field of interface Γ_j

N_j The number of basis functions expanding the fields of the interface Γ_j

a^j_n The outgoing wave at interface Γ_j b^j_n The ingoing wave at interface Γ_j

~g_n^j The magnetic basis function n at interface Γ_j

~f_n^j The electric basis function n at interface Γ_j

H~_n^{j ex} The magnetic field induced in the domain by the electric function ~f_n^j at one ot interfaces in the domain.

~

n_j Normal vector to interface Γ_j, is always directed inwards with aspect to the current sub domain.

pn Port n.

4.3 Decomposition of Computational Volume

The decomposition of the initial volume into sub volumes may be directly correlated to the physical division with the initial system in physical parts. However if the resulting subvolumes are complicated, the method also allows a further decomposition of the subvolumes in order to achieve subvolumes of equal complexity. Figure 4.1 displays an example of how the decomposition could be performed.

(26)

4.4. DETERMINING THE SCATTERING PARAMETERS OF A SUBDOMAIN 19

V₄ V₃ V₂

V₁

V₀

p1

p2

p3

p4

p5

p6

p7p8

Γ1 Γ2 Γ3

Γ4

Figure 4.1: A decomposition example for a horn antenna with a monopole feed placed in a cavity and a waveguide between the cavity and the horn.

The structure shows a horn antenna with a dipole feed, placed in a cavity, and a waveguide between the feed point and the horn. The initial volume Ω is decomposed into the sub volumes V₀-V₄. The horn, waveguide and cavity correspond to sub volumes V₁-V₃ while V₀ corresponds to the exterior region and V₄ to a fictive region representing the feeding of the monopole.

The volumes are connected by ports (p₁− p₈) at each of the interfaces Γ₁-Γ₄. The ports on the virtual interface Γ₄ correspond to a coaxial feed connected to the monopole.

The tangential electric (~E ) and magnetic field ( ~H ) of interface Γ_j are described according to eq. 1 in[6] as:

E~^j =

Nj

X

n=1

a^j_n+ b^j_n~f_n^j (4.1)

H~^j =

Nj

X

n=1

a^j_n− b^j_n~g^j_n (4.2)

Where N_j is the number of basis functions at interface Γ_j, a^j_n and b^j_n represent the amplitude of the inward and outward travelling waves respectively and ~f_n^j ~g_n^j the electric respectively the magnetic basis functions.

4.4 Determining the Scattering Parameters of a Subdomain

The scattering parameters for a specific subdomain can be determined from the admittance matrix of the same subdomain. The process of determining the admittance matrix for a subdomain is shown in Section 4.4.2 - 4.4.3. Thereafter follows the procedure of determining the scattering matrix based on the admittance matrix in Section 4.4.4. The method described follows the steps described by sections II.b - II.c in [6].

(27)

In order to visualise how the scattering parameters are determined an example involving the waveguide volume V₂, shown in Figure 4.1, will be used. This example will be referred to as the specific case. Each step of the procedure will first be described in a general case, followed by the same step for the specific case. The layout of the specific case is described in Section 4.4.1. The method directly follows the steps described by sections II.b - II.c in [6].

4.4.1 Specific Case Layout

The specific case consists of a waveguide which is a part of the system displayed in Figure 4.1. The waveguide is displayed in detail in Figure 4.2. The walls of the waveguide are perfect conductors. At the ends the waveguide has two interfaces, Γ₂ & Γ₃. The tangential electromagnetic fields at interface Γ₂ are expanded by three basis functions while interface Γ₃ is expanded by two basis functions. For this implementation waveguide modes have been chosen as basis functions according to [9].

xˆ yˆ zˆ

w

V2 h

Γ₂ Γ₃

Figure 4.2: A waveguide with two interfaces.

The example layout is defined in such a way that only transverse electric (T E) modes exist through the waveguide. Later in the implementation of the example, objects will be placed within the waveguide. These objects will however be placed in such a way that the TE modes are maintained and no transverse magnetic (TM) modes exist. The electric and magnetic modes inside the waveguide for the T E case are described by:

Em0,y(x) = sin

mπ w x

H_m0,x(x) = 1

Z_m0n_j× E_m0,y m = 1, 2, 3, ...

The impedance Z_m0 according to table 10.1 in [16] is given by:

Zm0 = η r

1 −^f_f^c²

(4.3)

(28)

Where η is the intrinsic impedance of the guide medium while f_c is the cut of frequency of the waveguide.

In the specific case the basis functions are chosen as:

~f_m^j = f_m^jˆy = E_m0,y (4.4)

~g_m^j = g_m^j x = Hˆ _m0,x (4.5) 4.4.2 Determining the Induced Magnetic Field

The process of determining the admittance matrix of a volume uses an expression which determines the magnetic field at interface Γ_i induced by the electric field in the volume V_j. This expression makes use of two different methods of determining the magnetic field which in turn are set to equal each other.

Determining the Magnetic Field, Method 1 General Case

The tangential magnetic field at interface Γ_jis given by Equation (4.2). This further provides that the surface current can be expressed as:

nˆ_j× ~H^j =

Nj

X

n=1

a^j_n− b^j_nnˆ_j× ~g_n^j (4.6)

Where ˆn is the normal vector to interface Γ_i inward to volume V_j. Specific Case

The surface current of the two interfaces, limiting the volume V₂, in the specific case can be determined by expanding the sum in Equation (4.6). For the interfaces Γ₂ and Γ₃ the resulting components are as follows:

Γ₂: ˆn₂× ~H² =a²₁− b²₁nˆ₂× ~g₁²+a²₂− b²₂ˆn₂× ~g₂²+a²₃− b²₃nˆ₂× ~g₃²

Γ₃ : ˆn₃× ~H³=a³₁− b³₁nˆ₃× ~g₁³+a³₂− b³₂nˆ₃× ~g₂³ Determining the Magnetic Field, Method 2

General Case

It would be desirable to find a method which expresses the magnetic field of interface Γ_i as a function of the electric field in the volume V_j. This can be achieved by a numerical calculation of the resulting H-field ( ~H_n^{j ex}) at Γ_iinduced by the respective

(29)

electrical basis function (~f_n^j), from Equation (4.1). Hence the magnetic field at interface Γ_i can be written as:

nˆ_j× ~H^j =^X

i Nj(i)

X

n=1

aⁱ_n+ bⁱ_nˆn_j × ~H_n^{i ex} (4.7) where i represents an iteration over all interfaces in the current volume. The equa- tion corresponds to the double sum formulation found at eq. 2 in [6].

Specific Case

The electric field of the two interfaces, based on Equation (4.1) are expressed as:

Γ₂ : ~E² =

N2=3

X

n=1

a²_n+ b²_n~f_n²

Γ₃ : ~E³ =

N3=2

X

n=1

a³_n+ b³_n~f_n³

By expanding equation (4.7) the two interfaces Γ₂ & Γ₃ can be expressed by their respective components:

Γ₂ : ˆn2× ~H² =a²₁+ b²₁nˆ2× ~H₁^{2 ex}+a²₂+ b²₂nˆ2× ~H₂^{2 ex}+a²₃+ b²₃nˆ2× ~H₃^{2 ex} +a³₁+ b³₁nˆ2× ~H₁^{3 ex}+a³₂+ b³₂nˆ2× ~H₂^{3 ex}

Γ₃ : ˆn3× ~H³ =a²₁+ b²₁nˆ3× ~H₁^{2 ex}+a²₂+ b²₂nˆ3× ~H₂^{2 ex}+a²₃+ b²₃nˆ3× ~H₃^{2 ex} +a³₁+ b³₁nˆ₃× ~H₁^{3 ex}+a³₂+ b³₂nˆ₃× ~H₂^{3 ex}

Equal the two Methods of Determining the Magnetic Field General Case

Once the fields have been expressed in components it is possible to summarise the contributions from each component in order to construct the total field applied on the interfaces. This procedure corresponds to the complete Equation 2 found in [6]:

ˆn_j × ~H^j =

Nj

X

n=1

a^j_n− b^j_n ~n_j× ~g_n^j = ^X

i Nj(i)

X

n=1

aⁱ_n+ bⁱ_n~n_j× ~H_n^{i ex} (4.8)

By the Equation above an expression has been found which describes the magnetic field at interface Γ_i induced by the electric field in the volume V_j.

(30)

Specific Case

For the specific case the resulting expressions for the interfaces Γ₂ & Γ₃ look as follows:

Γ₂: ˆn₂× ~H² =a²₁− b²₁nˆ₂× ~g₁²+a²₂− b²₂ˆn₂× ~g₂²+a²₃− b²₃ˆn₂× ~g₃²

=a²₁+ b²₁nˆ₂× ~H₁^{2 ex}+a²₂+ b²₂ˆn₂× ~H₂^{2 ex}+a²₃+ b²₃nˆ₂× ~H₃^{2 ex} +a³₁+ b³₁nˆ₂× ~H₁^{3 ex}+a³₂+ b³₂nˆ₂× ~H₂^{3 ex}

(4.9)

Γ₃: ˆn3× ~H³ =a³₁− b³₁nˆ3× ~g₁³+a³₂− b³₂ˆn3× ~g₂³

=a²₁+ b²₁nˆ3× ~H₁^{2 ex}+a²₂+ b²₂ˆn3× ~H₂^{2 ex}+a²₃+ b²₃nˆ3× ~H₃^{2 ex} +a³₁+ b³₁nˆ₃× ~H₁^{3 ex}+a³₂+ b³₂nˆ₃× ~H₂^{3 ex}

(4.10) 4.4.3 Introducing the Test Function and Scalar Product

When the magnetic field has been described the next step, according to [6], is to evaluate the resulting field by applying a family of test function ~t_m^j = ~f_m^j on each of the interfaces Γ_j. These are introduced together with the scalar product:

h ~a, ~b i = Z

Γj

~a · ~b^∗dS

As described by the Gallerkin procedure, the resulting field will be given by the weakly formulated scalar product when applying the test functions ~t_m^j = ~f_m^j = sin ^mπ_w x: to Equation (4.8).

Nj

X

n=1

a^j_n− b^j_nh ~n_j× ~g_n^j, ~f_m^j i = ^X

i Nj(i)

X

n=1

aⁱ_n+ bⁱ_nh ~n_j × ~H_n^{i ex}, ~f_m^j i (4.11) The formulation described in Equation (4.11) will be treated in two steps corresponding to the left and right side of the expression.

Treating the Left Hand Side of Equation (4.11) General Case

The general solution for each of the basis functions is given by:

a^j_n− b^j_nh ~n_j× ~g_n^j, ~f_m^j i = a^j_n− b^j_n Z Z

Γj

~ n_j × ~

g^j_n

~

f_m^j dS (4.12)

(31)

Specific Case

For the specific case each basis function can be treated in the following way

a^j_n− b^j_n Z Z

Γj

~ n_j× ~

g_n^j

· ~

f_m^j dS = a^j_n− b^j_n Z Z

Γj

ˆ

n_j× g_n^jˆx·f_m^j ˆy dS (4.13) For m 6= n the result of the integral in Equation (4.13) equals 0. For the case when m = n and j = 2 → ˆn_j× ˆx = ˆz × ˆx = ˆy the result is:

a^j_m− b^j_m Z Z

Γj

g_m^j · f_m^j dS = 1 Zm

a^j_m− b^j_m Z h

0

Z w 0

sin²

mπ w x

dxdy

=







s = ^mπ_w x dx = _mπ^w ds







= 1 Zm

a^j_m− b^j_mwh mπ

Z mπ 0

sin²(s) ds

= 1 Zm

a^j_m− b^j_m wh

2mπ[s − sin(s) cos(s)]^mπ₀

= 1 Zm

a^j_m− b^j_mwh 2

(4.14) Please observe that the solution makes use of the fact that the field direction in the example is constant since the layout is defined for T E wave propagation. This means that the resulting five cases in the specific case are given by :

a²₁− b²₁h ~n₂× ~g₁², ~f₁²i = a²₁− b²₁wh 2

a²₂− b²₂h ~n₂× ~g₂², ~f₂²i = a²₂− b²₂wh 2

a²₃− b²₃h ~n₂× ~g₃², ~f₃²i = a²₃− b²₃wh 2

a³₁− b³₁h ~n₃× ~g₁³, ~f₁³i = a³₁− b³₁wh 2

a³₂− b³₂h ~n₃× ~g₂³, ~f₂³i = a³₂− b³₂wh 2

These five terms correspond to the five diagonal terms when expressing the left hand side of Equation (4.11) using matrix notation.

Treating the Right Hand Side of Equation (4.11) General Case

Now, treating the right side of the expression given in Equation (4.11), a general solution for the different terms can be written as:

Domain Decomposition Methods for Finite Periodic Structures in Electromagnetics

Domain Decomposition Methods for Finite Periodic Structures in

Electromagnetics

Gustav Casselbrant

Domain Decomposition Methods for Finite Periodic Structures in Electromagnetics

Contents

Chapter 1

Introduction

1.1 Purpose

1.2 Methodology

Chapter 2

Computational Electromagnetics

2.1 The Method of Moments

2.2 The Finite Element Method

2.3 The Finite-Difference Time-Domain Method

Chapter 3

Domain Decomposition

3.1 Domain Decomposition and Electromagnetics

3.2 Selecting a Domain Decomposition Method

Chapter 4

FACTOPO - The Selected Method

4.1 Introduction to the FACTOPO Method

4.2 Notations and Indexes

4.3 Decomposition of Computational Volume

4.4 Determining the Scattering Parameters of a Subdomain