© H. Wang et al., Published by EDP Sciences, 2024

1 Introduction

In recent years, to accelerate the goal of achieving low-carbon and green aviation, countries have become increasingly eager to explore more/all-electric aircraft, hoping to develop more/all-electric aircraft with at least equivalent performance to traditional aircraft as soon as possible. The demand for artificial intelligence, new concept airborne equipment, and advanced power electronic technology in more/all-electric aircraft is constantly increasing, and they are developing towards high altitude, high speed, long endurance, and lightweight [1]. The required power of equipment carried by more/all-electric aircraft is increasing, and the thermal load of the aircraft is also growing. In addition, the heat dissipation limitations brought about by supersonic cruising result in an increase in the stagnation temperature of the external environment, and the heat sinks, such as fuel, are compressed [1]. These have put forward higher requirements for the thermal design of airborne power supply equipment such as the permanent magnet synchronous generator (PMSG). The PMSG is often used as an aircraft power supply due to its small size, high efficiency, and high power factor. In the trend of more/all electrified aircraft, the thermal design of PMSG is crucial for the overall thermal design of the aircraft system architecture.

Overheating of generator components can worsen the performance of high-efficiency permanent magnets, accelerate winding insulation aging, and decrease efficiency, shortening the service life of the PMSG [2]. Therefore, the close connection between the electromagnetic and thermal design of the PMSG cannot be ignored [3]. Suppose the electromagnetic system and mechanical components of the PMSG can be properly designed in the early stages of design, and the key parameters of the PMSG design can be adjusted based on the evaluation of thermal losses. In that case, it can significantly reduce the aircraft design time and material costs. The thermal design of the PMSG includes reducing heat generation and designing heat dissipation methods to manage the heat generated. The heat generation of PMSG is determined by the losses, and higher losses lead to additional heat. The heat dissipation of a PMSG is achieved through thermal convection, which transfers heat from the solid surface to the moving fluid and then to the heat sink components for fluid cooling. Therefore, in the early design stage, estimating the losses of the PMSG and properly designing the geometric structure and electromagnetic characteristics of the PMSG to reduce the heat accumulation is necessary. Simultaneously, suitable heat transfer media should be selected to implement the heat dissipation components.

Currently, many researchers have researched the thermal design of motors. Cai et al. proposed a technique for the thermal modeling of flux-switched permanent magnet (FSPM) motors that considers anisotropic thermal conductivity, contact thermal resistance, and convective heat transfer in air gaps but only studied from the perspective of the thermal discipline [2]. Nacu et al. only conducted numerical finite element simulation on the thermal system of a low-speed external rotor PMSG for wind power generation, and they have not been studied in conjunction with other systems in an integrated manner [3]. Tatevosyan et al. only established a thermal balance equation based on the magnetic system of the motor and conducted thermal calculations, which were verified in ANSYS software [4]. Iden et al. proposed a model that combines the thermal characteristics of aviation generators with electromagnetic-mechanical features and then conducted multi-disciplinary design optimization (MDO), but the optimization improvement was not significant [5]. Considering that the thermal design of the PMSG involves the interaction of thermal, electromagnetic, and mechanical and has strong discipline coupling, the traditional conceptual design method is less effective [3]. Therefore, the PMSG thermal design needs a multi-disciplinary conceptual design method with high coupling. However, the information required for multi-disciplinary system-level conceptual design has dramatically increased, posing significant challenges to multi-disciplinary motor design. Fortunately, multi-discipline optimization(MDO) has been extended to these system disciplines, which can ensure the interconnection of systems during conceptual optimization design and avoid unforeseen defects in the results of conceptual design in the later design stage [5,6].

This paper proposes an optimization design method for aviation PMSG based on MDO. Firstly, based on the coupling relationship between the thermal and electromagnetic mechanical performance of PMSGs, an MDO optimization model is established. Then, a multi-objective particle swarm optimization (MOPSO) based on a non-dominated solution is designed with the optimization goal of high efficiency and better heat dissipation performance as the main objective, combined with the optimization requirements of small size and light weight, and the critical parameters related to thermal, electromagnetic and mechanical of the generator as the decision variables. The MOPSO algorithm is used to optimize critical parameters of PMSG thermal design, solve multi-objective optimization results, and support multi-disciplinary optimization of the PMSG thermal design. Finally, a conceptual design method of PMSG with high discipline coupling and physical fidelity based on thermal modeling is presented.

2 Establishment of PMSG multi-disciplinary optimization analysis model

Establishing the MDO model of a PMSG requires decomposing the generator model into electromagnetic, mechanical, and thermal parts, clarifying the direct or indirect relationship between critical parameters and heat, and establishing each system's mathematical model.

2.1 Mathematical modeling of the electromagnetic part

The electromagnetic design of a PMSG is mainly to determine the electromagnetic size of the motor, that is, the value of the line load A and the maximum air gap flux density B. The electromagnetic size of the PMSG can be expressed by

$$AB=\frac{6.1\times P\text{'}}{{\alpha }_p^{\text{'}}{K}_{nm}{k}_{dp}{D}^2{L}_{ef}n}$$(1)

where n represents the motor speed; ${\alpha }_p^{\text{'}}$ represents the calculated pole arc coefficient, ${\alpha }_p^{\text{'}}=B_{ev}/B$ and B_av represents the average magnetic density of the air gap; K_nm represents the waveform coefficient of the air gap magnetic field, and K_nm= 1.11 when the air gap magnetic field waveform is a sine waveform; k_dp represents the armature winding coefficient, and usually equals to the fundamental winding coefficient; the calculation formula of the line load A is shown in equation (2), where m represents the number of phases of the stator winding, I represents the phase current of the stator winding and N represents the number of turns per phase of the stator winding

$$A=\frac{2mNI}{\pi D}.$$(2)

P ^′ represents the calculated power of the PMSG, and its formula is shown in equation (3), in which K_v represents the ratio of the induced electromotive force to the terminal voltage of the PMSG with rated load, S represents the rated apparent power, and η represents the efficiency under rated load

$$P^{\prime }=\frac{K_{v}S}{\eta }.$$(3)

D represents the diameter of the stator core, L_ef  represents the effective length of the core, D²L_ef can approximately represent the volume of the effective part of the stator, and the volume of the effective part of the rotor is also included. Therefore, D²L_ef  represents the main size of the motor. The numerical relationship between the two can be determined by the main size ratio formula (4), where p represents the number of pole pairs.

$$\lambda =\frac{2p{L}_{ef}}{\pi D}.$$(4)

From the analysis of equation (1), it can be concluded that the product of line load A and the maximum air gap flux density B will affect the volume of the motor. The larger the product AB, the smaller the volume of the PMSG. This is favorable for the generator design of aircraft. When a larger line load A is selected, the volume of the motor decreases, and the material requirement and the manufacturing cost decrease. However, a large line load will reduce the flux per pole. To ensure that the induced electromotive force is constant, the number of turns of the stator winding needs to increase, thereby increasing the amount of copper wire and copper consumption [7–9].

Because the core loss of the PMSG is proportional to the square of the magnetic density, to improve the efficiency of the PMSG, the internal magnetic density distribution of the PMSG should be reasonable. The maximum air gap magnetic density B should not be too large. Meanwhile, excessive air gap magnetic density will also increase the magnetic circuit saturation. Therefore, in the design of the PMSG, a lower maximum air gap flux density B is usually selected, and the higher line load A is established under the condition that the copper consumption is within a reasonable range to reduce the thermal load of the PMSG. Therefore, there is a tradeoff between the electromagnetic design and the thermal design of the PMSG, and the value of the electromagnetic size needs to be determined in combination with the thermal design.

2.2 Mathematical modeling of mechanical parts

The mechanical design of a PMSG mainly includes determining the maximum centrifugal stress that the rotor can withstand and determining the critical shaft diameter of the rotor.

2.2.1 Determination of maximum centrifugal stress of rotor

When the PMSG rotates, the force analysis of the rotor is shown in Figure 1. dm represents a mass infinitesimal of the rotor. The centrifugal force F on the rotor surface is:

$$F=\frac{M\rho D_z^2{\omega }^2}{4}$$(5)

where M represents the cross-sectional area of the rotor, ρ represents the rotor density, Dz represents the outer diameter of the rotor, ω represents the rotor angular velocity. The centrifugal stress γ of the PMSG is

$$\gamma =\frac{F}{M}=\frac{\rho D_z^2{\omega }^2}{4}=\rho v^2$$(6)

where υ represents the linear velocity on the outer surface of the rotor. Centrifugal stress should meet the following restriction:

$$\gamma =\rho v^2\le \frac{\left[\gamma \right]}{Q}$$(7)

where [γ] represents the permissible stress of the rotor material, Q represents the safety factor of the material. According to equation (6), the maximum centrifugal force γ_max in the PMSG rotor is proportional to the rotor density ρ and the square of the linear velocity v on the outer surface of the rotor. The rotor density ρ is determined by the material used in the rotor. In the case of the angular velocity range, the linear rotor velocity v is determined by the outer rotor diameter D_z [10]. From the analysis of equation (1), it can be seen that the volume of the rotor, that is, its diameter and length, are related to the electromagnetic size. Therefore, it can be inferred that the maximum centrifugal stress of the rotor is also related to the electromagnetic size. Meanwhile, it is also restricted by the thermal effect of the PMSG.

Fig. 1 Force analysis of micro mass units of the rotor.

2.2.2 Determination of the critical shaft diameter of the rotor

In PMSG design, the torsional stiffness condition is usually used to determine the critical shaft diameter of the shaft. When calculating the stiffness of the rotor shaft, the short circuit or braking of the PMSG must be taken into account. Therefore, the appropriate safety factor K_a should be selected. In general, the safety factor K_a= 2–3. The critical shaft diameter calculation formula of the shaft is as follows:

$$D_{s\_ctit}=\sqrt[4]{\frac{32K_a{T}_n\times 180}{G{\pi }^2\left[\theta \right]}}$$(8)

where T_n represents the rated torque of the PMSG, T_n = 9550P/n; G represents the shear elastic modulus of the shaft, which is determined by the material properties of the shaft.

2.3 Mathematical modeling of the heat part

According to the analysis of the first section, the thermal design of the PMSG includes heat generation analysis and heat dissipation design.

2.3.1 Heat generation analysis

PMSG generates losses during energy conversion, which dissipate in the form of thermal energy and generate heat. Energy loss includes winding copper, core, eddy current, and wind friction loss.

1) Winding copper loss

Since the rotor of the PMSG is a permanent magnet, no excitation winding is required. Only the stator core has windings. Therefore, only the copper loss in the stator needs to be calculated. Copper losses are losses caused by current through the armature winding, including DC losses, also known as I²R losses, and additional eddy current losses [4].

(1) I²R losses

The following formula gives I²R losses

$$P_{I^{2}R}=m{I}^{2}R$$(9)

where m represents the phase number of the stator winding, I represents the phase current of the stator winding, and R represents the winding DC resistance. The winding resistance will change with the temperature change. The higher the temperature is, the greater the resistance is. Usually, the winding resistance at 75 °C is taken as the standard. The resistance calculation formula at this temperature is as follows:

$$R_{75}=\frac{{\rho }_{75}al}{b}$$(10)

where ρ₇₅ represents the resistivity of the armature winding at 75 °C, a represents the number of conductors in series per phase, l represents the length of each coil winding, and b represents the number of copper wires in each conductor.

(2) Additional eddy current loss

The additional eddy current loss is derived from the skin effect of the conductor itself caused by the influence of AC around the conductor when AC is out of the conductor itself, as well as the proximity effect produced between conductors that are close together at this time. The calculation formula is as follows:

$$P_{addv}=P_{I^{2}R}\left(k_d-1\right)$$(11)

where k_d represents the average resistance coefficient, which is the ratio of effective AC resistance to DC resistance.

Analysis of equations (9) and (11) shows that the copper loss is temperature dependent, and as the temperature increases, the I²R loss will increase while the additional eddy current loss will decrease. When the PMSG operating frequency is lower than 1000 Hz, the average resistance coefficient is close to 1, and the additional eddy loss can be neglected. The additional eddy current loss must be considered only when the PMSG operating frequency is higher than 5000 Hz. Moreover, I²R loss and line load A are related to the stator winding. While I²R loss needs to be reduced, the line load A needs to be as large as possible; thus, the design requirements of these two are just opposite, and the coupling between these two must be considered for the PMSG design.

 2) Core loss

The stator core generates iron losses in an alternating magnetic field due to magnetization reversal and eddy currents heating. The iron losses are generally divided into hysteresis losses and eddy current losses [11].

(1)Hysteresis loss

The frequency of the alternating magnetic field and the flux density determine the hysteresis loss, and the hysteresis loss per unit volume of the iron core can be calculated as follows:

$$P_{mh}=f{K}_f{B}_{fe}^{\alpha }$$(12)

where f represents the winding current frequency, K_f represents the material constant, B_fe represents the maximum density of core flux, and the value α is related to the magnetic core density, ranging from 1.8 to 2.2.

(2) Eddy current loss

Eddy current loss is caused by the current induced by alternating magnetization in the core, and the following formula calculates the eddy current loss per unit volume:

$$P_{ec}=\frac{{\pi }^2{f}^2{t}^2{B}_{fe}^2}{{\rho }_z\beta }$$(13)

where t represents the material thickness, β represents the geometric coefficient, and ρ_z represents the resistivity of the material. The stator core is usually designed with thin silicon steel sheets stacked and insulated between each other, with eddy currents flow within each sheet to reduce eddy current losses.

The rotor flux density B_fe is determined by the flux Φ per pole and the corresponding magnetic circuit area M_fe. The calculation formula is as follows, where F_fe represents the wave amplitude factor reflecting the effect of magnetic circuit saturation:

$$B_{fe}=F_{fe}\frac{\Phi }{M_{fe}}.$$(14)

The following equation determines the magnitude of the magnetic flux per pole:

$$\Phi =B{\alpha }_p^{\text{'}}\tau L_{ef}.$$(15)

The rotor flux density varies with air gap flux density. A smaller air gap magnetic density in electromagnetic design helps reduce the core loss.

3) Rotor eddy current loss

The permanent magnet eddy current losses arise from the higher harmonics of the air gap flux, which is calculated as shown below [12]:

$$P_{zec}=\frac{1}{T_e}{\int }_0^{T_e}{\displaystyle \sum _{i=1}^k}{J}_i^2\Delta J_i$$(16)

where J_i represents the surface eddy current density, Δ represents the rotor surface unit area, T_e represents the rotor surface area, i =1∼k represents the number of units of rotor eddy current density. Usually, the rotor eddy current loss is negligible and only considered when the PMSG is running at high speed.

4) Wind friction loss

The wind friction loss is the energy lost in the gas due to the relative motion of the gas flowing between the rotor and the stator, and its contribution to heat at high speeds cannot be ignored [13], which can be calculated as follows:

$$P_{wind}=c{K}_m\pi {\rho }_f{\omega }^3{D}_z^4{L}_{ef}$$(17)

where c represents the rotor surface roughness, K_m represents the rotor surface friction coefficient, and ρ_f represents the fluid density.

5) Calculation of temperature rise due to losses

According to the principle of heat transfer, it is known that the relationship between the maximum temperature rise T of the PMSG and the copper loss, iron loss, and wind friction loss is [14]:

$$\begin{array}{l}{P}_{I^{2}R}+P_{addv}+{\lambda }_1(P_{mh}+P_{ec}+P_{zec})+{\lambda }_2{P}_{wind}\\ =\pi \sigma D_z{L}_{ef}\left({\lambda }_3+\frac{{\lambda }_4}{{\lambda }_s}\right)\Delta T_{emp}\end{array}$$(18)

where λ₁and λ₂ respectively represent the weighting factor of the effect of iron loss and wind friction loss on the PMSG temperature rise; σ represents the heat dissipation coefficient of the heat dissipation surface; λ₃ and λ₄ respectively represent the weighting factor considering the contribution of the rotor and stator surfaces to heat dissipation; λ_s represents the ratio of the outer diameter of the rotor and stator, generally known as the split ratio [15].

2.3.2 Heat dissipation design

Standard cooling methods include natural air cooling, forced air cooling, liquid cooling, etc. Different cooling methods are chosen for different cooling requirements of different equipments. Liquid cooling is usually selected for high power and high heat density equipment such as aviation PMSG. Because the heat transfer coefficient of liquid cooling (Fig. 2) [16] is much higher than the other two heat dissipation methods, the operating performance is preferable.

Liquid cooling often uses water, oil, or glycol aqueous solution as the liquid medium because oil can be directly poured into the interior of the PMSG and the rotor and stator windings of the PMSG for more efficient heat exchange. Oil's heat dissipation performance is better, the heat dissipation efficiency is higher, and the aircraft fuel tank is in a low-temperature environment. Thus oil can quickly dissipate heat to speed up the cooling cycle without adding extra weight, so aircraft generators usually use oil-cooled heat dissipation.

According to the motor heat transfer theory, heat will be transferred from high- to low-temperature place. Assuming that the object's state in the heat transfer process is constant and the heat transfer process is linear, the heat transfer process can be simplified to one-dimensional, steady-state heat transfer model. Then the temperature transformation of the PMSG and coolant in the heat dissipation process can be obtained as follows:

$$C_g{G}_g(T_e-T_r)={\rho }_o{C}_o{q}_{vo}{Q}_{in}{\frac{d{T}_e}{dt}}_{emp}$$(19)

where G_g represents the mass of PMSG, C_g represents the specific heat capacity of PMSG, Q_in represents the heat absorbed by the coolant, ρ_o represents the density of cooling oil, C_o represents the specific heat capacity of cooling oil, q_vo represents the volume flow of cooling oil, T_e represents the initial temperature of the generator, T_r represents the temperature of the generator after heat dissipation, assuming T_e − T_r = dT_cold, so that the maximum temperature rise of the PMSG can be changed to:

$$\Delta T_{temp}^{*}=\Delta T_{temp}-\Delta T_{cold}.$$(20)

The volume of the cooling oil can be determined according to equation (20). It is also necessary to determine the cooling and heat dissipation area of the PMSG:

$$M_{rg}=\frac{Q_{in}}{K_r\cdot \Delta T_r}$$(21)

where K_r represents the radiator heat dissipation coefficient and ΔT_r represents the average temperature difference between the cooling oil and the ambient air. From equations (19) to (21), the relevant parameters required for heat dissipation design can be determined.

Fig. 2 Comparison of heat transfer coefficients of different cooling methods [16].

3 Design of MOPSO algorithm for generator

The optimization algorithm is used to solve the critical parameters in the multi-disciplinary optimization model of the aeronautical PMSG. Since the conceptual design of the PMSG involves optimizing thermal, electromagnetic, and mechanical systems, the multi-objective optimization of the generator must be a multi-objective problem. It can be seen from the analysis in Section 2 that there are certain contradictions between these optimization objectives, so it is necessary to choose an appropriate optimization algorithm to coordinate multiple objectives and achieve the overall optimization.

The multi-objective optimization algorithm can be divided into normalization algorithm and non-normalization algorithm. The normalization algorithm transforms the multi-objective problem into a single objective problem by the weighting methods and then solves it by the traditional optimization algorithm. However, it is easy to be affected by the subjective judgment of designers, and the optimization effect and efficiency are not good. Non-normalization algorithm uses the Pareto mechanism to deal with multiple objectives of the optimization technology directly and has the advantages of easy operation, simplicity, and versatility. It can get a whole Pareto frontier at one time, rather than the traditional optimization algorithm to solve one by one, improving optimization efficiency. Particle swarm optimization (PSO) algorithm, as a typical non-normalized algorithm, has fewer parameters than other algorithms, solid global searchability and search efficiency, and fast computing speed. It is good at dealing with global optimization problems of multi-extremum functions, which is suitable for solving the designed multi-disciplinary optimization model.

3.1 The basic principles of PSO

PSO initializes a group of random particles. At each iteration, the particle updates itself by tracking two “extreme values”. One is the individual extreme value p_B, which is the optimal solution searched by the particle itself; the other is the global extreme value g_B, which is the optimal solution searched by the whole population at the current moment. After two extreme values are searched, the particle updates its speed and position according to the following update formula:

$$\begin{array}{l}{v}_{p_i}^{k+1}=w{v}_{p_i}^k+c_{1}rand\left(0\text{and}1\right)\left(p_{B_i}^k-x_{p_i}^k\right)\\ +c_{2}rand\left(0and1\right)\left(g_{B_i}^k-x_{p_i}^k\right)v_p\end{array}.$$(22)

$$x_{p_i}^{k+1}=x_{p_i}^k+v_{p_i}^{k\text{+1}}.$$(23)

$$w=w_{ini}*w_{damp}.$$(24)

In equations (22)–(24), v_p represents the velocity of particles, w represents the inertia weight, w_ini represents the initial weight, w_damp represents the weight attenuation coefficient. The earlier weight is generally larger, and the algorithm tends to globally search and determine the range of the optimal value; the later weight is smaller, and the algorithm tends to the optimal range of local search and determines the Pareto optimal solution. Then, x_p represents the position of the current particle p_B and g_b represent extreme individual value and extreme global value, respectively; rand () is a random number between (0,1), c₁, c₂ represent acceleration factors, representing the acceleration weight of particles advancing towards g_Best and p_Best, whose value ranges from 0 to 4, usually c₁= c₂ =2, i =1,2,3 ... pop, pop represents population size; k represents the current iteration number.

3.2 Non-dominant solution

Since it is difficult to ensure an optimal set of solutions to a multi-objective problem, but there will be a set of solutions that are better than other solutions in all directions, this group of solutions is the non-dominant solution, also called Pareto optimal solution. A non-dominant solution is defined as a solution, such as r1, that is not dominated by any other solution in a set. This means that for any other solution, r2, there is at least one target where r1 is not superior to r2. In other words, there is no solution that outperforms r1 in all targets simultaneously. If no other solution in the set dominates r1, then r1 is considered non-dominant [17]. The non-dominated solution can accurately evaluate the optimality of the solution of the multi-objective problem, so this paper will design the multi-objective PSO algorithm based on the non-dominated solution.

3.3 Fitness function design

According to the design content of electromagnetic, thermal, and mechanical parts of the generator, with low temperature rise, low volume, and large electromagnetic size as optimization objectives, and with each one-dimensional particle as the decision variable shown in Table 1, the fitness function of multi-objective optimization of the generator is shown in the following equations:

$$\begin{array}{c}fit=\min (-fi{t}_D;fi{t}_V;fi{t}_T)\\ -fi{t}_D=-AB\\ fi{t}_V=D^2{L}_{ef}\text{+}{G}_g\\ fi{t}_T=T_{emp}^{*}\end{array}$$(25)

where fit_T fit_V,-fit_D respectively, represents the low-temperature rise, low volume, and large electromagnetic size optimization goals.

At the same time, there are many equality and inequality relations in the calculation process of the above objective, which will become the constraints limiting the values of decision variables. Therefore, there are many constraints in the multi-objective optimization problem. The direct or penalty function methods can solve the constrained multi-objective problem [18]. The advantage of the direct method is simple, but it requires a large amount of calculation, and the convergence rate is slow, so it is unsuitable for the PMSG model [19]. The penalty function method can transform constrained problems into unconstrained problems, and the convergence speed is fast, which can realize efficient optimization design.

The penalty function method transforms the constraint conditions into penalty function P, which is added to the objective function. For the inequality constraint g_c ≤ 0, the penalty term G_p is:

$$G_p=\max (0,-g_c)$$(26)

For the equality constraint h_c, the penalty term H_p is:

$$H_p=\max (0,\left|h_c\right|-\xi )$$(27)

ξ represents the minimum deviation of equality constraint, |h_c| − ξ ≤ 0.

The penalty function P is the sum of each constraint penalty term:

$$P=\vartheta \left[{\displaystyle \sum _{k_g=1}^{l_g}}{G}_p^{k_g}+{\displaystyle \sum _{k_h=1}^{l_h}}{H}_p^{k_h}\right]$$(28)

where ϑ represents penalty factor, k_g= 1,2,3...l_g, l_g represents the number of inequality constraints, k_h = 1,2,3...l_h, l_h represents the number of equality constraints.

Considering the different importance of different constraints, the weight of corresponding importance should be given to penalty terms. For inequality constraints, its weight $W_g^{k_g}$ is treated as follows:

$$\begin{array}{c}{W}_g^{k_g}=\frac{{\displaystyle \sum _{i=1}^{pop}}{G}_p^{k_g}\left(x_{p_i}\right)}{{\displaystyle \sum _{k_g=1}^{l_g}}{\displaystyle \sum _{i=1}^{pop}}{G}_p^{k_g}\left(x_{p_i}\right)}\\ {\displaystyle \sum _{k_g=1}^{l_g}}{W}_g^{k_g}=1\end{array}.$$(29)

For equality constraints, its weight $W_g^{k_g}$ is treated as follows:

$$\begin{array}{c}{W}_h^{k_h}=\frac{{\displaystyle \sum _{i=1}^{pop}}{H}_p^{k_h}\left(x_{p_i}\right)}{{\displaystyle \sum _{k_h=1}^{l_h}}{\displaystyle \sum _{i=1}^{pop}}{H}_p^{k_h}\left(x_{p_i}\right)}\\ {\displaystyle \sum _{k_h=1}^{l_h}}{W}_h^{k_h}=1\end{array}.$$(30)

The numerator in equations (29) and (30) represents the sum of the penalty terms of all particles under this constraint; the denominator represents the sum of penalty terms of all particles under all inequality or equality constraints. The final penalty function P is as follows:

$P=\vartheta \left[{\displaystyle \sum _{k_g=1}^{l_g}}{W}_g^{k_g}{G}_p^{k_g}+{\displaystyle \sum _{k_h=1}^{l_h}}{W}_h^{k_h}{H}_p^{k_h}\right].$(31)

The objective function can be rewritten as:

$$fit=\min (-fi{t}_D;fi{t}_V;fi{t}_T;P).$$(32)

3.4 The design process of the MOPSO algorithm

The flowchart of the MOPSO algorithm is shown in Figure 3, the specific process is as follows:

• Initialize the MOPSO, set the number of iterations T, population size pop, population size rpop of the non-dominant solution, particle dimension D, the range of particles in each dimension [popmin(j), popmax(j)], and the velocity range of particles [Vmin(j), Vmax(j)], where j represents particle dimension.

• Initialize the speed and position of each particle, and initialize the current position as the historically optimal position of the population.

• Calculate the three fitness values of each particle fit_D[x_p], fit_v[x_p], fit_T[x_p], and penalty function value as the initial value of algorithm iteration, and preliminarily screen the non-dominant solution set rep.

• Select any particle position from the non-dominated solution set as the optimal position for the population. Update the velocity of each particle according to the velocity update formula. Then the position of the particle is updated according to the position formula.

• Calculate each particle's fitness and penalty function values.

• Judge whether the particle is dominated; otherwise, add the non-dominated solution set rep. Then judge whether there is a dominant solution in the non-dominated solution set, and delete it if there is. Finally, judge whether the non-dominated solution reaches the maximum population number of non-dominated solutions, and delete the excess particles if it exceeds.

• Update weight coefficients. Judge whether the algorithm's operation reaches the maximum number of iterations; otherwise, return to step (2) to continue iteration. If yes, the algorithm's process is terminated, and several groups of non-dominated solution sets and corresponding fitness function values are output.

Fig. 3 Flowchart of MOPSO algorithm.

4 Optimization results of the MOPSO algorithm

The non-dominated population rpop of the MOPSO algorithm is set to pop/2, and the maximum number of iterations T is set to 400; the initial weight coefficient w_ini = 0.7, and the weight attenuation coefficient w_damp = 0.98; the acceleration factors c₁ = 1.4995 and c₂ = 1.4995. Before the conceptual optimization design of the PMSG, it is necessary to confirm the number of populations that the MOPSO algorithm can achieve global convergence. The Pareto frontier is useful to determine a suitable particle population size pop when the maximum number of iterations T is selected, then the MOPSO algorithm can achieve global convergence accurately and rapidly. The Pareto frontier curves under different population sizes are shown in Figure 4. In the figure, V represents the volume of the PMSG; G represents the weight of the PMSG; T_r represents the temperature rise of the PMSG; and AB represents the electromagnetic size. A comparison of the Pareto frontier curves for the population size of 200, 300, and 400 shows that when the population size increases from 200 to 300, the position of the Pareto frontier curve changes greatly, indicating that it is only locally convergent at this time. When the population size varies from 300 to 400, the position of the Pareto frontier curve is almost unchanged. Comparative results show that the algorithm can achieve global convergence when the population size is over 300. Considering that an increase in population size will slow down the calculation speed, a population size of 300 was ultimately chosen.

As shown in Figure 4, the volume and weight of PMSGs typically increase as the temperature rise increases. In contrast, the electromagnetic size decreases as the temperature rise increases. Considering that the conceptual design method is based on thermal modeling, a small portion of the air gap magnetic density will be sacrificed during the algorithm optimization process to reduce the heat source and requirements for the radiator.

The previous analysis shows that the comprehensive optimality of each solution on the Pareto curve is the same for the three objectives. Therefore, any data set is selected to compare with the PMSG performance parameters before optimization and calculate the corresponding rate of change. Following the MOPSO algorithm in Figure 3 based on the electromagnetic, mechanical, and thermal models in Section 2, the critical parameters optimized by MOPSO in this group are shown in Table 2. The optimized performance parameters obtained based on these critical parameters are shown in Table 3. The selected PMSG-rated parameters for comparison are a rated output power of 500kW, rotate speed of 25000 r/min, a current of 722.53 A, and a voltage of 400 V.

From the comparison of the data in the table, it can be seen that the optimized PMSG efficiency has significantly improved, and the temperature rise has also significantly decreased. The volume radiator area of the PMSG and the volume flow rate of oil have slightly increased, but the change is tiny. It can be seen that under similar PMSG volume conditions, the thermal performance of the PMSG has been greatly improved.

Although the electromagnetic size of the PMSG has decreased to a certain extent compared to the original, which will lead to an increase in volume and weight, it can reduce copper and iron losses, which is beneficial for improving thermal performance. The permissible stress of the mechanical design parameters of the PMSG increases with the increase in volume, and the requirements for the mechanical endurance of the material may increase. The critical diameter of the shaft changes slightly only because the autogenous material and torque limit it, but not because of the improved thermal properties.

Fig. 4 Comparison of Pareto frontier curves under different population numbers: (a) electromagnetic size and temperature rise and (b) volume + weight and temperature rise.

5 Conclusion

This paper establishes a multi-disciplinary optimization analysis model for the PMSG with coupled electromagnetic, mechanical, and thermal systems. An efficient and intelligent multi-objective optimization algorithm is presented, forming a conceptual design method for aviation PMSGs based on MDO. The convergence of the optimization algorithm is verified, and the optimization results are compared and analyzed. The following conclusions can be drawn:

• According to the needs of the optimal design of the PMSG, The critical parameters of the electromagnetic, mechanical, and thermal systems of the PMSG are determined. The coupling in the electromagnetic, mechanical, and thermal systems is then analyzed to derive the relationship between each critical parameter and thermal parameters. Finally, the MDO model of the motor is built.

• Combining non-dominated solutions and penalty functions, The MOPSO algorithm is designed based on the MDO model. The algorithm simultaneously considers the interactions between various parameters and disciplines, while balancing the conflict between PMSG thermal and other requirements. Ultimately, a reliable and coordinated conceptual design method for aviation PMSGs based on thermal modeling is formed.

• The proposed conceptual design method can improve the thermal performance of aviation PMSG. It effectively supports the multi-disciplinary optimization of PMSG thermal design and can reduce the design cost of aerospace generators and shorten the development cycle.

Funding

This research is from free exploration research.

Conflicts of interest

The authors declared no potential conflict of interest with respect to the research, authorship, and/or publication of this article.

Data availability statement

The authors confirm that the data supporting the findings of this study are available within the article.

Author contribution statement

He Wang conceived of the presented idea. Fengming Ai developed the theory and performed the computations. Linke He and Zhongzheng Zhou verified the analytical methods. Yangang Wang supervised the research.

References

All Tables

All Figures

	Fig. 1 Force analysis of micro mass units of the rotor.
In the text

	Fig. 2 Comparison of heat transfer coefficients of different cooling methods [16].
In the text

	Fig. 3 Flowchart of MOPSO algorithm.
In the text

	Fig. 4 Comparison of Pareto frontier curves under different population numbers: (a) electromagnetic size and temperature rise and (b) volume + weight and temperature rise.
In the text