Multi-objective free-form shape optimization of a synchronous reluctance machine

Abstract

Purpose

This paper aims to deal with the design optimization of a synchronous reluctance machine to be used in an X-ray tube, where the goal is to maximize the torque while keeping low the amount of material used, by means of gradient-based free-form shape optimization.

Design/methodology/approach

The presented approach is based on the mathematical concept of shape derivatives and allows to obtain new motor designs without the need to introduce a geometric parametrization. This paper presents an extension of a standard gradient-based free-form shape optimization algorithm to the case of multiple objective functions by determining updates, which represent a descent of all involved criteria. Moreover, this paper illustrates a way to obtain an approximate Pareto front.

Findings

The presented method allows to obtain optimal designs of arbitrary, non-parametric shape with very low computational cost. This paper validates the results by comparing them to a parametric geometry optimization in JMAG by means of a stochastic optimization algorithm. While the obtained designs are of similar shape, the computational time used by the gradient-based algorithm is in the order of minutes, compared to several hours taken by the stochastic optimization algorithm.

Originality/value

This paper applies the presented gradient-based multi-objective optimization algorithm in the context of free-form shape optimization using the mathematical concept of shape derivatives. The authors obtain a set of Pareto-optimal designs, each of which is a shape that is not represented by a fixed set of parameters. To the best of the authors’ knowledge, this approach to multi-objective free-form shape optimization is novel in the context of electric machines.

Citation

Gangl, P., Köthe, S., Mellak, C., Cesarano, A. and Mütze, A. (2022), "Multi-objective free-form shape optimization of a synchronous reluctance machine", COMPEL - The international journal for computation and mathematics in electrical and electronic engineering, Vol. 41 No. 5, pp. 1849-1864. https://doi.org/10.1108/COMPEL-02-2021-0063

1. Introduction

In many industrial applications, the design of electric machines has to be tailored to the application at hand because off-the-shelf solutions are not available. The design of electric machines is usually based on engineering knowledge and is sometimes refined by geometric optimization. The most widely used approach is to introduce geometric parameters and optimize these, either using stochastic optimization algorithms or derivative-based methods, see Bramerdorfer et al. (2018) for an overview article. While derivative-based optimization algorithms successively improve a given initial geometry by means of gradient information and are known to converge to a local optimum rather fast, stochastic algorithms include random effects and are less prone to getting stuck in local optima. In practice, one is usually confronted with several conflicting objective functions, thus making multi-objective optimization capabilities for finding a Pareto optimal set of designs important. The extension to a multi-objective setting is more straightforward in the case of many stochastic optimization algorithms; however, it can also be achieved in the case of derivative-based methods (Doganay et al., 2019).

For a thorough overview over the field of multi-objective optimization, including both evolutionary and derivative-based algorithms, we refer the reader to the monographs (Deb, 2001; Miettinnen., 1998). The most widely used approaches of gradient-based multi-objective optimization consist in a way to scalarize vector-valued optimization problems. For instance, the weighted average method minimizes convex combinations of the given cost functions for a range of different convex parameters, or the ϵ-constraint method chooses one cost function to be minimized and poses constraints on all other cost functions. While the former method has the shortcoming that non-convex regions of the Pareto front cannot be recovered, this is possible using the latter approach. For an overview over other, more sophisticated aggregation methods, we refer the reader to Deb (2001). Besides the mentioned scalarization approaches, there also exist methods that work directly with the vector-valued cost function and seek a descent direction that is common to all cost functions. We mention the approach introduced in Désidéri (2012) and the method introduced in Fliege and Svaiter (2000). Both approaches have been extended to the case of parametric shape optimization, see Giacomini et al. (2014) and Doganay et al. (2019) and Bolten et al. (2021), respectively. These methods are able to recover solutions lying on non-convex regions of the Pareto front. A further important criterion for a successful construction of a Pareto front is a method’s ability to create evenly spaced points on the Pareto front. This is usually a difficult issue that can be addressed, e.g. by adaptive choices of weights (Eichfelder, 2009) or a predictor-corrector scheme (Schmidt and Schulz, 2008).

In recent years, non-parametric shape optimization methods based on the mathematical concept of shape derivatives (Delfour and Zolésio, 2011) (often referred to as free-form shape optimization approaches) have become a more and more popular tool for the design optimization of electric machines, see (Gangl et al., 2015; Kuci, 2018; Putek et al., 2016) for approaches using the finite element method or the recent work (Merkel et al., 2021) in the context of isogeometric analysis. In these approaches, the geometry is not parametrized by a finite number of scalar values, but the design variable is a set, e.g. the set of points occupied by ferromagnetic material in the rotor of an electric machine. Starting out from a given initial design, the design is updated by the action of a smooth vector field, thus allowing for any kind of design that is topologically equivalent to the initial design. This way, often new and innovative designs can be obtained.

The purpose of this paper is twofold: on the one hand, we extend the gradient-based multi-objective optimization method introduced in the case of a parametrized geometry in Doganay et al. (2019) to the case of free-form shape optimization. This is done by replacing the notion of gradients of functions defined on ℝn by the mathematical concept of shape derivatives of shape functions defined on a set of shapes. This allows to exploit the flexibility of free-form shape optimization methods, as well as their fast convergence properties also in the practically important case of multiple competing objective functions. On the other hand, we use this method on both, the more standard single-objective case and the case of two objective functions, to find (Pareto-)optimal designs of a synchronous reluctance machine. Comparing our results with the results obtained by a stochastic parameter optimization confirms the higher degree of flexibility and computational efficiency of our approach compared to parametric design optimization.

The rest of this paper is organized as follows: in Section 2, we introduce the problem at hand and state the mathematical model. We recall the main ingredients for a free-form shape optimization method and apply the algorithm to our problem in Section 3. In Section 4, we show an extension of the gradient-based free-form shape optimization algorithm to the case of multiple objective functions before concluding in Section 5.

2. Problem description

2.1 Physical model

We consider the design optimization of a synchronous reluctance machine (SynRM), i.e. a motor that is based solely on the reluctance principle. This motor generates torque exclusively by a difference of reluctance between two axes, namely, the d-axis and the q-axis (the location of the axes is defined by the number of poles of the machine). Thus, torque generation is not based on any transient behavior or quantity and a static magnetic field analysis is sufficient. The machine under investigation is intended for the use in an X-ray tube for medical applications. The considered rotor will be operated in a vacuum and therefore must be built of solid pieces of metal (as opposed to the commonly used steel sheet structure). Additionally, the air gap of the motor is unusually large (e.g. 10 mm with an outer stator diameter of 130 mm) decreasing the torque capability of the machine. Furthermore, the rotor has to withstand temperatures of up to 450^∘C (Mellak et al., 2018).

The synchronous reluctance machine is particularly suitable for such an application mainly because of its ruggedness and construction simplicity and the absence of rotor windings (Xu et al., 1991). As per the operation mode of the machine, quick acceleration and subsequent braking of a tungsten disk are required. Typically, this sequence takes at maximum 10 s. Figure 1 shows the machine under investigation. The stator is a three-phase stator with one pole pair, the rotor consists of alternating magnetically conducting (blue) and non-conducting layers (gray). The reference design parameters of the machine are stated in Table 1.

Figure 2 shows the simplified vector diagram of a synchronous reluctance machine. The d-axis of the machine is the path with least reluctance, the q-axis is the path with the highest reluctance. In the d-q axis theory, the torque is expressed as:

T=3Np2(λdIq−λqId),

where N_p denotes the number of pole pairs, λ_d and λ_q are the magnetic flux linkages and I_d and I_q are the currents in d-axis and q-axis direction, respectively. Alternatively, using the inductances L_d and L_q as well as the stator current I_s and current angle β, the torque is expressed as:

(1) T=3Np4(Ld−Lq)Is2sin(2β) .

Evidently, as per equation (1), assuming linear lossless behavior and a fixed stator current I_s, the maximum torque can be achieved with a machine current angle β (angle between current vector and d-axis of the machine, Figure 2) of 45^∘ (Spargo, 2013).

2.2 Optimization goal

A static analysis is chosen to calculate the reluctance torque. Therefore, a current is impressed on the windings according to Table 2. Subsequently, the rotor is rotated and fixed clockwise to create the optimal current angle β of 45^∘. The objective is to increase the torque with the given stator at a constant current and air gap length at the optimum current angle β. The number of conducting and non-conducting layers remains unchanged. Solely the shape of each individual layer is subject to the optimization as to increase the d-axis inductance L_d while, ideally, decreasing the q-axis inductance L_q at the same time.

2.3 Mathematical model

We consider a two-dimensional (2D) cross-section of the machine in the setting of 2D magnetostatics, i.e. B = curlA, where the magnetic vector potential is of the form A = (0,0,u(x₁, x₂))^T. Let D⊂ℝ2 denote the computational domain, which comprises the 2D cross section of the machine as well as a surrounding air region, and let Ω ⊂ D denote the ferromagnetic parts of the machine. The mathematical design optimization problem reads:

(2) maxΩ∈A T(u)

(3) s.t.  −div(νΩ(x,|∇u|)∇u)=Ji, x∈D,u=0,  x∈∂D,

where T represents the torque for the considered rotor position, A is a set of admissible shapes, J_i represents the impressed current density and the magnetic reluctivity is defined piecewise as:

νΩ(x,s)={ν^(s)x∈Ω,ν0x∈D∖Ω¯.

Here, ν^ is a non-linear function, which represents the magnetic reluctivity of the ferromagnetic material, and ν₀ corresponds to the magnetic reluctivity of air. The partial differential equation (PDE) constraint (3) admits a unique solution under natural assumptions on the non-linear function ν^ (Pechstein and Jüttler, 2006). Note that the torque T depends on the shape Ω of the ferromagnetic components via the solution to the PDE constraint (3). Denoting the unique solution to equation (3) for given Ω∈A by u_Ω, we define the reduced cost function T(Ω):=T(uΩ).

3. Free-form shape optimization

We propose a free-form shape optimization algorithm based on the mathematical concept of shape derivatives, which is capable of improving the shape of a given initial geometry without the need of defining geometric parameters. We will outline the main ingredients to the method in the following. We introduce the theory for a general cost function J and will choose J:=−T later in Section 3.3.

3.1 Shape derivative

The shape derivative of a general shape function J=J(Ω) represents the sensitivity of J when the domain Ω is perturbed by the action of a given vector field V. Given a smooth vector field V, which is defined on D, let Ω_t = (id+tV)(Ω) denote the perturbed domain for t > 0, see Figure 3 for an abstract illustration. The shape derivative of J in the direction given by V is defined as:

(4) dJ(Ω;V):=limt↘0 J(Ωt)−J(Ω)t,

provided that this limit exists and the mapping dJ(Ω;V):=limt↘0 J(Ωt) − J(Ω)t, is linear and continuous (Delfour and Zolésio, 2011).

The shape derivative for problem (2)–(3) in the direction given by a vector field V can be derived in an analogous way as it was done in Gangl et al. (2015).

We briefly sketch the main steps to be taken and refer the reader to Gangl et al. (2015) and Gangl et al. (2021) for more details. The shape derivative can be computed as the derivative of a perturbed Lagrangian with respect to the shape perturbation parameter t:

(5) dJ(Ω;V)=∂tG(0,u,p)

where G:[0,τ)×E×F→ℝ is defined as:

G(t,φ,ψ):=L(Ωt,φ∘Φt−1,ψ∘Φt−1)

and, for our problem (2)–(3),

L(Ωt,φ,ψ):=−T(φ)+∫DνΩt(x,|∇φ|)∇φ·∇ψ dx−∫DJiψ dx.

Here, τ > 0 is a small parameter, and E and F denote the function spaces for state and adjoint variable, respectively. Moreover, Φ_t(x) = x + tV (x) denotes the shape perturbation given by a vector field V, and u and p are the state and adjoint state.

For a vector field V that is only supported on the rotor, the shape derivative of problem (2)–(3) reads (Gangl et al., 2015):

(6) dJ(Ω;V)=∫DνΩ(x,|∇u|)((div V)I−∂VT−∂V)∇u·∇p dx−∫D∂sνΩ(x,|∇u|)|∇u|(∂VT∇u·∇u)(∇u·∇p) dx.

Here, p denotes the solution to the adjoint equation, which for the case of the maximization of the torque reads in its strong form:

(7) −div(AΩ(u)∇p)=∂T∂u, x∈D,p=0,     x∈∂D.

with

AΩ(u):=νΩ(x,|∇u|)I+∂sνΩ(x,|∇u|)|∇u|∇u⊗∇u.

3.2 Descent direction

Given a closed formula for the shape derivative, a descent vector field V can be obtained by solving an auxiliary boundary value problem as follows. Let X be a Hilbert space and b:X×X→ℝ a symmetric and positive definite bilinear form. Then, the solution W ∈ X to the variational problem:

(8) b(W,V)=−dJ(Ω;V) ∀V∈X

is a descent direction since it satisfies by construction:

dJ(Ω;W)=−b(W,W)<0.

Thus, it follows from the definition in equation (4) that perturbing Ω a small distance into the direction W will yield a decrease of the cost function J.

The user has some degrees of freedom in the choice of the bilinear form b(·,·) as well as the space X. Common choices include X=H1(D,ℝ2) and b(W,V)=∫D∂W:∂V+W·Vdx or b(W,V)=∫DCϵ(W):ϵ(V)+W·Vdx where ϵ(V)=12(∂V+∂V⊤) and C is a fourth-order elasticity tensor. The latter choice is known to preserve mesh quality better compared to other choices of b(·,·) (Iglesias et al., 2018). An alternative strategy for extracting a descent direction, which also allows for the extension to multiple objective functions, will be discussed in Section 4.2.

3.3 Numerical results

The procedure outlined in Sections 3.1 and 3.2 constitutes the following free-form shape optimization algorithm for minimization of shape function J=J(Ω):

Algorithm 1. Given initial design Ω0, cost function 

J, tolerance tol, k = 0.1) Solve state equation (3) and adjoint equation (7)2) Compute shape derivative 

dJ(Ωk;V) given in (6)3) Compute shape gradient W as solution to (8)4) If ‖W‖ < tol then stop       else set ‖W‖ < tol where 

t=max{1,12,14,18,…}        such that 

J(Ωk+1)<J(Ωk).5) k ← k + 1 and go back to 1)

Note that the shape gradient W computed in Step 3 is a vector field that is defined on the full rotor domain. In Step 4, the new geometry is obtained from the old geometry Ω_k+1 = (id + tW)(Ω_k) = {x + tW(x): x ∈ Ω_k} by just moving every mesh node a distance t into the direction given by W(x). Here, the parameter t is chosen by a line search to guarantee a descent of the cost function J.

We applied Algorithm 1 to problems (2)–(3), i.e. we chose to minimize J(Ω):=−T(Ω), using the finite element software package NGSolve (Schöberl, 2014). In particular, we used the automated shape differentiation capabilities provided by NGSolve, which enables the automated computation of the shape derivative dJ(Ω;V) for a large class of PDE-constrained shape optimization problems (Gangl et al., 2021).

For the space X in (8), we chose the space of all vector-valued H¹ functions defined on the rotor of the machine whose normal component vanishes on the top and bottom boundary parts of the rotor and which vanish at the left and right boundary parts. For the bilinear form b(·,·), we chose the H¹ inner product:

b(V,W):=∫Drot∂V:∂W+V·W dx,

where D_rot denotes the union of the five iron and four air layers, as depicted in Figure 4. The initial design along with the descent vector field W in the first iteration and the results obtained after 70 iterations of Algorithm 1 are depicted in Figure 4. The torque was increased by about 26% from 1.007 to 1.270 Nm. The computational time to obtain the optimized design was about 10 min on a single core.

3.4 Validation

In this section, we validate the results obtained in our numerical experiments by comparing them to an optimization run in JMAG (JSOL-Corporation, 2022). Motivated by the results of the gradient-based optimization, see Figure 4, we parametrized our rotor geometry by means of 13 geometric parameters under symmetry conditions, see Figure 5, and ran a genetic algorithm, which is built into JMAG to maximize the torque. We started with a population size of 300 and ran the algorithm for 50 generations, allowing for 60 children in each generation. The JMAG model used a mesh consisting of 25,109 nodes and 43,776 elements and polynomial degree 1 (compared to 14,063 nodes, 27,310 elements and degree 1 for the NGSolve model). The computational time for evaluating one design in JMAG is about 8 s (compared to 4 s for the coarser NGSolve model). The total computation time of the genetic algorithm consists of the evaluation of the single designs as well as some additional time spent on reading/writing data and on determining the next generation of designs.

The computational time used by the genetic algorithm was about 19 h and a total of 300 + 50 (60 + 13 + 1) = 4,000 designs were examined. The four designs with the highest torque values are depicted in Figure 6. It can be seen that the best designs are similar to the design we obtained by the gradient-based algorithm (Figure 4), but also that the torque values were not quite reached. While one might be tempted to explain such a discrepancy by the fact that different simulation tools were used, we mention that the calculated torques in the two simulation softwares (NGSolve and JMAG) showed a good match for the initial geometry. Thus, it seems like the design in Figure 4 is superior to those obtained by the genetic algorithm in JMAG because more general geometries can be obtained. Of course, the computation time of 19 h could be reduced by reducing the parameters of the genetic algorithm; however, the general order of magnitude remains. Finally note that, as the choice of the geometric parameters was inspired by Figure 4, the designs in Figure 6 would have been unlikely to be found without the knowledge provided by the free-form shape optimization algorithm.

4. Multi-objective shape optimization

In this section, we will consider an extension of the gradient-based free-form shape optimization method presented in Section 3 to the setting of multiple objective functions. We will show how to compute a descent vector field W that assures a descent with respect to several objective functions and will use this approach to obtain an approximation of the Pareto front. We will apply the method to the bi-objective free-form shape optimization problem:

minΩ(J1(Ω)J2(Ω))

where J1(Ω):=−T(Ω) corresponds to the negative of the torque related to Ω and J2(Ω):=Vol(Ω) denotes the volume of the ferromagnetic subdomains of the machine.

4.1 The steepest descent methods for multicriteria optimization in ℝn

We briefly recall the method introduced in Fliege and Svaiter (2000) for optimization problems in ℝn before extending it to the case of multi-objective (free-form) shape optimization. We explain the method for the case of two functions f1,f2:ℝ2→ℝ; however, a general number of objective functions n ≥ 1 is possible.

The main idea of the method is the following: one chooses an initial guess (x₀, y₀) and in each iteration k computes a bi-descent direction, i.e. a vector W(k)∈ℝ2 such that ∇f₁(x_k, y_k)·W^(k) < 0 and ∇f₂(x_k, y_k)·W^(k) < 0. Then, one chooses a small enough step size t and sets (x_k+1, y_k+1) = (x_k, y_k) + tW^(k). By construction, when t is chosen small enough, one obtains a descent of both objective functions, i.e. f_i(x_k+1, y_k+1) < f_i(x_k, y_k) for i = 1, 2.

According to Fliege and Svaiter (2000), the bi-descent direction W = (W₁, W₂) in iteration k can be obtained by solving the following auxiliary optimization problem:

(9) minρ,W1,W2ρ+12(W12+W22)subject to   ∇f1(xk,yk)·W≤ρ,∇f2(xk,yk)·W≤ρ.

As (ρ, W₁, W₂) = (0,0,0) is a feasible point of (9), the solution of (9) will satisfy ρ ≤ 0 and thus ∇f_i(x_k, y_k)·W ≤ 0 for i = 1, 2. The second term in the cost function of (9) is meant to keep the norm of W bounded.

4.2 Multi-objective descent direction for shape derivatives

Given two shape functions J1, J2 and their corresponding shape derivatives dJi(Ω;V), i = 1, 2, we want to find a vector field W such that:

dJ1(Ω;W)<0 and dJ2(Ω;W)<0.

We extend the idea explained in Section 4.1, which has already been used in the framework of parametric shape optimization in (Doganay et al., 2019), to the setting of free-form shape optimization. For that purpose, we consider a finite element discretization using piecewise linear and globally continuous finite elements on a triangular mesh. Denoting the corresponding hat basis functions by φ₁,….,φ_n where n is the number of mesh points and Φ_i = (φ_i,0)^T, and Φ_n+i = (0,φ_i)^T, i = 1,…,n, we have that:

{Φ1,…,Φ2n}

is a basis for the set of all 2D vector fields on the mesh. Thus, after discretization, each vector field W_h can be written as Wh=∑i=12nWiΦi with the coefficient vector W¯:=(W1,…W2n)⊤. Note that we can identify the finite element function W_h with its coefficient vector W¯. To obtain a discrete bi-descent direction W_h, we solve the auxiliary optimization problem to find (ρ,W¯)∈ℝ×ℝ2n:

(10) minρ,W¯ρ+12∑i=12nWi2,s.t.  dJ1(Ω;Wh)≤ρ,dJ2(Ω;Wh)≤ρ.

Due to the linearity of the shape derivatives dJi(Ω;Wh) with respect to W_h, the solution (ρ,W¯)=(0,0)∈ℝ×ℝ2n is a feasible point of (10). Therefore, it follows that the solution (ρ,W¯) to (10) satisfies dJi(Ω;Wh)≤ρ≤0, i = 1, 2, thus giving a bi-descent direction W_h whenever the optimal ρ is negative. Again, the second term in the cost function of (10) is meant to keep the norm of W_h bounded. Of course, an extension of this approach to account for more than two cost functions J1,…,JN is straightforward.

4.3 Obtaining a Pareto front

Proceeding as described in Section 4.2 allows to obtain a bi-descent direction W_h. Thus, starting out from an initial design, iteratively computing a bi-descent vector field and moving the interface a small distance in the direction given by this vector field constitutes a gradient-based free-form shape optimization algorithm for two cost functions. When no further decrease can be obtained, a Pareto optimal point is found.

To obtain many Pareto optimal points, one could start with many different initial designs. However, for the case of shape optimization, it turns out to be more convenient to proceed as follows: consider different scalings of the two objective functions, i.e. apply the gradient-based bi-objective descent algorithm for the two objective functions J1 and wJ2 with different values of the weight w (Doganay et al., 2019). Each choice of the weight w corresponds to a run of the bi-objective descent algorithm and will yield a point on the Pareto front.

4.4 Discussion of multi-objective optimization method

We stress that, in contrast to the widely used weighted sum method, the approach discussed above is also feasible for finding solutions on non-convex regions of a Pareto front, as also noted in Doganay et al. (2019). To illustrate this feature of the method, we apply it to the simple numerical example involving the two cost functions f1,f2:ℝ2→ℝ:

(11) f1(x1,x2)=x1f2(x1,x2)=1+x22−x1−0.1   sin(3πx1)

to be minimized. For this problem, the Pareto front is known to be non-convex [3, Section 3.1.4]. Figure 7 shows the Pareto front of the problem as well as 18 points on it which were obtained by means of the bi-descent algorithm of Section 4.1. Here, to obtain different points on the Pareto front, six different starting values were chosen. In addition, three different weights of f₂ were used as discussed in Section 4.3.

While problem (11) is only academic and is not a shape optimization problem, it illustrates that the procedure that is based on finding descent directions for both cost functions can very well find designs on non-convex parts of the Pareto front. This is an advantage of the method over some aggregation methods such as the well-known weighted sum method.

4.5 Numerical results for synchronous reluctance machine

The proposed algorithm to obtain an approximation of a Pareto front consists in a loop over different weights w where each iteration uses an algorithm similar to Algorithm 1 to obtain an optimized design. In contrast to Algorithm 1, however, here the descent direction is obtained by solving the auxiliary optimization problem (10) rather than an auxiliary boundary value problem of the form (8). The algorithm reads as follows:

Algorithm 2. Given initial design Ω0, cost functions 

J1,J2, tolerance tol, set of weights {w1,…,wM}.For j = 1,…,M:1) If j > M then stop else set 

J1˜←J1, J2˜←wjJ2.2) Set 

k←0, Ω0(j)←Ω03) For k = 0,1,2,…(i) Solve state equation (3) and adjoint equation (7) (ii) Compute shape derivatives 

dJ1˜(Ωk(j);V), dJ2˜(Ωk(j);V) (iii) Compute bi-objective descent direction Wh as solution to (10) with 

dJ1˜(Ωk(j);·), dJ2˜(Ωk(j);·) (iv) If ‖Wk‖ < tol then j ← j + 1 and go to 1)else set 

Ωk+1(j)=(id+tW)(Ωk(j))‖Wk‖ < tol where 

t=max{1,12,14,18,…} such that 

Ji(Ωk+1(j))<Ji(Ωk(j)), i = 1, 2.

In our implementation, we solved the quadratic optimization problem involving linear inequality constraints (10) by means of a sequential least squares programming optimization algorithm using the functionality scipy.optimize(…). To reduce computation time, we restricted problem (10) to the degrees of freedom on the material interfaces that are subject to optimization and neglected the interior degrees of freedom. This is motivated by the fact that a movement of points inside a subdomain does not alter the shape. Proceeding like this, we obtain a deformation vector field that is only supported on the material interfaces and vanishes on all interior mesh nodes. To avoid intersection of the mesh when updating the geometry, we extend the vector field from the interfaces to the whole rotor domain by harmonic extension, i.e. by solving an elliptic PDE. As additional constraints, we imposed the linear equality constraints that the normal component of the vector field on the boundary on the rotor domain vanishes, i.e. W_x(z)n_x(z) + W_y(z)n_y(z) = 0 for all mesh points z ∈ ∂D_rot. These constraints ensure that the radius of the rotor remains unchanged.

Figure 8 (left) shows the results of the bi-objective descent algorithm for minimizing the negative torque and w times the volume, J1(Ω)=−T(Ω) and wJ2(Ω)=wVol(Ω), for different choices of the weighting factor w. The right picture of Figure 8 depicts a zoom on the obtained Pareto optimal points. The computational effort for obtaining one Pareto optimal design is comparable to the cost of one single-objective optimization run (Section 3.3), amounting to a computational time of about 2 h on a single core to obtain the depicted Pareto front. The Pareto optimal designs corresponding to three different choices of w can be seen in Figure 9. The Pareto front appears to be convex in this problem. However, as discussed in Section 4.4, this is not a necessary requirement for the bi-descent algorithm to deliver satisfactory results.

5. Conclusion and outlook

We addressed the problem of finding the optimal shape of the rotor of a synchronous reluctance machine as used in an X-ray tube by means of a gradient-based free-form shape optimization method, which is based on the shape derivative. This approach allowed to obtain an optimized shape, which exhibits an increase of the torque by 26% within only several minutes of computation time. The results are confirmed by a geometric parameter optimization in JMAG where the parametrization is motivated by the design obtained by free-form optimization. Moreover, we introduced an extension to the setting of multi-objective shape optimization and showed a way to obtain an approximate Pareto front while significantly decreasing the computation time when compared to evolutionary algorithms. While we applied the methods to a specific motor design featuring a large air gap, the shape optimization approaches are universal and can easily be applied to more common rotor geometries, see (Merkel et al., 2021) for an application to a permanent magnet-assisted synchronous machine.

One aspect to further improve the proposed method is the incorporation of a mechanism to create evenly spaced solutions on the Pareto front. Here, on the one hand, adaptive choices of the initial designs or weighting parameters discussed (Section 4.3) in the spirit of Eichfelder (2009) seems promising, or, on the other hand, a predictor-corrector scheme as used in Schmidt and Schulz (2008). However, the extension of these approaches to the case of free-form shape optimization requires a more detailed analysis and is subject to future work. Moreover, in this paper, we only considered shape optimization approaches that cannot alter the connectivity of the initial design. A next step would be to consider topology optimization methods to additionally allow for changing topologies, in particular in the context of multi-objective optimization. While this was beyond the scope of this paper, also this is subject of future work.

Figures

Figure 1.

Upper half of synchronous reluctance machine with a three-phase, two-pole stator

Figure 2.

Vector diagram of a synchronous reluctance machine for the simplified model in d-q reference frame

Figure 3.

Illustration of domain perturbation by the action of a smooth vector field, Ω_t = (id + tV)(Ω)

Figure 4.

Top: initial design of rotor, T = 1.007 Nm. Center: initial design together with descent vector field W. Bottom: optimized design obtained after 70 iterations of Algorithm 1, T = 1.270 Nm

Figure 5.

Geometric parameters used for genetic algorithm

Figure 6.

Best results obtained by genetic algorithm in JMAG based on geometric parametrization of Figure 5 after 300 generations

Figure 7.

Pareto front of problem (11) and points obtained by bi-descent algorithm for six different initial values and three different weighting factors of f₂

Figure 8.

Left: Values of different designs obtained in the course of gradient based two-objective optimization algorithm for different weights w. Right: Zoom on approximated Pareto front

Figure 9.

Three different designs obtained on the approximated Pareto front by using different weights w for J2

Table 1.

Example case machine design parameters

Parameter	Value
Stator
Inner radius	26.5 mm
Outer radius	47.5 mm
Number of slots	24
Number of phases	3
Number of poles	2
Axial length	50 mm
Winding type	Single-layer distributed
No. of turns per slot	64
Phase resistance Rs.20°C	7.1 Ω
Rated voltage Ueff	230 Vac/400 Vac
Connection	Star
Rotor
Outer radius	18.5 mm

Table 2.

The current values for each winding

U-phase	V-phase	W-phase
12 A	−6 A	−6 A