Geometric design of electric motors using adjoint-based shape optimization

Interest in electric aircraft has increased due to developments in electric propulsion technology and concerns regarding aircraft carbon emissions. The emerging urban air mobility industry aims to provide convenient short-range air travel using electric aircraft. An important factor in the design of electric aircraft is the modeling and design of electric motors. The many degrees of freedom in electric motor design make it a complex design problem. To mitigate this complexity, we have developed an adjoint-based electric motor design methodology using free-form deformation for geometry parameterization, PDE-based mesh warping with exact derivatives for mesh manipulation, and ﬁnite element analysis for electromagnetic modeling. This paper highlights the approach and details of the proposed method and presents results from its applicationtoarepresentativemotordesignproblem.Theoptimizationresultsina35% decreaseinmotormassanda3%increaseinefﬁciencycomparedtothebaselinedesign. Theseresultsdemonstratetheefﬁcacyofanadjoint-basedoptimizationapproachwith exactanalyticalderivativesforelectricmotordesign.


Introduction
Interest in electric aircraft has increased over the years.This is heavily influenced by the current trends in aircraft emissions.The aviation industry is a significant contributor to carbon emissions worldwide, with jet fuel usage accounting for about 7% of world oil consumption and just under 3% of overall carbon emissions for energy use (Schäfer et al. 2019).It is expected that by 2050, at the current unchecked rate, the aviation industry will be responsible for roughly 25% of global carbon emissions (Graver et al. 2018).Advancements in electric propulsion technology have also prompted interest in aircraft electrification.Electric powertrain technologies, such as batteries and electric motors, have seen steady improvement over the years (Koh and Magee 2008;Crabtree et al. 2015).These technological improvements have brought realism to the idea of electric aircraft.
Urban air mobility (UAM) is emerging as a novel form of transportation, with the promise of shortening commutes and providing convenience for short-range air travel using electric aircraft.The vision for UAM is to provide a safe and reliable air transportation system for people, goods, and services (Cohen et al. 2021).By utilizing air transportation, it is anticipated that commuting times can be significantly reduced (Holden and Goel 2016).
Interest in aircraft electrification motivates research into the modeling and design of new, aviation-class electric motors.The electric motor design problem involves many geometric and operational parameters.To deal with this high dimensional parameter space, optimization is commonly utilized for electric motor design.A variety of motor design optimization approaches have been proposed, and these vary in application and scope.System-level optimization approaches often focus on modeling the influence of high-level motor parameters, such as outer radius and yoke flux, to capture general trends.System-level electric motor design processes, as seen for an electric aircraft (Ragot et al. 2006) and for an electric vehicle (Ahn et al. 2015), utilize analytical and simplified relationships to model the electric motor.However, these design processes focus on the system level, therefore limiting the modeling scope of the motor.
The majority of motor design optimization approaches in the literature focus on the detailed design of the motor.These methods contain both operational variables, like current and control parameters, and geometric variables, like rotor radius and air-gap thickness.Discrete parameters, such as the number of pole pairs, are also considered.Objective functions are selected to be any combination of high-level outputs significant in overall motor operation and design.Examples include torque (Parasiliti et al. 2012), output power (Di Noia et al. 2020), efficiency (Idir et al. 1998), total power losses (Grenier et al. 2021) and motor mass (Parasiliti et al. 2012;Di Noia et al. 2020).Coupling between these parameters also motivates multiobjective design problems, aimed at optimizing the geometry for a combination of these states (Parasiliti et al. 2012;Di Noia et al. 2020).
Some motor design optimization approaches focus on a subset of geometric variables to target particular adverse phenomena of motor operation.As a result, the optimal motor design is determined by a low-level aspect of motor operation, as opposed to a high-level parameter like torque or mass.For example, motor optimization approaches in Chai et al. (2016), Yamazaki and Ishigami (2010), Pellegrino and Cupertino (2010) and Shin et al. (2007) isolate the motor design problem to focus on the rotor structure and magnet shape.These parameters heavily dictate the harmonic behavior of motors, and minimization of harmonic effects is important to reduce iron losses and torque ripple.
Most of the existing motor design optimization approaches utilize gradient-free algorithms (Parasiliti et al. 2012;Idir et al. 1998;Chai et al. 2016;Yamazaki and Ishigami 2010;Pellegrino and Cupertino 2010;Di Noia et al. 2020) or incorporate derivatives in a simplified and inefficient manner (Grenier et al. 2021;Shin et al. 2007).The derivative computation cost of high-dimensional problems, typically with tens of variables, exhibits poor scaling with the number of design variables (Martins and Ning 2022); this cost can be decreased by reducing the number of design variables, but this sacrifices the optimization problem scope.A gradient-based approach is required to make large-scale optimization problems tractable to solve.
Although uncommon, there have been cases of solving motor design optimization problems using a gradient-based approach.The sequence of studies in Gangl et al. (2015Gangl et al. ( , 2016) ) utilize shape derivatives and topological derivatives to optimally design the rotor and magnet structures to minimize total harmonic distortion.However, a gradient-based formulation alone is not enough to improve the scaling of computation time relative to the number of design considerations.Despite the abundance of existing electric aircraft motor design approaches, the current methods are inefficient in modeling the effects of the full motor design space or are too limited in scope.
Gradient-based optimization approaches using the adjoint method overcome these limitations.With the adjoint method, the cost of derivative computation scales with the number of outputs rather than the number of design variables; for optimization problems with a single objective, updating the design variables using the adjoint method is a more efficient approach (Martins and Ning 2022;Martins and Hwang 2013;Gray et al. 2019).Adjoint-based sensitivity analysis has been applied to efficiently optimize systems with many degrees of freedom in various engineering problems, such as topology optimization (Yan et al. 2022), aerodynamic shape optimization (Chauhan and Martins 2021) and multidisciplinary aircraft design optimization (Sarojini et al. 2023).
We have developed a motor design optimization methodology for high-dimensional design problems using a collection of geometry parameterization, mesh manipulation and motor analysis methods centered around adjoint-based sensitivity analysis.In the proposed methodology, the geometry is parameterized using free-form deformation and the mesh is updated using PDE-based mesh warping.A semi-adaptive load-stepping approach is applied to ensure mesh warping convergence, with the geometry parameterized in a manner that is amenable to exact analytical derivatives.Finite element analysis is utilized to solve Maxwell's equations over a full motor mesh.This paper highlights the details and approach for the presented motor design methodology.We demonstrate the methodology by performing design optimization of a three-phase permanent magnet synchronous motor for a range of operating conditions.We compare results from a baseline case and a comprehensive geometric design case to quantitatively assess the efficacy of using adjoint-based optimization for electric motor design.
Recently, Babcock et al. independently developed a similar methodology using adjoint-based derivative computation for motor design optimization with electrothermal coupling (Babcock et al. 2023).There are differences between their methodology and the methodology proposed here, such as the electromagnetic modeling approach, the material model and the electromagnetic output computation method.
Most notably, the work from Babcock et al. utilized the finite-difference method to capture sensitivities with respect to changes in the mesh and geometry, while the methodology proposed here uses exact analytical derivatives coming from the geometry parameterization and mesh warping methods.The current methodology consists of a different set of methods with end-to-end analytical derivatives to solve the highdimensional motor design problem.
The structure of the paper is as follows.Section 2 discusses electric motor design background, including motor types, and various motor analysis fidelity levels.Section 3 outlines the motor design methodology details, divided into theoretical and computational aspects.Section 4 provides results, including verification of the electromagnetic solver, a grid independence study, and a system-level demonstration of the presented methodology.Finally, Sect. 5 concludes the paper and discusses future work.

Background
This section provides background on high-level motor design and analysis choices.The demands of electric aircraft propulsion place strict limitations on the type of electric motor.As a result, these discrete choices must be made early on in the design stage.First, we address the benefits and drawbacks of various electric motor types in electric aircraft operating conditions.Next, we compare different motor flux orientations.Finally, we discuss the capabilities of various motor analysis approaches across a range of fidelity levels.

Comparison of motor types
There are strict requirements for aviation-grade electric motors; these include compactness, high power density, efficient operation at high speed conditions, and minimal noise footprint.There are numerous standard motor topologies used in various applications; these include DC (brushed or brushless), induction, switched reluctance, and permanent magnet motors.Of the existing motor types, we need to select one that is well suited for the performance requirements of aircraft operation.
Studies conducted in Bolam et al. (2020), Tenconi et al. (2014), Yildirim et al. (2014), Hashemnia and Asaei (2008) and Zeraoulia et al. (2006) highlight the benefits and drawbacks of each motor type.DC motors, for example, provide simple control and field weakening capabilities, but construction and maintenance of the motor structure pose difficulties with usage; these motors serve well in low-power applications (Hashemnia and Asaei 2008).Induction motors are the most mature motors running on AC power, and are favored in general applications for their reliability and robustness.However, this class of motors suffers from lower efficiency at high speed due to the large rotor winding power loss, and is not well suited for the high-speed demands of aircraft operation (Hashemnia and Asaei 2008).Switched reluctance motors are rigid and operate reliably under a wide constant power region.However, this motor type suffers from high torque ripple and electromagnetic interference (Hashemnia and Asaei 2008;Yildirim et al. 2014).Permanent magnet synchronous motors (PMSMs) are known for their high power density and efficient operation in high-speed conditions; however, the reliance on permanent magnets in the rotor puts them at risk of demagnetization at high temperatures.Compared to other motors, these motors contain better thermal management systems than induction motors to mitigate excessive joule losses (Bolam et al. 2020;Tenconi et al. 2014;Zeraoulia et al. 2006) and less harmonic and vibration effects compared to SRMs due to less torque ripple (Hashemnia and Asaei 2008;Tenconi et al. 2014;Yildirim et al. 2014;Zeraoulia et al. 2006).Based on the studies mentioned above, we demonstrate our design methodology using a PMSM topology; among the available motor types, it is best suited for the demands of electric aircraft operation.

Motor flux orientation comparison
Another significant design choice that impacts the performance and application of electric motors is the magnetic flux orientation.The two types are radial flux, in which the magnetic flux between the rotor and stator is oriented radially, and axial flux, where the flux is oriented along the length of the motor.Selecting a flux orientation for an electric aircraft motor is not as straightforward as the motor choice, as both types have their benefits and tradeoffs.
For example, studies done in Cavagnino et al. (2002) and Zhang et al. (2014) show that, when motor size is heavily restricted, axial flux motors have higher output torque and torque densities than radial flux motors.This suggests that axial flux motors are well suited for electric aircraft applications.However, other attributes of the two flux orientations need to be considered for aircraft operation as well.A performance comparison done by Parviainen et al. (2005) suggests that, in operating cases of identical current density and electric loading, radial flux motors operate with higher efficiency than axial flux motors.For an axial flux motor to match the efficiency of a radial flux motor, a significant volume increase is required.Additionally, it is impractical to produce axial flux motors at short lengths due to the space of the end connections.Work done in Pippuri et al. (2013) compares the torque density at a baseline power and speed condition between axial and radial flux motors; it was found that, for the maximum efficiency case and most compact case, the torque density of the radial flux motor exceeded that of the axial flux motor.
Another important aspect of aircraft design is noise; initial studies have shown that radial flux motors have reduced noise profiles compared to axial flux motors.Wei et al. (2022) compared radial and axial flux motors with identical power demand and varying pole-slot combinations to assess the effects on noise; it was found that the axial flux machines possess more complex noise profiles with higher dominant harmonic modes than radial flux motors.However, limited comparisons of noise profiles between radial and axial flux motors have been performed.Further studies must be conducted to conclude which flux orientation has lower harmonic effects.
Although aviation-class electric motors require high output power and torque densities, the practicality and operation of motors for both flux directions must be taken into account.Based on the need for compactness and high efficiency across a wide range of operating conditions, as shown in the aforementioned studies, we have chosen to demonstrate our motor design optimization methodology using a radial flux PMSM.

Motor modeling approaches
Electric motor modeling approaches span a wide array of methods and levels of fidelity.The level of desired accuracy is dictated by many factors, such as problem scope and application.Electric motor models can be divided into three main categories: simple analytical models, equivalent circuit models, and finite-element-based models.
The lowest level of fidelity of motor models uses simple analytical equations or surrogate models; these approaches simplify the motor representation into its most fundamental aspects and are utilized in the design and optimization of systems with many variables.McDonald (2013) developed a method to model motor performance based on data on motor efficiency and power loss; the approach uses polynomial fitting to model efficiency and power loss as a function of torque and speed.The work conducted in Dantsker et al. (2019) describes an optimization approach for the design of an unmanned solar-powered aircraft, modeling the electric motor using Ohm's law and standard analytical equations to calculate back electromotive force (emf).At a similar level of fidelity, Falck et al. (2017) performed trajectory optimization of the NASA X-57 Maxwell aircraft and modeled motor efficiency using data interpolation.Analytical relationships relate the efficiency to the motor power losses, and the heat dissipation and temperature are determined to assess the thermal constraints.
Although the analytical and surrogate approaches are sufficient for initial motor modeling and design, they are heavily limited by the lack of generality and extensibility to different motor geometries.Analytical approaches rely on assumptions about the motor operation, and behavior; as a result, higher-order effects are often ignored and only the general trends are captured.Additionally, these models do not extend well to general motor analysis, as the modeling assumptions are specific to motor geometries.
Magnetic equivalent circuit (MEC) models represent an improvement in fidelity level from the analytical approach.MEC models benefit from a short computation time and reduced computational complexity while still providing reasonable accuracy; as a result, these methods are desirable for the initial rapid design and analysis of existing motor topologies.MEC methods discretize different zones of the motor geometry into circuit components, and a basic circuit analysis approach is applied to model the magnetic flux.Source terms within standard equivalent circuit methods represent motor flux sources, and magnetic reluctances denote field hindrance in areas like the air gap.Sebastian et al. (1986) devised an MEC approach modeling two PMSMs with different magnet inset layouts within the rotor.The analytical MEC results showed promising qualitative correlation with the experimental torque data across a sweep of current lead angles, while the inductance errors were within 2%; the results also highlight the versatility of the MEC approach to represent motor geometry.Further development of the MEC approach expanded the considerations within the equivalent circuit to capture more complex phenomena; work done in Ba et al. (2022), Wijenayake and Schmidt (1997) and Sheikh-Ghalavand et al. (2010) highlights the significance of considering additional equivalent circuit components to estimate saturation and core losses, both of which play a significant role in determining overall power loss and motor performance.Upon introducing additional components to model more complex motor phenomena, validation with results from finite element analysis and experimental testing shows an improvement in the accuracy of the MEC approach.
Despite the extensibility of MEC to capture additional pertinent physical phenomena of motor operation, this approach models the motor in a nodal fashion.Additionally, parameters must be averaged over regions of the motor geometry to adhere to the nodal nature of the MEC method.As a result, the method is unable to capture the full flux field effects, and information about the flux fields and motor behavior is lost, leading to inaccuracies in parameters like core loss and saturation.This limits the MEC approach when it comes to analyzing completely new motor designs (Yilmaz and Krein 2008).For a full in-depth analysis and understanding of the flux field, a finite element analysis approach is required to model the comprehensive behavior of the motor.
High-fidelity analysis is conducted using the finite element analysis (FEA) approach by solving a form of Maxwell's equations over a mesh representation of a motor.The solution field is directly used to calculate important outputs like flux linkage, core losses, and torque.An additional benefit of the FEA approach is the dependence of the solution field on the geometry; this captures detailed effects of the flux field and motor operation and provides a more meaningful comparison between the behavior of different motor topologies (Yilmaz and Krein 2008).As a result, this approach provides more insight when conducting design sweeps of different motor topologies for discrete parameters such as pole pairs, and stator teeth, or assessing new motor structures.Work done in Güemes et al. (2011), Powell et al. (2003) and Hebala et al. (2021) explores these parameters using the FEA approach, highlighting the benefits and drawbacks of different motor designs.
Additionally, studies from Grenier et al. (2021), Powell et al. (2003), Lin et al. (2018), Hebala et al. (2021Hebala et al. ( , 2022) ) simulated electric motors in aircraft operating conditions using a high-fidelity approach; these studies are aimed towards motor designs of high power density and efficiency, key requirements for electric motors in aircraft operating conditions.The FEA approach provides a comprehensive output of the flux fields; as a result, many critical phenomena in motor operation can be analyzed, and the motor design approach can focus on one or more of these phenomena.These include the minimization of AC losses due to the high-frequency nature of motors (Hebala et al. 2022) and the analysis and mitigation of cogging torque or torque ripple, a key output for understanding the motor harmonics and vibrations (Hebala et al. 2020(Hebala et al. , 2021;;Lin et al. 2018;Güemes et al. 2011).
The accuracy and computational efficiency across motor modeling approaches of varying fidelity levels can be leveraged in a multi-fidelity approach.Ragot et al. (2006) demonstrate this in the optimization of a solar airplane; a simplified MEC model determines the air gap field strength, while other motor parameters are found using analytical relationships.The FEA approach can also be used to inform low-fidelity models.Elsherbiny et al. (2022) utilize results from FEA simulations as a pre-processing calibration step to generate maps for flux linkage, torque, and inductance as a function of current components; these maps are utilized as an inverse function in a MEC model to estimate torque and power losses.Methods that directly couple MEC and FEA also exist, where both approaches are utilized simultaneously to leverage the computational benefits of MEC with the accuracy of FEA in different zones of the motor geometry (Nedjar et al. 2012).However, calibration and coupling in multi-fidelity approaches are often too dependent on either the operating case or the geometry; thus, the approach is not modular enough for general motor design.In addition, the use of coupling introduces implementation complexity with potentially little gain in accuracy.
For our motor design optimization methodology, we are interested in an accurate motor analysis method that is extensible to various motor geometries.Based on the above studies comparing different fidelity levels of various motor analysis models, we select the high-fidelity FEA approach for our motor analysis model.

Methodology
The details of the novel motor design optimization methodology are explained in this section.We first outline the theoretical approach, and then discuss the implementation and computational aspects of the workflow.We then combine these details to present the complete motor design methodology.

Mathematical model
In this subsection, we outline the theory behind the motor design and shape optimization methodology.First, we discuss free-form deformation, a differentiable approach to parameterizing the motor geometry.Then, we formulate the governing electromagnetic modeling PDE from the standard Maxwell's equations and outline the approach for modeling magnetic permeability in a differentiable manner.Finally, we summarize the post-processing formulation used to compute high-level outputs from the FEA solution.

Free-form deformation
To parameterize and deform the subdomains of our geometry in a computationally efficient manner, we use a technique to manipulate solid geometries called free-form deformation (FFD).This technique is significantly more cost effective than physically motivated approaches, such as linear spring models, or finite-element methods, and is easier to implement (Moore and Molloy 2007;Gibson and Mirtich 1997).For simple geometric movements, where the underlying mechanics of the geometric deformation are not necessary, FFD provides low computational cost and high control for altering the geometric structure (Moore and Molloy 2007;Gibson and Mirtich 1997).
FFD applies geometric deformations using parametric curve representations, such as B-splines or Bézier curves, to parameterize and deform a grid containing the geometry (Sederberg and Parry 1986;Kenway et al. 2010;Gibson and Mirtich 1997).In practice, FFD is applied by surrounding a geometry with a pseudo-structure defined by control points; the geometry is embedded within the structure and is defined paramet-Fig. 1 We shown an example of a two-dimensional geometry (in black) embedded within an FFD face (in red).When embedded into the red FFD face, the corners of the black rectangle V i can be represented parametrically in the (u, v) space.Movement of the FFD face control points P jk will warp V i based on the undeformed (u, v) coordinates.(Color figure online) rically as a function of the control points.Movement of the control points reshapes the physical form of the geometry to maintain the parametric coordinates of the embedded geometry.In two-dimensional space, we refer to these pseudo-structures as FFD faces.A simplified two-dimensional example is illustrated in Fig. 1, where mesh entities are colored in black and FFD entities used for parameterization are shaded red; V 0 , V 1 , V 2 , and V 3 represent the defining vertices of the mesh.We can formulate this parametrization in the x-y coordinate frame as p(x, y) = P 00 (x, y) + u 0 0 v (P 11 (x, y) − P 00 (x, y)), (1) where P 00 and P 11 represent the FFD control points, and u, v are the parametric coordinates of the local coordinate frame defined in red.Nodes within the FFD face contain a physical coordinate within the x-y plane, along with a set of constant parametric coordinates in the u-v plane.As the FFD control points P 00 and P 11 move, the embedded geometry defined by V i will deform in physical space in a way that maintains its parametric coordinates within the FFD face.In addition, the FFD parameterization is prescribed within the geometry changes; thus, we can assign design variables corresponding to the changes of the geometry rather than directly defining the geometry itself.

Electromagnetic modeling
We model the motor by solving Maxwell's equations: where Ampere's circuital law (2) governs the field behavior based on current sources J, and the permeability μ relates the magnetic flux density B and the magnetizing field H (4). We rely on two assumptions to simplify the governing equations.First, the motor frequency is within the kilohertz range, low enough to disregard the transient effects of the displacement field D. Second, the end effects of the motor are negligible, and the operating states are invariant along the motor length; this allows the problem to be solved in the cross-sectional plane of the motor.Gauss's law of magnetism (3) states that the magnetic flux density is divergence-free, so the vector field B can be represented as a curl of a vector field: and this vector field is called the magnetic vector potential A. The variation of the magnetic flux density is purely in-plane, so (5) implies the magnetic vector potential A is one-dimensional normal to the geometry.Substituting ( 4) and ( 5) into (2), we obtain a governing equation in terms of the magnetic vector potential: where each term contains a component normal to the solution plane.Utilizing the gauge freedom of the vector potential field, we set ∇ •A z = 0 to formulate a governing equation over a two-dimensional field for A z represented as a Poisson equation: where the source term J z represents the field excitation.This term contains two distinct contributions: the input signal excitation in the stator windings J w and the flux density distribution from the magnets J m .The equivalent source term from the magnets J m can be represented as using the magnetic coercivity H c (Weiss et al. 1984).
A model for the permeability μ is required to relate the magnetic vector potential field to the source terms.In ferromagnetic materials, such as iron, permeability is a highly nonlinear property.Discrete approaches to permeability and nonlinear material modeling are commonly used in electric motor analysis due to the fitting complexity; however, these are infeasible for our application as gradient-based optimization requires smooth and continuous functions.We model this effect with smoothed piecewise fits based on material data from Ansys Maxwell (Ansys Accessed 03/20/2021).Each discrete function is selected such that the overall permeability fit qualitatively matches regular permeability curve trends (Jaafar et al. 2004;Hao et al. 2020;Gmyrek and Cavagnino 2021).To guarantee differentiability, we use a cubic function to fit the center region of the data.This allows for matching gradients at the bounds of each section, ensuring compatibility with gradient-based optimization.The resultant fit is shown in Fig. 2.

Post-processing
The post-processing step uses the outputs from FFD and FEA to determine high-level motor states such as torque, efficiency, and power loss and geometric parameters such as area.We consider two different coordinate frames in motor operation: the static abc coordinate frame and the rotating dq coordinate frame, which rotates with the rotor rotation and input signal; we use a direct-quadrature-zero (DQ0) transformation to convert between the two frames.The FEA solution is computed in the abc coordinate frame; we introduce the dq coordinate system to simplify the post-processing steps, circumventing the complexities of analyzing a rotating system in a static frame.
We use the FEA solution to calculate flux linkage, which is a measure of the flux through a surface in space found by integrating the magnetic flux density field B. The integration surface is represented by a pseudo-surface in the air gap connecting two stator winding slots of the same phase, as shown in Fig. 3; this surface is extended along the length of the motor.We reformulate the flux linkage in terms of the magnetic vector potential field A by applying Stokes' Theorem, as follows: where S and l represent a surface and a line defining the pseudo-surface, respectively.
Given that the magnetic vector potential acts normal to the motor plane, (9) can be simplified further in terms of parameters from Fig. 3 to where the subscripts s, s + 1 describe adjacent slots of the same phase.The nodal flux linkage approach (10) has been used in both the educational and research fields (Lowther and Silvester 2012;Kang et al. 2000;Salon 1995;Bianchi 2005).The total flux linkage per phase is then computed by summing each λ per pole across the entire Fig. 3 We model the air-gap with a pseudo-surface representation.The flux linkage can be represented using the magnetic vector potential A z rather than the magnetic flux density B motor as follows: where the subscript i represents one of three phases a, b, c, subscript j represents the per-pole flux linkage from (10) and N s is the number of stator slots.The per-phase flux linkage is converted from the abc frame to the dq frame with the DQ0 transform; the dq flux linkage components are used to calculate voltage, as shown below: where ω is the rotational speed and I d , I q are the dq components of the current.The flux linkage is also used to calculate electromagnetic torque, represented as: where p is the number of pole pairs.Power losses and efficiency can be found from this output torque.Four sources of power loss are considered: copper loss, core loss, windage loss, and stray loss.Copper loss is characterized as heat loss due to the current excitation in the windings, and is represented as where m is the number of phases, R is the per-phase resistance, and I is the rootmean-square current.Core losses arise from the periodically changing magnetic field within the motor core, and can be divided into two distinct factors: hysteresis loss and eddy current loss.We model these two power losses following Mi et al. (2005) taking advantage of the periodic behavior of motor operation to calculate the overall core loss.The hysteresis loss is represented as where f is the input signal frequency, k h is the hysteresis coefficient, and β is the Steinmetz coefficient, set to 1.7.The eddy current loss can be modeled similarly as where k ec is the eddy current coefficient.The windage loss represents losses due to air resistance within the air gap and can be represented as where k r is the motor surface roughness, f is the friction coefficient, ρ is the density of air, ω is the rotational speed, r r is the rotor radius, and l is the motor axial length.
The friction coefficient is a measurement of the fluid interaction with the motor core, and is dependent on the Reynolds number in the air gap: where ρ is air density, ūθ is the average tangential velocity in the air gap, δ is the air gap thickness and ν is the kinematic viscosity of air.This simplified approach to model windage loss, however, is unable to capture the effects of the complex vortex interaction and fluid instability due to the stator geometry, leading to under-prediction (Saari 1998).Stray loss P s in the motor is also considered to account for variations in the load during motor operation; we assume this to be one percent of the output power (Tong 2014).After calculating the power loss components, the effective output load torque can be determined.We use a simple torque balance model, where the torque loss is determined by the mechanical power losses and rotational speed: where τ l represents the load torque produced by the motor.Note that the effects of copper loss are neglected in the torque balance model.Although copper loss reduces the overall output power, it does not hinder the mechanical behavior and rotation of the motor; as a result, it is not considered in the torque balance.The efficiency can then be computed as a ratio between the output power and input power, where the Fig. 4 A visual representation of the shape parameter types is shown.A constant shape parameter applies uniform movement to the geometry in one direction.The linear shape parameter divides the prescribed deformation linearly across the geometry in a direction; this will maintain the symmetry of the deformation output power is the product of rotational speed and effective output torque and the input power is equal to the output power and the sum of the power losses: η = τ l ω τ l ω + P cu + P h + P ec + P w + P s (21)

Computational approach
In this subsection, we discuss the computational tools and aspects of the motor design methodology.These details of the model build upon the theoretical foundation described above in Sect.3.1.

Geometry parameterization and free-form deformation
To perform FFD and mesh warping, we need a robust mesh generator that provides detailed mesh information for geometry parameterization.We use GMSH (Geuzaine and Remacle 2009), a versatile three-dimensional finite element mesh generator with a wide range of capabilities, such as classical geometry generation, boolean operations, and object assignment.We utilize this software to generate the motor mesh and extract detailed nodal information for the initial geometry parameterization.
As outlined in Sect.3.1.1,FFD faces are overlaid onto the geometry around areas of interest for shape optimization.To connect design variables to the geometry parameterization, we assign shape parameters to the FFD faces; we have defined two types of shape parameters: constant and linear.A constant-type shape parameter applies a constant shift to the FFD face; this operation corresponds to a fixed translation in space.A linear-type shape parameter, on the other hand, divides the transformation across the FFD face in a linear fashion, acting symmetrically across the selected dimension; this is synonymous with stretching or compressing the geometry.Figure 4 shows a visual representation of the two shape parameter types, where bold variables denote shape parameters and the ˆsymbol denotes the geometric change of that shape parameter.
The deformation of the FFD face changes the geometry; a simple parameterization example is visualized in Fig. 5. First, the movement is propagated onto the major vertices of the geometry subdomains V 0 , V 1 , V 2 , V 3 ; in the parametric coordinate Fig. 5 We visually represent the FFD vertex and edge parameterization steps 4 and 5.The black represents the geometry, and the blue represents the geometric aspect that is being parameterized.At step 4, the parametric coordinates of the vertices defining the geometry are found.Then during step 5, a one-dimensional parameterization for nodes along the geometry edges is conducted.This maps the changes in the geometric vertices to the boundary nodes on the rest of the geometry frame, these vertices have constant u-v coordinates.The movement of the control points P 00 , P 01 , P 10 , P 11 will alter the physical coordinate of the subdomain vertices V 0 , V 1 , V 2 , V 3 to maintain their original u-v coordinates in the parametric frame, as shown on the left in Fig. 5.The subsequent movement of the subdomain vertices is then applied to the nodes along each of the subdomain edges, parameterized similarly.As shown on the right, the movement of the vertices V 1 , V 2 will alter the position the nodes e 10 , e 11 , e 12 , e 13 of edge E 1 whilst maintaining a constant parametric coordinate w along the edge.We can consider FFD as a mapping ψ from the original geometry to an altered geometry, where ψ contains the parameterization of the subdomains relative to the FFD faces.For all subdomain edges, the new edge nodes can be defined as ēi j = ψ(e i j ) and the changes in the geometry can be defined as δ e = ēi j − e i j = ψ(e i j ) − e i j .
The full set of computational steps for FFD is summarized below.
1.A region of the geometry is embedded within an FFD face.A combination of shape parameters is set to define the allowable deformations of the FFD face.2. Geometric design variables related to the embedded geometry are set; these design variables define significant parameters of the geometry (for example, the rotor radius in a motor).3. Analytical relationships between the design variable and a set of shape parameters are defined; these relationships between the design variables and shape parameters define a mapping that represents the change in the variable to the representation of this change on the geometry, altering the location of the FFD control points.4. The movement of the FFD control points warps the FFD face; the deformation of the FFD face shifts the major vertices of the embedded surface V i to maintain the parametric coordinates.5.The movement of the major vertices is mapped onto the edges E i by moving the edge nodes e i j that define the embedded surface; the edge nodes will shift accordingly to maintain their parametric representation along the curve.
The approach presented here allows changes in geometric variables of interest in motor design to be mapped to the boundaries defining the structure of the motor geometry in an efficient and differentiable manner.

PDE-based mesh warping
We utilize a PDE-based approach for mapping the discretization of the major features (e.g., magnet shape, outer boundary) to the mesh for the electromagnetic analysis.
For simple subdomain movements, as shown in the previous section, non-physical techniques like splines and FFD are simple and cost-effective.However, for complex mesh deformations, non-physical techniques are not as robust and prone to generating skewed mesh elements (Moore and Molloy 2007).
Interpolation-based methods, such as inverse distance weighting (IDW) and radial basis functions (RBFs), lack a physical analogy but are commonly used and computationally robust.However, these approximation methods require tuning parameters and functions beforehand, and have their own shortcomings related to error and feasibility.Solving the system of equations generated with RBFs for problems with more than a few thousand data points becomes very costly and memory intensive (Selim and Koomullil 2016;Amidror 2002).IDW interpolation avoids solving a system of equations but still requires summing over all of the nodes, and the accuracy suffers near the control points due to the vanishing derivatives (Luke et al. 2012;Amidror 2002).
PDE-based approaches using finite element methods (FEM), compared to linear spring models, are more flexible with complex meshes, boundary conditions and material types; they are also less prone to matrix stiffness, solver instability and localized interactions, which can all lead to poor mesh element shapes (Moore and Molloy 2007;Gibson and Mirtich 1997).Although the computational expense of FEM approaches is significant, we opt for this method given the importance of a high quality mesh.
As explained in the previous section, FFD is implemented to model the movement of subdomain boundaries comprising the mesh outline.The mesh warping step updates the mesh based on the geometric changes computed during FFD.The geometric deformation computed from FFD is utilized as the boundary condition for the mesh warping.The state variable in the mesh warping technique is the mesh node displacements, denoted as u and defined as where X represents the original position of the mesh nodes and φ is a mapping to the new mesh node locations.The finite element domain is denoted by , and indicates the boundaries of the finite element domain.We formulate the mesh warping problem as a hyperelasticity problem, represented as where the boundary condition term u bc = δ e is the boundary movement given by FFD, and f represents a body force, which is equal to zero.Additionally, P is the first Piola-Kirchhoff stress tensor, defined as where F is the deformation gradient given as and I is the identity matrix.The St. Venant-Kirchhoff constitutive model is used for the material model, in which the strain energy function is defined as and the Green-Lagrange strain function E is defined as Thus, the first Piola-Kirchhoff stress tensor in (24) can be reformulated as in terms of the Green-Lagrange strain function.The mesh warping problem is solved using FEniCSx (Baratta et al. 2023), an open-source platform for solving PDEs; FEniCSx contains derivative computation capabilities, a necessity for gradient-based optimization.We interface to FEniCSx using femo,1 a general framework for solving finite element-based optimization problems.The hyperelasticity approach is highly nonlinear and presents difficulties related to solver robustness and convergence.We use a method known as load-stepping, in which the boundary condition data for the mesh movement is divided into small steps.As a result, the mesh warping algorithm solves a sequence of nonlinear problems with smaller deformations.The reduced step size makes the problem more linear, and the solver is less susceptible to divergence.The number of steps is computed based on mesh element qualities and maximum FFD boundary movement.The hyperelastic nature of the approach also makes each load step independent of the solution history.The load-stepping approach has been applied adaptively; upon solver divergence, the mesh reverts to the previous successful solution.The step size is then recalculated based on updated mesh characteristics and the remaining boundary movement.This allows for automated correction during the load-stepping process.This process is visualized in Fig. 6.Fig. 6 A visual representation of the adaptive load stepping approach is shown here.The total load is divided into a discrete number of steps.Upon solver divergence, a backward step is taken to recover the last successful solve, and a new step size is calculated.We refer to this process as the adaptive recalculation (AR) step.This process is repeated for "i" instances until the total load is applied, indicated by the dashed green line.To avoid recurring solver divergence, it should be noted that Treatment of the nonlinear solver is important for convergence, as this influences the accuracy of both model evaluation and derivative computation.Divergence signifies that the mesh did not converge to the correct geometry based on the boundary information.As a result, the post-processing computations that are reliant on the geometry are both compromised.However, the influence of solver convergence on the model derivatives plays an arguably more significant role.The derivatives provided by the FEniCSx solution depend on the convergence status of the mesh warping solver; this influences the trajectory of the design variables, causing the design to divert from the optimal solution.Furthermore, critical failure of the optimization can occur upon the existence of infs and NaNs.

Electromagnetic modeling
Similar to the PDE-based mesh warping, FEniCSx is utilized to solve the electromagnetic problem over the motor mesh.We use (7) to formulate the electromagnetic finite element problem as where we set g = 0 to represent dissipation of the magnetic vector potential A z at the boundaries.This boundary condition requires special attention due to the boundary movement from FFD.We format (29) into the weak form by multiplying each side with a test function v and integrating over the domain.We now have a partial-differential equation in the weak form represented as: which can be implemented in FEniCSx.Nitsche's method is applied to handle the treatment of A z at the mesh boundaries due to the boundary movement from FFD.First, the PDE on the undeformed mesh is considered.By moving both terms of (30) to the same side and applying Nitsche's method for the Dirichlet boundary condition, the PDE on the undeformed mesh can be written as in the form of a residual function, where C pen is a scaling factor, h is the cell diameter computed using the cell circumradius, and n represents the facet normal direction.The nonlinear solver drives this residual to zero to satisfy the boundary conditions due to changes in the geometry and mesh.The approach to determine the variational form for the deformed mesh uses a coupling method that incorporates u from ( 23) into ( 31) by modifying the differential operators using Nanson's formula; the details of this are approach are further outlined in Bazilevs et al. (2008).

Optimization framework software
The motor design optimization framework is developed using two open-source software packages.The model is built using the Computational System Design Language (CSDL) (Gandarillas et al. 2022), a Python-based embedded domain-specific language designed to facilitate the modeling and optimization of complex multidisciplinary systems.CSDL contains a standard library of operations with automated derivatives; as a result, the complexity of model implementation and derivative computation in CSDL is low.CSDL also automates adjoint-based derivative computation.Any models generated with the CSDL front end will automatically compute derivatives using the adjoint method during optimization.To connect the workflow in CSDL to an optimizer, we use modOpt,2 an environment developed to interface with and design optimization algorithms.The complexity and demand of this workflow require a robust optimizer capable of gradient-based optimization.We utilize SNOPT (Gill et al. 2005), an optimizer based on the SQP method.

Implementation details
The mathematical models from Sect.3.1 and details of the computational approach from Sect.3.2 are combined to assemble the general motor design framework; each discipline is shown in Fig. 7, along with the inputs and outputs that connect each discipline.The geometry pre-processing step generates the mesh and FFD parameterization data.In the current approach, the rotor torque and rotational speed inputs are constant, Fig. 7 The design workflow for the motor methodology is shown.The pre-processing steps in orange provide information for the four main evaluation models, shown in green.We denote the outputs from the rotor performance model using the subscript r to refer to the rotor torque τ r and rotor rotational speed ω r .We use δ g as shorthand for a general set of geometric design variables and δ e to represent the changes in subdomain edge position given from FFD.The terms f , c represent the objective and constraints, respectively.(Color figure online) so rotor-motor coupling is neglected.Therefore, we incorporate a rotor performance model outside of the optimization loop to evaluate the operating torque.In general, however, we can easily consider this coupling by including a rotor performance model within the optimization loop.

Results
In this section, we demonstrate the proposed motor design optimization methodology.First, we present verification results and quantify any existing errors in the proposed method.Next, we conduct a grid independence study to determine reasonable element sizes for the mesh.Finally, we demonstrate design optimization results for a motor within an electric aircraft system.All results are generated for a 12-pole, 36-slot radial flux permanent magnet synchronous motor; an example of the geometry is shown in Fig. 8.

Model validation
The FEniCSx electromagnetic solver is verified using the Ansys Maxwell software (Ansys Accessed 03/20/2021).We begin with a qualitative comparison of the output magnetic flux density fields, shown in Fig. 9.The two cases considered here are the no-load (no current) case and the 280 A current amplitude case.We see qualitative similarities in the no-load flux density field distributions between the FEniCSx solver in Fig. 9a and Ansys Maxwell in Fig. 9b.The flux field dissipates in a similar fashion around the boundaries and between stator teeth, and the flux density concentration around the magnets indicates that the modeling approach using (8) for the magnets is accurate.However, higher saturation values of the flux density exist around the magnets in the Ansys Maxwell results.Comparing Fig. 9c and d, we see similar trends for the 280 A load case, with larger discrepancies in the saturation areas.
To assess the error in the magnetic flux density fields shown in Fig. 9, the output torque from the FEniCSx solution and the Ansys Maxwell solution are compared.A no-load case produces no torque; therefore, we compare the output torque using an additional load case of 100 A. The results are shown in Table 1.We see a larger torque error as the current amplitude increases; the error does not exceed 14% within the tested cases.We believe the sources of error in these results originate from the nonlinear permeability fitting within the motor core and the nodal evaluation approach for calculating the torque from (10); this is further discussed in Sect. 5 under future work.For the current study, we accept this error as reasonable to demonstrate the design methodology.

Grid independence study
We conduct a grid independence study to assess the influence of mesh element sizes on the model outputs.We utilize these results to select an ideal grid size that leverages accuracy and computational cost.The prescribed element sizes are partitioned into four distinct values, varying radially; the element sizes are not controlled globally across the mesh, but rather around points that define the structure and subdomains of the motor geometry.This is shown in Fig. 10a, where the prescribed element sizes l 1 , l 2 , l 3 , l 4 vary radially based on the boundaries of the motor structure.The values of the element sizes are determined by the dimensions and proximity of the surrounding subdomains.This necessitates smaller element sizes in regions like the air gap and stator teeth, where at least one dimension of the subdomain is small.Larger element sizes are allowable around the mesh boundaries in the stator and rotor core, where there are larger gaps between different subdomains.
To assess the influence of mesh element size, we analyze the electromagnetic torque and core losses; these parameters depend directly on the magnetic flux field solutions.The torque is computed directly from the vector potential solution in the air gap, which assesses the validity of element sizes l 3 and l 4 .Core loss, which is the sum of hysteresis and eddy current loss, is calculated by integrating the solution field over Fig. 10 vary the mesh element sizes for the grid independence study across three levels.Each mesh corresponds to a value for l i , dictating an element size in a particular region of the geometry.The different meshes had element counts of 32,274, 11,866, and 6122, respectively Table 2 We compared the torque and core loss error percentage for the medium and coarse meshes relative to the fine mesh

Output
Medium mesh (%) Coarse mesh (%) Current load cases of I = 30 A and I = 100 A were analyzed the motor core; this output is influenced mainly by element sizes l 1 and l 2 .The grid independence study is conducted for load cases of 30 A and 100 A.
The grid convergence error results shown in Table 2 are relative to the fine mesh; the comparison for each output is done by varying the element sizes with the greatest influence in the respective regions.For both load cases, the error for both parameters increases as the mesh coarsens.The core losses exhibit noticeable errors for changes in l 1 and l 2 .The rotor and stator yokes require sufficiently fine element sizes to resolve the vector potential distribution due to the zero-boundary condition.Large mesh elements misrepresent the drop in the vector potential around the boundaries, especially with higher magnitude loads.As the load increases, more elements are required to integrate over the motor core with reasonable accuracy.The torque, on the other hand, remains relatively unchanged for variance in l 3 and l 4 .This is because the air gap is already thin; having at least one radially intermediate node is sufficient to capture the flux movement effects between the rotor and stator.Although the error percentages in the coarse mesh are acceptable, the medium mesh will be utilized; the use of small element sizes benefits the performance of the mesh warping solver, and this level of fineness can balance accuracy and computation time.

Optimization problem background
We demonstrate the motor design optimization methodology by optimizing the motor design considering metrics at the aircraft system level.The reference vehicle is the Lift + Cruise aircraft concept from NASA (Silva and Johnson 2021).We aim to optimize the pusher motor.
The influence of the electric motor in the full system-level aircraft design can be captured using the electric range equation.The electric range equation relates basic aircraft operation to electric motor properties in the cruise condition under the assumption of steady-level flight and is given by where W represents weight, R is the flight range, g is the gravitational constant, η p is the propulsive efficiency, E/m batt is the energy density and L/D is the lift-todrag ratio.For this demonstration, the range, battery energy density, and lift-to-drag ratio are all fixed.The propulsive efficiency η p is a product of the efficiency of the inverter η inv , motor η m and rotor η r ; we assume constant inverter and rotor efficiency.
We break down the aircraft into three main components and represent the total aircraft mass as where the airframe mass m air f rame refers to the remaining parts of the aircraft.The motor design has two significant roles in this approach.First, the geometric design of the motor directly alters the motor mass m motor .Second, the battery mass is influenced by motor efficiency.By rearranging (32) and converting the weight on the left-hand side to mass, the battery mass can be represented as where C is the parameter defined in (32).The dependency of motor efficiency on the battery mass is seen through the overall propulsive efficiency η p defined in (32).

Optimization problem
The goal is to optimize a baseline motor design with a comprehensive set of motor design properties; this includes geometric and operational variables and constraints.The objective function in the optimization problem is the total aircraft mass m total .The problem is summarized in Table 3.The geometric variables, shown relative to the motor geometry in Fig. 11, span all aspects of the motor cross-section.The motor geometry determines the mass and the size contributes to the torque equality constraint.The current and motor length are also introduced as design variables; these variables have a significant role in the motor operation.For a fixed geometry, the current has inverse effects on the output torque and efficiency.Thus, the current is a key design variable in this optimization problem because there is a need to match the input torque constraint whilst maintaining high efficiency for minimal battery mass.The motor length dictates another opposing relationship between torque and mass.As the motor Subject to τ l = τ r Torque balance 1 We aim to minimize total aircraft mass with respect to geometric and operational design variables, subject to a torque balance and voltage limit length increases, both the torque and motor mass increase; this is beneficial to ensure that torque equality is satisfied, but negatively impacts the overall objective.
A parameter sweep is conducted to optimize the motor across a set of design conditions.The airframe mass is selected as the sweep parameter for two reasons.First, the battery mass scales with the airframe mass, as seen in (34).Second, the airframe consists of most of the aircraft mass, in comparison to the motor and battery; thus, the operating RPM and torque necessary for steady-level flight are mostly dictated by the airframe mass.These values are calculated as a pre-processing step before optimization using a blade element momentum-based ideal loading design approach from Ruh and Hwang (2021); given a required thrust, this method determines the torque, RPM, and blade shape that minimizes aerodynamic losses.Under the assumptions of steady-level flight and a constant lift-to-drag ratio, the required thrust for each case in the parameter sweep can be calculated based on total aircraft weight.Although the total aircraft mass during optimization will change, the magnitude of the changes in battery and motor mass is small compared to the total aircraft mass; the expected variations in torque and RPM are small and negligible.Therefore, it is sensible to assume a constant torque and RPM for each case in the parameter sweep.We use a gear ratio of four.

Comprehensive versus baseline optimization comparison
The parameter sweep was conducted under a set of six different airframe masses, to optimize the motor design by minimizing total aircraft mass under a variety of operating conditions.We performed a set of baseline and comprehensive design optimizations for each case in the parameter sweep.The comprehensive motor design optimization includes the detailed design of the motor cross-section.The baseline case neglects changes in cross-sectional motor shape and only considers current and motor length as design variables.Thus, the only variable affecting motor mass is the motor length, a variable that also scales proportionally with the torque.This allows us to assess the influence of the motor geometry on performance.Table 4 summarizes the parameter sweep inputs and the results for the comprehensive design optimization, and Fig. 12 compares the baseline and comprehensive shape optimization results.
The comprehensive design optimization cases produce an overall improved design compared to the baseline optimization case; the benefits are amplified as the airframe mass increases.Minimal improvements are seen in the total aircraft mass; however, this is not unexpected, given that the airframe mass is constant and significantly larger than the other mass components.The battery mass exhibits larger differences as airframe mass increases, primarily due to the improvement in efficiency.Referring back to (34), the battery mass depends on the efficiency, motor mass, and airframe mass.In general, the motor mass is significantly smaller than the airframe mass.As a result, the changes in battery mass are influenced mainly by the motor efficiency.
We can formulate a more informative comparison by analyzing the influence of motor geometry.In general, torque is proportional to volume; more specifically, torque  12 A comparison of the baseline and comprehensive design optimization shows an improvement in the overall design and results when geometric optimization is taken into account scales with the geometry as such: where D and L represent the motor diameter and length, respectively.The freedom for both the motor geometry and length to change allows the diameter to shrink to reduce mass, whilst placing more emphasis on the motor length to meet the torque constraint.This can also be leveraged to produce torque from the motor length rather than the current amplitude, reducing copper losses and increasing efficiency.This explains the improvement in efficiency between the baseline and shape optimization results in Fig. 12.

Comparison of optimization results across parameter sweep cases
We can also draw comparisons between different cases in the parameter sweep.Each case operates under a different torque and RPM, so different designs are expected across the parameter sweep.The aircraft and battery mass exhibit similar trends; referring back to Fig. 12, these two outputs display a monotonically increasing trend with the airframe mass.The battery mass is driven by the airframe mass and efficiency; the latter is heavily determined by the current amplitude, the driving force in attaining the torque balance constraint.Copper losses scale quadratically with the current, so a decrease in efficiency is expected with an increase in required output torque, and thus the airframe mass.
Characterizing the variation of the optimized motor design across the parameter sweep is not as straightforward.The motor mass exhibits a different trend relative to the other components.The motor mass remains relatively constant for the lower airframe masses and begins to increase halfway through the parameter sweep.The larger airframe mass assumes a larger torque, so an increase in motor mass is expected to produce additional torque without compromising the efficiency.Increases in current Fig. 13 The variation of a subset of design variables is visualized in the parallel coordinates plot.Each case in the parameter sweep converges to a different optimized design, and the magnitude of variation based on airframe mass can be seen here optimization case; this is due to the modeling limitations and the physics within the methodology, discussed further in Sect. 5.
The optimized geometries relative to the initial motor dimensions are shown in Fig. 14.Each parameter sweep case began with identical motor dimensions; this allows for a common reference when comparing optimized designs.We see a reduction in overall size for each parameter sweep case; this is expected considering the objective of minimizing overall aircraft mass.For the high torque cases, additional changes in the optimal layout are seen around the magnets and stator slots, as the geometry tries to compensate for the higher current needed to satisfy the higher load.
To more clearly identify the geometric differences across the parameter sweep, Fig. 15 shows a one-sixth slice of each motor with the geometric design variabes overlaid onto the geometry.Each slice corresponds to a parameter sweep case from Table 5, identified by the large enclosed number near the inner radius.In the latter parameter sweep cases (four to six), we start to see significant geometric changes as a result of the increased torque demand.The motor structure is manipulated to simultaneously produce torque and minimize the core and copper losses coming from the high current necessary to produce satisfactory torque and maximize efficiency for the battery.
The limitations of the current methodology are centered around the tight bounds on the geometric design variables, motivated by the difficulties of our mesh warping technique.Mesh warping is prone to generating highly skewed mesh elements, which can lead to poor solver performance; an example of this is shown in Fig. 16.Under a constant mesh parameterization, the elements around the magnet and air gap become highly skewed from large changes in the magnet width.This is a consequence of utilizing a constant mesh parameterization throughout the optimization.Although re-meshing and reparameterization of the mesh during optimization would solve this Fig. 14 A comparison of the original and optimized motor geometries for each case shows a general reduction in size.The latter cases (4-6) differ more due to the higher torque demand, resulting in the further detailed design of the stator and magnet dimensions issue, these are discrete operations and are not feasible for gradient-based optimization.As a result, tight geometric design variable bounds have been set to avoid this issue.A caveat to this decision is the freedom of the motor geometry to change is limited, and resultant optimized motor designs are falsely constrained by design variable bounds rather than physics.Figure 14 shows side-by-side comparisons of initial and optimized geometries for each parameter sweep case.Although the internal motor structure shows variation, many of the geometric variables approach the same value due to the tight bounds.An example of this is the shaft diameter, which is dictated by the maximum shear stress of the shaft material; however, the shaft diameter bound was reduced to prevent significant skewing of the rotor mesh elements.An alternative approach is necessary to allow for larger geometric changes without compromising the mesh quality.

Conclusion
In this paper, we have presented an adjoint-based motor design optimization methodology for high-dimensional design problems.This methodology is motivated by the poor scaling of computational time with respect to the number of design variables, a roadblock prevalent in existing motor design approaches.Free-form deformation is incorporated to parametrize the motor geometry in a differentiable manner.A PDEbased mesh warping algorithm is utilized to manipulate the mesh; the mesh boundary  5.Some variables like t slot experience significant changes between cases, while others like r sha f t do not Fig.16 Large geometric changes can lead to highly skewed mesh elements.The significant width reduction of the magnet causes skewing in the mesh elements and makes the solvers more prone to divergence conditions are applied using fractional load-stepping, and the updated mesh is found by solving a sequence of nonlinear problems with small deformations.The motor analysis is conducted using high-fidelity finite element analysis by solving Maxwell's equations.The resultant output field is utilized to calculate high-level parameters, such as power losses, torque, and efficiency.Derivatives are computed analytically for each method.As a result, even though model evaluation is expensive, the derivative computation cost is reduced using the adjoint method, and updating the motor design is efficient.
The methodology was demonstrated using an aircraft system-level design problem, in which it searches for a motor design that minimizes aircraft mass across a series of operating conditions.A baseline optimization case with a fixed motor shape is analyzed for comparison.Results show marginal improvements in aircraft mass for each case in the parameter sweep; this is due to the small mass ratio between the motor and battery relative to the entire aircraft.Improvements of up to 35% in motor mass and 3% in motor efficiency are found.The results indicate the feasibility and success of an adjoint-based motor design optimization approach; we see improvements in the overall design, and the shape optimization results show that efficient exploration of the full design space is possible using adjoint-based optimization.Future development of this work has two distinct avenues: improving the modeling of the physics and addressing solver ill-conditioning due to large geometric deformations.
The first area of future work addresses the shortcomings of the motor analysis model.The torque error shown in Table 1 and flux field qualitative differences in Fig. 9 can be attributed to the permeability model and torque computation approach.Alteration to the permeability model is required, as the curve fit visualized in Fig. 2 deviates from the data in the flux-density-magnitude range expected within the motor air-gap.In addition, we assume the direction of the air-gap flux density is purely radial.A higher fidelity method, such as virtual work or the Maxwell stress tensor, is more suitable to calculate torque using the FEA field solution.The current methodology also does not utilize periodicity, a key property of electric motor operation.Periodicity is incorporated by modeling a fraction of the motor geometry swept by one pole pair.This concept is commonly used for motor analysis (Ahn et al. 2015;Yamazaki and Ishigami 2010;Pellegrino and Cupertino 2010;Lin et al. 2018) as a means to reduce computational expense.In addition, coupled multiphysics models are required to accurately capture the influences from other disciplines.This includes thermal constraints to model demagnetization and temperature dependence of power loss, and structural considerations related to harmonics and maximum stresses.Feedback coupling of the temperature distribution is necessary to capture the proper thermal effects in the electromagnetic analysis (Babcock et al. 2023).
The second major area of future work addresses the issue of solver ill-conditioning due to significant skewing of mesh elements.In the current work, small geometric design variable bounds are applied to avoid this issue.Improvements to the current mesh adjustment approach are required to allow for larger geometric changes in the motor design.The current method suffers from severe mesh element warping due to the constant mesh parameterization.To adjust the mesh information during optimization, we can integrate an outer-loop geometry refinement step to reset the mesh and FFD parameterization.The outer loop can be automatically called based on mesh quality Fig. 17 The outer-loop geometry reparameterization would be integrated in the manner shown here.The geometric design variables δ g update the base geometry and regenerate the parameterization, and the optimization continues from the updated geometry metrics like skewness and aspect ratio; during this step, the mesh of the new motor geometry is regenerated with the most recent geometric deltas.The inner optimization loop is then restarted and initialized with the design variables of the previous optimization run.This process is shown visually in Fig. 17.
This concept has been utilized by He et al. (2019) to conduct aerodynamic shape optimization, showing a circle transforming into an airfoil; the FFD parameterization is adaptively updated by exiting the optimization loop to adjust the FFD control points using an adaptation metric.To expand the geometric freedom of the shape optimization approach, this is a promising avenue for future work.

Fig. 2
Fig. 2 A set of smooth functions is fit to discrete permeability data.A joiner function (in green) is used to connect the linear and exponential sections.The joining function is chosen to be a cubic function to allow for matching values of permeability and derivatives at the section bounds.(Color figure online)

Fig. 9
Fig. 9 Magnetic flux density output fields from FEniCSx and Ansys Maxwell for no-load and 280 A current cases are shown here.A qualitative comparison shows agreement in overall field trends, with numerical discrepancies at higher load cases.Despite qualitative agreement between the 280 A cases between FEniCSx and Ansys Maxwell, we see a 14% error in the output torque

Fig. 11 A
Fig. 11 A visual representation of the geometric design variables relative to the motor geometry is shown here.The design variables represent changes in the geometric aspect

Fig. 15
Fig. 15We compare the optimized design from each parameter sweep case with the geometric design variables overlaid on the geometry.The large numbers within the circles correspond to the parameter sweep case indices shown in Table5.Some variables like t slot experience significant changes between cases, while others like r sha f t do not

Table 1 A
torque comparison between Ansys Maxwell and FEniCSx solvers shows increasing error with the applied current load

Table 3
We have described the optimization problem here

Table 4
A comprehensive summary of the otor design optimization parameter sweep inputs and outputs is shown here