On torque computation in electric machine simulation by harmonic mortar methods

The use of trigonometric polynomials as Lagrange multipliers in the harmonic mortar method enables an efficient and elegant treatment of relative motion in the stator-rotor coupling of electric machine simulation. Explicit formulas for the torque computation are derived by energetic considerations, and their realization by harmonic mortar finite element and isogeometric analysis discretizations is discussed. Numerical tests are presented to illustrate the theoretical results and demonstrate the potential of harmonic mortar methods for the evaluation of torque ripples.


Introduction
A particular challenge for electric machine simulation is the relative motion of stator and rotor and the computation of quantities of interest depending on the rotation angle, e.g., the torque as a measure for the magneto-mechanic energy conversion. Various approaches have been proposed to tackle the coupling of subdomain problems on moving geometries, e.g., moving band, mortar and Lagrange multiplier or discontinuous Galerkin methods [1,7,13,16], all with their specific advantages and shortcomings. In addition, different strategies have been discussed for the numerical computation of the torque, e.g., via virtual displacements or the Maxwell stress tensor; see [12] for an overview.
In this paper, we consider harmonic mortar finite element and isogeometric analysis methods proposed in [3,10], which are based on finite element or isogeometric analysis approximations of the magnetic problems in the stator and rotor subdomains, coupled by a Lagrange multiplier technique using trigonometric functions. These methods have several interesting properties: (a) They are based on the Galerkin approximation of a certain weak formulation of the problem using Lagrange multipliers. This will allow us to give an unambiguous definition of the torque resulting in an exact energy conservation principle also on the discrete level. (b) The computed torque, as a function of the rotation angle, can be shown to be smooth, i.e., infinitely differentiable, and therefore no spurious torque ripples are introduced by non-smoothness of the numerical approximation.
(c) The use of trigonometric functions for the Lagrange multipliers additionally allows for an efficient assembling of the coupling matrices for multiple rotation angles, avoiding the repeated integration of interface terms. The resulting methods therefore are particularly well-suited for multiquery applications, e.g., the construction of performance maps, the setup of reduced order models, coupled frequency or time domain simulations. This will be illustrated by computation of cogging torque in the numerical tests.
The remainder of the manuscript is organized as follows: Sect. 2 introduces the magnetostatic setting for electric machines in two dimensions. The strong and weak formulations of the problem are presented and the energy balance is used to derive the torque. The discretization of the problem is given in Sect. 3 using the Galerkin approach. Section 4 discusses the harmonic mortar method and its implementation. Numerical results for the torque computation using the harmonic mortar method are demonstrated in Sect. 5. Finally the paper is concluded with a discussion in Sect. 6.

Problem statement
For ease of presentation, we consider a linear magnetostatic setting in two space dimensions. The geometric setup represents the cross section of a cylindrical device and consists of a stator domain S completely enclosing the rotor domain R with common interface = ∂ S ∩ ∂ R = ∂ R ; see Fig. 1. Motivated by the application context, we assume that is a circle placed within the air gap and centered around the origin. Using a standard vector potential formulation [1,10], the magnetic fields in the two subdomains are described by where a S , a R are the z-component of the magnetic vector potential, j e is the impressed electric source density, and M ⊥ = (m y , -m x ) the rotated magnetization vector. The equations are formulated in two different coordinate systems attached to the respective subdomains, such that the position of material continuities stays fixed under relative motion. We further assume that a S = 0 on (3) Figure 1 Sketch of the geometry. x θ denotes a point in the fixed coordinate system (red) of the stator domain S , whereasx θ denotes a point in the rotated coordinate system (blue) of the rotor domain R . They are related by x θ =x θ -α where α is the angle of rotation on the outer boundary = ∂ S \ . The continuity of the magnetic field and the magnetic vector potential across the interface are encoded in where ρ α describes the rotation of a point by angle α around the origin and (f • g)(x) = f (g(x)) is the concatenation of functions. Furthermore, n denotes the unit normal vector on pointing from R to S . The second coupling condition allows us to introduce which amounts to the tangential trace of the magnetic field at the interface.
Dependence on α By (4)-(5) the solutions (a S , a R , λ) of (1)-(6) implicitly depend on the rotation angle α, even if the sources j e and M are independent of α, which we assume in the sequel. We write u(·|α) to highlight the dependence of a function u on the angle α and denote by the derivative of such a parameter-dependent function u with respect to the parameter α. By differentiation of (4) with respect to α, we can see that where x θ-α = ρ -α (x θ ) is the point in the rotated coordinate system of ∂ R corresponding to a point x θ = (r cos θ , r sin θ ) on the fixed interface ; see Fig. 1. Furthermore, d dθ u(x θ ) = r n ⊥ · ∇u(x θ ) amounts to the tangential derivative of u along and r is the radius of .
Weak formulation Using standard arguments, sufficiently regular solutions of (1)- (5) can be seen to satisfy the variational identities for all smooth test functions v S , v R , and μ defined on S , R , and , respectively, with v S = 0 on . Using the considerations outlined above, the coupling condition (11) can also be rephrased as which will be used in the following considerations.
Energy balance and torque computation As a next step, we introduce the magnetic energy of the system. For ease of notation, we assume a unit length L of the device in axial direction. Then Let us note that, via the fields a S , a R , the energy also implicitly depends on the rotation angle α. By elementary calculations, we may then compute where we used the variational identities (9)- (10) with v S = a S and v R = a R in the second step. Using the differential form (12) of the coupling condition, we arrive at the energy balance where the first term denotes the electric work required to maintain the electric current j e in the stator coils (neglecting Ohmic losses), and the second term is the mechanic work required for an infinitesimal rotation of the system. These energy-based considerations immediately give rise to the following definition of the torque Note that a R and λ and hence T are functions of the rotation angle α. In general, the energy and torque have to be scaled by the axial length L.
Remark 1 The above formula can be derived in various ways, see [4,9,12]. Since the stator and rotor were both assumed to be rigid bodies, the torque can be computed here directly without resorting to the Maxwell stress tensor. Using rotational symmetry as well as integration by parts along the interface , we can express the torque alternatively as which may be more convenient depending on the particular setting.

Galerkin approximation
For discretization of the problem (1)-(5), we consider Galerkin approximations of the variational equations (9)-(11) using appropriate finite-dimensional sub-spaces The discrete problem for a fixed angle α then reads as follows.
By choosing bases for the approximation spaces V h,S , V h,R and M N , the discrete variational problem can be turned into a linear system of equations for the corresponding coordinate vectors; see Sect. 4 for details. The above variational form is however more convenient for the analysis of the discretization.
Well-posedness For the following considerations, we assume that V h,S , V h,R , and M N are finite-dimensional and that ν S , ν R are uniformly positive and bounded from above. The source currents j e and magnetization vector M, as well as the domains S , R are assumed to be sufficiently regular. We can then establish the well-posedness of the discretization scheme as follows.

Lemma 1 Assume that
Then the discrete variational problem (16)-(18) is uniquely solvable and The assertion is a direct consequence of Brezzi's saddle-point theory [6]; also see [7,11] for application in the current context. Under the stability condition (19), the Galerkin approximations obtained by Problem 1 thus always lead to quasi-optimal error estimates with respect to the chosen approximation spaces. Two particular examples of Galerkin approximations will be discussed below.
Discrete torque computation and energy balance Mimicking the formula (14) for the torque obtained on the continuous level, it seems natural to define the discretized approximation of the torque as The torque computation for a given angle α thus only requires the solution of one magnetostatic interface problem (16)-(18) for this specific angle. With the same reasoning as on the continuous level, we obtain the following result.
with discrete torque T h,N defined by the formula (20).
This assertion illustrates that the electro-magneto-mechanic energy balance holds exactly also on the discrete level.
Alternative representation and smooth dependence on α If the space M N of discrete Lagrange multipliers consists of sufficiently smooth functions, we may express the discrete torque equivalently as By a recursive argument, one can then show the following result.

Lemma 3
Let M N consist of k -1 times continuously differentiable functions. Then the solution (a S , a R , λ N ) is k times differentiable with respect to α.
A proof of this assertion follows from the algebraic structure of the discretized problem; see the end of the next section.

Harmonic mortar discretizations
As a particular choice of a Galerkin approximation, we now consider the following choice of spaces: V h,S and V h,R are standard H 1 -conforming finite element spaces [5] over appropriate triangulations of the stator and rotor domain S and R . Note that no conformity of the mesh across the interface is required. For approximation of the Lagrange multipliers, we consider the space of trigonometric polynomials of degree ≤ N .
Remark 2 As shown in [11], the inf-sup condition (19) holds true, if N is chosen not too large. As a rule of thumb, N must not exceed the number of degrees of freedom in V h,S and V h,R associated to the interface . For a finite element approximation on a quasi-uniform mesh and with fixed polynomial degree p, stability can thus be guaranteed for N h -1 p; see [11,Theorem 1]. Similar results hold for isogeometric methods for the stator and rotor domains [3,8]. Following the considerations of the previous section, both methods allow for a torque computation consistent with a discrete energy conservation principle.
Implementation of the harmonic mortar finite element method In the following, we discuss some details concerning the implementation of the harmonic mortar finite element method outlined above. As we will see, the use of trigonometric polynomials for the discrete Lagrange multipliers brings advantages also for the numerical realization.
After choosing appropriate basis functions for the finite element subspaces V h,S , V h,N and the trigonometric polynomials in M N , the discrete variational problem (16)-(18) can be turned into a linear system The source vectors J e and J M and the matrices K S , K R and B S are independent of the angle α, whereas B R and hence also the solution vectors a S , a R and λ depend on the rotation angle. For the harmonic mortar, the entries of the coupling matrix B R (α) can be computed by the following expressions where φ j are the basis functions for the Galerkin approximation on R .
Algebraic energy balance and torque computation On the algebraic level, the energy is Hence the discrete torque can be simply computed by the algebraic formula Note that B R (α) can be expressed as B R (α) = R (α)B R (0) by using (27). Moreover, the matrix R (α) is block diagonal with 2 × 2 blocks whose entries can be computed analytically; see above.
Derivatives with respect to α and proof of Lemma 3 By formal differentiation of the above equations, we can see that Together with (28)-(30), one can see that the solution of (24)-(26) is differentiable with respect to the angle α. By recursion, we may further conclude that the torque T h,N of the harmonic mortar method is a smooth function of the angle α. The same argument provides a proof of Lemma 3.

Numerical results
In the following, we demonstrate the application of the harmonic mortar finite element and isogeometric analysis method for computing the torque in a six-pole permanent magnet synchronous machine shown in Fig. 2. The geometry consists of a rotor with three pole pairs with embedded permanent magnets and a stator with 36 slots filled with copper windings. The length of the actual machine is L = 0.1 m, so the energy and torque have to be scaled by L as outlined above. Some information about the geometry and material parameters is summarized in Table 1; for details see [2, Chapter V.A].
Simulation setup In our numerical tests, we compute the cogging torque via formula (22) in the absence of excitation currents j e = 0. For the discretization of the magnetostatic problems in the stator and rotor domains, we used open source package GeoPDEs   Table 2 Relevant and irrelevant Fourier modes of the Lagrange multipliers of an IGA simulation result for α = 7°, C n = c 2 n + d 2 n is the total amplitude n / ∈N C n 9.2621 × 10 3 1.6651 × 10 3 5.5848 × 10 3 4.9720 × 10 -6 [17]. The computations are carried out with B-splines of degree p = 1 and p = 2 with maximal continuity as basis functions for the solution; this corresponds to a standard finite element method (FEM) and a multipatch spline discretization with patchwise C 1continuous splines of second order (IGA). To accommodate for material discontinuities, only C 0 -continuity is enforced across patch borders; see [3] for details. At the external boundaries of the stator and rotor, the magnetic vector potential a is set to zero. Symmetry of the solution A snapshot of the magnetic vector potential and its isolines, which correspond to the magnetic flux lines, is depicted in Fig. 2. Note that by symmetry of the rotor and stator geometry, the solution components (a S , a R , λ) have a 60 • antisymmetry. In particular, for the Fourier expansion of the computed Lagrange multiplier only the coefficients c n , d n with index n ∈ N = {3, 9, 15, . . .} are non-zero. By careful construction of the meshes, the geometric symmetries are preserved on the discrete level and the same behavior is observed for the numerical solution, as can be seen from Table 2.
In practice, the computational cost can be reduced substantially by solving the problem only for one pole of the machine and using anti-periodic boundary conditions; see [3] for details. In the present case the number of degrees of freedom can be reduced by a factor of six. Using this prior knowledge, most of the Lagrange multipliers can be eliminated already during assembling leading to more efficient computations, in particular in connection with the Schur complement technique outlined in the Appendix.  Torque computation In Fig. 3, we display the torque computed via formula (14) with FEM (p = 1) and IGA (p = 2) for different orders N of harmonic basis functions. By the 10 • symmetry of the stator slots and the corresponding symmetry of the mesh, the computed torque can be shown to exhibit a periodic behavior with a period length of 10 • . Moreover, the torque can be seen to be anti-symmetric with respect to half of the period. We further observe that the torque converges towards that of the reference solution when increasing the number N of harmonic basis functions. As indicated in Remark 2, the computation becomes unstable when the dimension of the Lagrange multiplier space M N exceeds the number of degrees of freedom of V h,S and V h,R associated to the interface . Due to periodicity, the torque T h,N can be approximated by a Fourier series  Table 3.

Discussion
We have discussed the harmonic mortar method using trigonometric polynomials as Lagrange multipliers for the efficient treatment of rotation in electric machine computation. The energy balance has been used to derive an explicit formulation for the computation of the torque in the continuous and discretized settings. It has been shown that the electromagneto-mechanic energy balance holds exactly also on the discrete level. Furthermore, if the functions in the Lagrange multiplier space are smooth then this is inherited by the torque. The derived formulation has been applied to a typical permanent magnet synchronous machine model using isogeometric analysis and a lowest order finite element method with an exact geometry mapping. The cogging torque has been computed with the proposed method and the influence of the number of Lagrange multiplier basis functions on the solution has been investigated.