# T-Ω Formulation with Higher-Order Hierarchical Basis Functions for Nonsimply Connected Conductors

Authors: Ahmed Khebir; Paweł Dłotko; Bernard Kapidani; Ammar Kouki; Ruben Specogna
Journal: Mathematical Problems in Engineering
 (2018)

---

## Abstract

This paper presents in detail the extension of theT-Ω formulation for eddy currents based on higher-order hierarchical basis functions so that it canautomatically deal with conductors of arbitrary topology. To this aim, we supplement the classical hierarchical basis functions with nonlocal basis functions spanning the first de Rham cohomology group of the insulating region. Such nonlocal basis functions may be efficiently and automatically found in negligible time with the recently introduced Dłotko–Specogna (DS) algorithm. The approach presented in this paper merges techniques together which to some extent already existed in literature but they were never grouped together and tested as a single unit.


---

## Body

## 1. Introduction

This paper considers the problem of computing the so-called eddy currents [1] that come up in an electrically conducting body when it is immersed in a slowly varying magnetic field h. Electromagnetic phenomena are governed by Maxwell’s equations [2] and material constitutive relations. For slowly time-varying fields whose change in magnetic field energy is dominant and electromagnetic wave propagation can be ignored, it is wise to neglect the displacement current density in the Ampère–Maxwell’s equation [2]. Thus, the equations to be solved in the time-harmonic case [1] are(1)curl⁡h=σe,(2)curl⁡e+iωμh=0,where e is the electric field, σ is the electric conductivity (that vanishes inside the insulators), μ is the magnetic permeability, ω is the angular frequency, and i is the imaginary unit.In this paper we focus on the widely usedT-Ω Finite Element formulation [1] for solving eddy currents problems. This formulation is hardly the only choice for the solution of this particular physical problem; nevertheless it is one of the most popular and computationally attractive solutions given that it reduces the number of unknowns by exploiting the fact that the magnetic field is irrotational in the insulating region (where σ vanishes).Let us assume that the conductors form a combinatorial three-dimensional manifold [3] Kc and that K=Kc∪Ka, where Ka represents the insulating region and K the whole computational domain. To make the presentation easier, we assume that K is simply connected. This is the case in the great majority of practically relevant problems. Furthermore this assumption does not affect the generality of the results of this paper, given that it may be removed as described in [4].If one wants to solve the eddy current problem with the magnetic scalar potentialΩ in configurations containing topologically nontrivial conductors (e.g., conductors with “handles” like a torus), representatives of first cohomology group generators of insulator are the most efficient mean to make the problem well defined. They are typically addressed ascuts in the community.A common misconception is that the problem of finding cuts has been solved many years ago. Yet, there seems to be not a single proposed solution available in literature which is fully automatic (in the sense that the user needs not to provide any homological information as input), provably correct, and acceptably fast. To back up this statement let us review the approaches available in the literature so far.One pioneering attempt to solve this issue for the lowest orderT-Ω Finite Element formulation was proposed in [5, 6]. In these papers Kotiuga proposed various algorithms to produce a set ofcuts, whose elements are representatives of generators of a second relative homology group basis [3] of Ka realized as discrete surfaces in Ka. Yet, those are restricted to non-self-intersecting surfaces, which is a computationally expensive choice. As a result, the algorithms, although a considerable progress, are not to be considered practical for industrial size engineering problems. Incidentally, the drawbacks of these approaches is one of the main reasons why the authors chose to use cohomology [3] theory instead.In principle, one can use rigorous methods from standard algebraic topology, like the ones presented in [7], to compute representatives of a cohomology basis. These approaches are based on reduction of the complex Ka followed by a Smith Normal Form (SNF) computation on the boundary matrix of the reduced complex. Unfortunately, this route to finding cohomology generators easily leads to the solution of the topological problem to be the bottleneck of the whole simulation (see [8] for a recent comparison); therefore such methods to solve this issue are not used in electromagnetic solvers.On the other hand, there is a large body of literature in which the problem of dealing with nonsimply connected conductors is claimed to be solved by using heuristic approaches. Most of these papers present a general formulation for nonsimply connected conductors, but the required topological preprocessing is not actually completely dealt with. For example in [9] the unknown source fields which have to be precomputed to use the formulation must be still provided by hand. In fact, in the vast majority of papers, an algorithm to produce these unknown source fields is omitted, and even a precise mathematical definition of what these fields are is totally missing. Some papers like [10, 11] propose a partial solution, in the sense that the computation of the unknown source fields is addressed. Yet, they still expect to have some topological information as input. In particular, the algorithm of both papers assumes to have one current distribution as input for each unknown source field to be computed. It is easily seen from nuclear fusion engineering applications, in which the conductor is a n-fold torus with n being O103, that finding n of such current distributions is a very hard problem by itself which is not addressed in any paper. Addressing this issue is a computationally nontrivial task: since we assume that the meshes of the conductor and of the insulator are such that the sum of them is topologically trivial, the so-called exact sequence [12] of the pair of homology groups of insulating and conducting complexes yields to a one-to-one correspondence between the topology of the conductor and the topology of the insulator. In other words, finding independent currents to compute the unknown source fields is to some extent as difficult as finding the unknown source fields. In simple academic cases (like a simple coil, i.e., a torus) those independent currents can be found by hand or by solving a direct-current (DC) problem, but for more complicated cases (as the aforementioned n-fold torus) a dedicated algorithm is needed.In this paper, we propose to perform the topological preprocessing with the DS algorithm (see [4, 13]) which is, to the best of the authors knowledge, the first and so far the only fully automatic, provably correct, and fast way of obtaining cohomology generators. Yet, the set of representatives provided by the DS algorithm is alazy cohomology basis [4, 13]: the provided set of representatives span the needed cohomology group but contains additional, dependent elements. The size of the lazy basis is no more than twice the size of a standard cohomology basis and with moderate effort one may produce a standard cohomology basis (see [4] or thelean cohomology computation described in [14]) given a lazy one. However, it has been verified that these techniques which produce a standard cohomology basis do not provide any speedup in the solution of the electromagnetic problem while giving exactly the same solution in terms of induced currents up to linear solver tolerance. We remark that the linear system matrix is singular when a set of lazy generators is used and therefore a linear solver able to deal with singular systems is required. We used the Intel MKL Pardiso which is able to solve general complex and symmetric systems (by setting thematrix type input parameter to 6).The DS algorithm was originally developed in the framework of the Discrete Geometric Approach (DGA) described in [15]. Yet, both the DGA and the standard Finite Element counterpart based on Whitney edge elements [16] provide a current density which is only uniform inside each mesh element. Higher-order basis functions are therefore very attractive since they yield greater accuracy for a given computational cost and smoother current density vector field. There are various possibilities of obtaining a high order of convergence. The hierarchical basis functions introduced in [17, 18] are particularly appealing. They allow having a good control over the distribution of degrees of freedom (dofs) given that different orders can coexist on the same mesh. Moreover, these different orders may be set by hand or, as in p-adaptivity, automatically by a mesh refinement scheme. Furthermore, these hierarchical higher-order basis functions have been applied to the solution of eddy currents problems using the T-Ω Finite Element formulation in [19]. Yet, the authors of [19] assume that the conducting region is simply connected.In this contribution we will supplement the basis functions introduced in [19] with nonlocal basis functions spanning thefirst de Rham cohomology groupHdR1(Ka) [3] of the insulating region.We were able to find in literature the closest attempt to solve the problem reported in [20]. In particular, the main contribution of [20] and our paper is exactly to combine the use of cohomology theory [3] to the well-known T-Ω Finite Element formulation using hierarchical higher-order basis functions. Yet, in Section 2.4, we point out several shortcomings for the solution proposed in [20], which is based on the algorithm introduced in [16] and also used in [21], to demonstrate that the solution of the problem addressed in this paper is not yet documented in literature. Different from [20], in this paper, we propose to perform the topological preprocessing with the DS algorithm [4], whose typical computational complexity is linear and its correctness is proved for any possible input.The paper, which details the material presented in [22], is organized as follows. In Section 2 we introduce the novel T-Ω formulation employing a lazy cohomology basis together with higher-order hierarchical basis functions to solve eddy current problems with conductors that have arbitrary topology. In the same section, we also survey the algorithm used to extract the lazy cohomology basis. Section 3 presents the results obtained when solving two Testing Electromagnetic Analysis Methods (TEAM) problem benchmarks. Finally, in Section 4, the conclusions are drawn.


## 2.T-Ω Formulation with Higher-Order Hierarchical Basis Functions and Nonlocal Lazy Cohomology Basis Functions

### 2.1. De Rham Cohomology

IfKa is simply connected (or topologically trivial), as in Figure 1(a), a curl-free vector field, like the magnetic field h in an insulator, can be expressed as the gradient of the magnetic scalar potential Ω.Figure 1
The insulating regionKa is obtained as the complement with respect to a box (not represented in the picture for clarity) of (a) a trivial conductor, (b) a solid torus, (c) a twofold solid torus, that is, a conductor with two handles.
(a)
(b)
(c)Topology starts playing a role in theT-Ω formulation when Ka is not simply connected. Let us consider an example where Ka is the complement of a conductive torus Kc with respect to a box as in Figure 1(b). The magnetic field h is curl-free in Ka but its circulation on the loop c of Figure 1(b), because of Ampère’s law, has to match the electric current ic that flows around Kc. The consequence is that the gradient of the magnetic scalar potential is not enough to span the space where the magnetic field h lives.In the following we use concepts of algebraic topology that due to the limited space cannot be reproduced here. Please consult [23] for a formal introduction or [4, 6, 15, 24] for an informal one.In algebraic topology, thefirst de Rham cohomology groupHdR1(Ka) of Ka, is exactly the space of curl-free vector fields in Ka that are not gradients [3]. Therefore, by its very definition, HdR1(Ka) is the space we have to add to the space generated by the well-known hierarchical basis functions described in [25]. It is known that the dimension of this space is equal to the first Betti number β1(Ka) of Ka [3].Focusing again on the example in Figure1(b), β1(Ka)=1, thus the basis for the HdR1(Ka) space is composed by just one element calledgenerator. Let us construct the generator g∈HdR1(Ka) as a curl-free field that has a circulation 1 on a loop linking the conductor as c in Figure 1(b). If one now considers the product icg∈HdR1(Ka), by varying the independent current ic one is able to span the whole HdR1(Ka) space. Therefore, ic is a new degree of freedom that has to be added as an additional unknown of the eddy current problem.This property is generalizable to more complicated topologies; see [4, 26] for a more formal explanation. Here we informally show what happens in a more complicated example, that is, let us consider a twofold solid torus as Kc; see Figure 1(c). Here, β1(Ka)=2 (the conductor has two handles); thus there are two independent currents ic1 and ic2. The HdR1(Ka) space is generated by ic1g1+ic2g2, where the two generators g1,g2 are such that(3)∫cigj·t^dl=δij, where t^ is the tangent vector of the loop ci and δij is the Kronecker delta.In general, the magnetic fieldh is represented as(4)h=∇Ω+∑k=1β1Kaickgk, where the space of gradients of the magnetic scalar potential Ω has been already defined in [25] and each independent current has to be included as an additional degree of freedom to span the whole space of the eddy current problem solution.


### 2.2. Construction ofHdR1(Ka) Generators

Let us try to replicate the properties ofg in the example of Figure 1(b) at the discrete level. Let us call g the discrete counterpart of g. First, the circulation of g on a loop c made of edges that links the conductor like c has to be 1 as for g. Therefore, it is natural to define g as a discrete field having integer coefficients assigned to mesh edges. Second, we want g to be curl-free as its continuous equivalent. This translates in the discrete case to the condition Cag=0, where Ca is the face-edge incidence matrix of the mesh restricted to Ka. That is, the discrete circulation of g over the edges in the boundary of each face of the mesh must be zero. It is possible to prove that such edge discrete fields may be interpreted as the representatives of the generators of the first cohomology group over integers H1(Ka,Z) [4, 26]. Then, to retrieve g, we just need to interpolate g with the standard Whitney edge elements basis functions [1]. Thus, it is evident that there is no need to use higher-order basis functions to span the de Rham space.H
1
(
K
a
,
Z
) are realized using thethick cut technique [28]; see Figure 2.Figure 2
(a) Let us consider a cross section of the example in Figure1(b). (b) The thick black edges belong to the thick cut g (i.e., to the support of the representative of the cohomology generator H1(Ka,Z)). The coefficients assigned to thick edges are such that Cag=0, where Ca is the restriction of C to Ka, and that the circulation on all discrete cycles linking the conductive torus like c is 1.
(a)
(b)The advantage in switching to cohomology over the integers is that the generators ofH1(Ka,Z) can be rigorously and efficiently computed. In this respect, it is important to note that the property of the basis functions of being hierarchical is of fundamental importance: thus, the term of h related to H1(Ka,Z) can be accounted for when using higher order by adding it to the scalar first-order dofs on the edges of Ka. In what follows we give a detailed description of an algorithm that performs this task.


### 2.3. Efficient Computation ofH1(Ka,Z) Generators

We find the required nonlocal basis functions spanning the first de Rham cohomology group space by exploiting the recently introduced idea oflazy cohomology generators [4]. In what follows we recall the DS algorithm that finds the lazy generators of the first cohomology group H1(Ka,Z):(1)
First the discrete surfaceS=Kc∩Ka (see Figure 3(a)) is extracted and its first cohomology group generators H1(S,Z) are computed with a combinatorial algorithm with linear worst-case complexity [4]. In Figure 3(a) the supports of the dual in S of two possible representatives c1 and c2 of the cohomology generators are shown.(2)
Thinned currents are found by premultiplying the generators of S by the incidence matrix Cc between face and edge pairs [4] restricted to Kc.(3)
Finally, a vectorialized version of the ESTT algorithm [4] is run on the whole complex K for all thinned currents at once. The ESTT algorithm is a general version of the Webb–Forghani iterative algorithm [25] to obtain a discrete field whose discrete curl has to match the thinned current vectors obtained at the previous step. Its typical complexity is linear, even though the worst-case complexity is cubical. The output of the ESTT restricted to Ka forms the required lazy cohomology generators.Figure 3
(a) The cycles on the dual complex ofS which represent the supports of the dual of the representatives c1 and c2 of one of the possible cohomology basis H1(S,Z). This particular basis is such that the dual cycle relative to c1 is bounding in the complement of K, whereas the one of c2 is bounding inside K. (b) The support of the dual of the representative of a “mixed” generator c3.
(a)
(b)The output lazy generators span theH1(Ka,Z) cohomology group, but they do not form a base given that they are dependent. In the example, the generator obtained by processing c1 span H1(Ka,Z), whereas the generator obtained with c2 is not useful given that the dual of c2 bounds a discrete surface inside Kc (therefore it is not topologically interesting but it is used anyway in the final system of equations).Disentangling boundary generators is far from being trivial especially because it is in general not enough to pick half of the generators to produce the suitable basis. This is because the generators of the basis may be “mixed,” meaning that one of them is a linear combination of the two asc3 in Figure 3(b). Thus, in general there are no easier ways to do it than the change of basis described in the appendix of [4].The novel idea of lazy cohomology generators [4] is exactly that all boundary generators produce a cut, even if half of them may be eliminated. This approach does not provide a sensible slowdown in the system solver time, whereas it provides great speedup in the topological preprocessing part. Concerning the quality of the solution, the results with a standard versus lazy basis in terms of eddy currents are the same up to solver tolerance.This technique is appealing first of all because the topological preprocessing requires mere seconds even on very complicated problems. Moreover, the DS algorithm is proved to be general; therefore it works for every possible input, no matter how complicated the geometry or the topology of the insulating region is.


### 2.4. Critical Survey on an Alternative Solution

In this section we explain why the technique proposed by He et al. [20] lacks a number of decisive details. Without them, it is impossible to assess its correctness nor to implement and use it. The algorithm presented in [20] consists of two steps, Step 1 and Step 2, summarized below.In Step1, the algorithm of [20] finds acut in the two-dimensional boundary of Kc by searching for a maximal collection of 2-simplices S∈∂Kc such that S does not have any nontrivial cycle. One can prove that by doing so thesum of supports of a first cohomology group basis of ∂Kc is the edges in the difference ∂Kc∖S. However, this technique does not specify how to isolate each cut from their union.Let us now suppose that we are able in some way to separate the each cut. In this case, let us assume that Step2 of [20] algorithm starts from the two canonical cohomology generators of a torus, like c1 and c2 of Figure 3(a). Then, the tetrahedra from Ka are iteratively added to the surface S. One can think about this as process of thickening the surface S in the direction of the insulator.The problem appears when this thickened surface collides with itself. A tetrahedron located in the spot of the collision will not be added to the thickenedS given that it would close a loop. This produces an additionalsprue to the output of Step 2 that has to be eliminated to obtain a correct output; see also [29, Figure  5]. No technique to remove it is discussed in [20]. Let us assume, however, that there is a way to remove all sprues. What we can get in this case is a sum of cut surfaces. Once again it is not discussed how the individual cuts may be obtained from their union. Note that obtaining separated cuts is essential for the electromagnetic computations. Even if there exists a way to separate those cuts, for some surfaces S, one simply cannot get the generators by using this procedure. Let us consider an example of a Sudanese Möbius band cut surface represented in Figure 4. It can be obtained when the complement of the surface S in the boundary of the conductor contains a noncanonical generator of the torus, as c3 in Figure 3(b). In this case, the “cut” one may obtain from it is a Möbius band; it is therefore nonorientable. Since it is not orientable, it cannot be used to fulfill the global Ampère’s law implicitly. In this case, again, the method as presented in [20] fails.Figure 4
Four views form a different angle of the so-called Sudanese Möbius band or Möbius snail (image taken from [27]).


## 2.1. De Rham Cohomology

IfKa is simply connected (or topologically trivial), as in Figure 1(a), a curl-free vector field, like the magnetic field h in an insulator, can be expressed as the gradient of the magnetic scalar potential Ω.Figure 1
The insulating regionKa is obtained as the complement with respect to a box (not represented in the picture for clarity) of (a) a trivial conductor, (b) a solid torus, (c) a twofold solid torus, that is, a conductor with two handles.
(a)
(b)
(c)Topology starts playing a role in theT-Ω formulation when Ka is not simply connected. Let us consider an example where Ka is the complement of a conductive torus Kc with respect to a box as in Figure 1(b). The magnetic field h is curl-free in Ka but its circulation on the loop c of Figure 1(b), because of Ampère’s law, has to match the electric current ic that flows around Kc. The consequence is that the gradient of the magnetic scalar potential is not enough to span the space where the magnetic field h lives.In the following we use concepts of algebraic topology that due to the limited space cannot be reproduced here. Please consult [23] for a formal introduction or [4, 6, 15, 24] for an informal one.In algebraic topology, thefirst de Rham cohomology groupHdR1(Ka) of Ka, is exactly the space of curl-free vector fields in Ka that are not gradients [3]. Therefore, by its very definition, HdR1(Ka) is the space we have to add to the space generated by the well-known hierarchical basis functions described in [25]. It is known that the dimension of this space is equal to the first Betti number β1(Ka) of Ka [3].Focusing again on the example in Figure1(b), β1(Ka)=1, thus the basis for the HdR1(Ka) space is composed by just one element calledgenerator. Let us construct the generator g∈HdR1(Ka) as a curl-free field that has a circulation 1 on a loop linking the conductor as c in Figure 1(b). If one now considers the product icg∈HdR1(Ka), by varying the independent current ic one is able to span the whole HdR1(Ka) space. Therefore, ic is a new degree of freedom that has to be added as an additional unknown of the eddy current problem.This property is generalizable to more complicated topologies; see [4, 26] for a more formal explanation. Here we informally show what happens in a more complicated example, that is, let us consider a twofold solid torus as Kc; see Figure 1(c). Here, β1(Ka)=2 (the conductor has two handles); thus there are two independent currents ic1 and ic2. The HdR1(Ka) space is generated by ic1g1+ic2g2, where the two generators g1,g2 are such that(3)∫cigj·t^dl=δij, where t^ is the tangent vector of the loop ci and δij is the Kronecker delta.In general, the magnetic fieldh is represented as(4)h=∇Ω+∑k=1β1Kaickgk, where the space of gradients of the magnetic scalar potential Ω has been already defined in [25] and each independent current has to be included as an additional degree of freedom to span the whole space of the eddy current problem solution.


## 2.2. Construction ofHdR1(Ka) Generators

Let us try to replicate the properties ofg in the example of Figure 1(b) at the discrete level. Let us call g the discrete counterpart of g. First, the circulation of g on a loop c made of edges that links the conductor like c has to be 1 as for g. Therefore, it is natural to define g as a discrete field having integer coefficients assigned to mesh edges. Second, we want g to be curl-free as its continuous equivalent. This translates in the discrete case to the condition Cag=0, where Ca is the face-edge incidence matrix of the mesh restricted to Ka. That is, the discrete circulation of g over the edges in the boundary of each face of the mesh must be zero. It is possible to prove that such edge discrete fields may be interpreted as the representatives of the generators of the first cohomology group over integers H1(Ka,Z) [4, 26]. Then, to retrieve g, we just need to interpolate g with the standard Whitney edge elements basis functions [1]. Thus, it is evident that there is no need to use higher-order basis functions to span the de Rham space.H
1
(
K
a
,
Z
) are realized using thethick cut technique [28]; see Figure 2.Figure 2
(a) Let us consider a cross section of the example in Figure1(b). (b) The thick black edges belong to the thick cut g (i.e., to the support of the representative of the cohomology generator H1(Ka,Z)). The coefficients assigned to thick edges are such that Cag=0, where Ca is the restriction of C to Ka, and that the circulation on all discrete cycles linking the conductive torus like c is 1.
(a)
(b)The advantage in switching to cohomology over the integers is that the generators ofH1(Ka,Z) can be rigorously and efficiently computed. In this respect, it is important to note that the property of the basis functions of being hierarchical is of fundamental importance: thus, the term of h related to H1(Ka,Z) can be accounted for when using higher order by adding it to the scalar first-order dofs on the edges of Ka. In what follows we give a detailed description of an algorithm that performs this task.


## 2.3. Efficient Computation ofH1(Ka,Z) Generators

We find the required nonlocal basis functions spanning the first de Rham cohomology group space by exploiting the recently introduced idea oflazy cohomology generators [4]. In what follows we recall the DS algorithm that finds the lazy generators of the first cohomology group H1(Ka,Z):(1)
First the discrete surfaceS=Kc∩Ka (see Figure 3(a)) is extracted and its first cohomology group generators H1(S,Z) are computed with a combinatorial algorithm with linear worst-case complexity [4]. In Figure 3(a) the supports of the dual in S of two possible representatives c1 and c2 of the cohomology generators are shown.(2)
Thinned currents are found by premultiplying the generators of S by the incidence matrix Cc between face and edge pairs [4] restricted to Kc.(3)
Finally, a vectorialized version of the ESTT algorithm [4] is run on the whole complex K for all thinned currents at once. The ESTT algorithm is a general version of the Webb–Forghani iterative algorithm [25] to obtain a discrete field whose discrete curl has to match the thinned current vectors obtained at the previous step. Its typical complexity is linear, even though the worst-case complexity is cubical. The output of the ESTT restricted to Ka forms the required lazy cohomology generators.Figure 3
(a) The cycles on the dual complex ofS which represent the supports of the dual of the representatives c1 and c2 of one of the possible cohomology basis H1(S,Z). This particular basis is such that the dual cycle relative to c1 is bounding in the complement of K, whereas the one of c2 is bounding inside K. (b) The support of the dual of the representative of a “mixed” generator c3.
(a)
(b)The output lazy generators span theH1(Ka,Z) cohomology group, but they do not form a base given that they are dependent. In the example, the generator obtained by processing c1 span H1(Ka,Z), whereas the generator obtained with c2 is not useful given that the dual of c2 bounds a discrete surface inside Kc (therefore it is not topologically interesting but it is used anyway in the final system of equations).Disentangling boundary generators is far from being trivial especially because it is in general not enough to pick half of the generators to produce the suitable basis. This is because the generators of the basis may be “mixed,” meaning that one of them is a linear combination of the two asc3 in Figure 3(b). Thus, in general there are no easier ways to do it than the change of basis described in the appendix of [4].The novel idea of lazy cohomology generators [4] is exactly that all boundary generators produce a cut, even if half of them may be eliminated. This approach does not provide a sensible slowdown in the system solver time, whereas it provides great speedup in the topological preprocessing part. Concerning the quality of the solution, the results with a standard versus lazy basis in terms of eddy currents are the same up to solver tolerance.This technique is appealing first of all because the topological preprocessing requires mere seconds even on very complicated problems. Moreover, the DS algorithm is proved to be general; therefore it works for every possible input, no matter how complicated the geometry or the topology of the insulating region is.


## 2.4. Critical Survey on an Alternative Solution

In this section we explain why the technique proposed by He et al. [20] lacks a number of decisive details. Without them, it is impossible to assess its correctness nor to implement and use it. The algorithm presented in [20] consists of two steps, Step 1 and Step 2, summarized below.In Step1, the algorithm of [20] finds acut in the two-dimensional boundary of Kc by searching for a maximal collection of 2-simplices S∈∂Kc such that S does not have any nontrivial cycle. One can prove that by doing so thesum of supports of a first cohomology group basis of ∂Kc is the edges in the difference ∂Kc∖S. However, this technique does not specify how to isolate each cut from their union.Let us now suppose that we are able in some way to separate the each cut. In this case, let us assume that Step2 of [20] algorithm starts from the two canonical cohomology generators of a torus, like c1 and c2 of Figure 3(a). Then, the tetrahedra from Ka are iteratively added to the surface S. One can think about this as process of thickening the surface S in the direction of the insulator.The problem appears when this thickened surface collides with itself. A tetrahedron located in the spot of the collision will not be added to the thickenedS given that it would close a loop. This produces an additionalsprue to the output of Step 2 that has to be eliminated to obtain a correct output; see also [29, Figure  5]. No technique to remove it is discussed in [20]. Let us assume, however, that there is a way to remove all sprues. What we can get in this case is a sum of cut surfaces. Once again it is not discussed how the individual cuts may be obtained from their union. Note that obtaining separated cuts is essential for the electromagnetic computations. Even if there exists a way to separate those cuts, for some surfaces S, one simply cannot get the generators by using this procedure. Let us consider an example of a Sudanese Möbius band cut surface represented in Figure 4. It can be obtained when the complement of the surface S in the boundary of the conductor contains a noncanonical generator of the torus, as c3 in Figure 3(b). In this case, the “cut” one may obtain from it is a Möbius band; it is therefore nonorientable. Since it is not orientable, it cannot be used to fulfill the global Ampère’s law implicitly. In this case, again, the method as presented in [20] fails.Figure 4
Four views form a different angle of the so-called Sudanese Möbius band or Möbius snail (image taken from [27]).


## 3. Numerical Results

The DS algorithm together with theT-Ω formulation with hierarchical basis functions has been implemented inside the EMS solver (https://www.emworks.com). We validated the software and assessed the performances of first and second-order approximation in solving various eddy current problems.As a first benchmark, we use the TEAM Workshop problem 7; see [30]. It consists of a racetrack shaped coil driven by a time-harmonic current (amplitude 2742AT, frequency f=200 Hz) over a square aluminum plate. Figure 5 represents the eddy currents resulting from the solution with first and second order. Table 1 contains the comparison of the computed heat losses P=∫Kρj2, where j is the current density vector field and ρ is the resistivity.Table 1
TEAM problem7 benchmark.
Number oftetrahedra
Losses, 1st order [W]
Losses, 2nd order [W]
52539
7.11
8.02
126615
7.59
8.83
225478
9.02
9.37
423841
9.32
9.60
755333
9.66
-Figure 5
Eddy currents obtained with the electromagnetic code EMS.
(a)
1st order
(b)
2nd orderAs a second benchmark, we consider the TEAM Workshop problem21, see [31]. Figure 6 contains the eddy currents resulting from the solution with first and second order in case of benchmark 21a-3. Table 2 contains the comparison of the computed heat losses for the cases 21a-1, 21a-2, and 21a-3.Table 2
Loss in Watt for the TEAM problem 21 benchmark.
Problem
Tetrahedra
Loss, 1st ord.
Loss, 2nd ord.
Measured
21a-1
45177
3.19
3.40
3.40
21a-2
29577
1.58
1.66
1.68
21a-3
199781
0.95
1.04
1.25Figure 6
Eddy currents obtained with the electromagnetic code EMS ((a) 1st order and (b) 2nd order) on the TEAM problem 21a-3.
(a)
(b)


## 4. Conclusions

Lazy cohomology generators technology has been integrated inside a commercial eddy current solver employing second-order hierarchical basis functions. As expected, the second-order formulation indeed provides more accurate results than its first-order counterpart. The inclusion of automatic mesh adaptivity to render the second order more efficient is left for further studies.



---

*Source: 8308643-2018-01-11.xml*