Matlab-Simulink Coupling to Finite Element Software for Design and Analysis of Electrical Machines

1. Introduction

Classical electric machine design processes normally do not take into account some specific aspects (PWM voltage supply, mechanical resonances, etc) of the application where the motor is integrated. Actually the electromagnetic properties of the machine can influence the overall performances of the system, or the performances of the motor can get worse due to the influence of other components that are connected to the motor. For example, at low speeds, the torque ripple produced by the machine could cause undesirable speed pulsations and inaccuracies in motion control. At the same time, machine performance could be negatively affected by the Power Converter. For instance, the large amount of harmonic components contained in a PWM wave can increase the iron losses. Furthermore, when long connection cables are used between motor and power converter, high dv/dt voltages, which are characteristic in PWM signals, might cause over voltages at motor terminals which could damage seriously the winding insulation.

A possible solution to identify these problems is to add a new step into the design process previously to the prototyping stage, where the behavior of the machine into the final application is analyzed, evaluating the interaction between different elements of the system.

In order to achieve this purpose, this chapter deals with the integration of the numerical magnetic field computation software FLUX and the system simulator MATLAB-SIMULINK into only one simulation tool. The electric machine is modeled with FLUX software, whereas the control, the electric components and the mechanical systems are implemented in MATLAB-SIMULINK. In addition, the Finite Element Method (FEM) for electrical motor design and analysis is described. On the other hand, MATLAB-SIMULINK coupling to FEM software is explained and some clues related to this issue are given. And finally two case studies are presented in order to illustrate the benefits of this multi-domain simulation method. In both cases a surface mounted permanent magnet synchronous motor is analyzed.

2. Co-simulation tools for electrical machine design

Typically conventional electric machine design process is divided into two stages: an initial analytic design is complemented by a more accurate analysis based on FEM Software. A disadvantage of this methodology is that some specific aspects of the application for which the machine is designed are not entirely taken into account during the design process. Consequently, after integrating the machine into the final application, some performance problems may appear due to the interaction between the motor and different elements of the system. A new design step where these phenomena are identified by multi-domain analysis may help in the task of achieving an optimal design.

An electric drive is a non-linear multi-domain system which involves electric, electromagnetic, thermal and mechanical domains. The analysis of the whole system requires a specific simulation tool for each domain. But this way it is not possible to evaluate the interaction between different subsystems. A possible solution could be a simulation tool capable to integrate all domains. However, it is well known that the more are the domains of the simulation tool, the less is the accuracy. Therefore, a solution may be the integration of different simulation tools in a multi-domain analysis. Few years ago, that solution was impossible due to limitations in computation capacity of the computers and incompatibilities between different simulation software. Nevertheless, nowadays the development of new powerful work stations and updated simulation tools oriented to co-simulations make possible to analyze relatively complex multi-domain systems by means of co-simulations.

Figure 1.

Proposed Electric Machine Design Methodology

In this section some useful co-simulation examples for electric machine design will be presented.

2.1. Coupling between control unit – Electric machine – Mechanical system

Normally the electromagnetic torque produced by an electric machine is not entirely smooth but it presents a certain ripple component. This oscillatory component is due to a non-sinusoidal spatial distribution of flux density in the air-gap, and due to the effect of slots known as cogging torque.

Regarding to the mechanical system, many times it may present resonances which can be excited by the drive. In that case important level of vibrations and noise can be produced.

Carrying out co-simulations concerning the control unit, the electric motor FEM model and the mechanical system model it is possible to evaluate the effect of torque ripple over the mechanical system (see Fig 2-A). Furthermore, this multi-domain simulation makes possible a pre-evaluation of control algorithms, analyzing the influence of different control strategies as well as the influence of regulator parameters upon system performances.

Figure 2.

Examples of different multi-domain simulations

2.2. Coupling between power electronics converter – Electric machine

Few decades ago, technological advances in power converters made possible the development of new variable speed drives, leading to important performance improvements in many applications. In this kind of drives, the speed of the motor is controlled varying the frequency of the supply voltage by means of an power converter, which generates a well known Pulse Width Modulated (PWM) voltage.

The high harmonics content of this non-sinusoidal supplying voltage causes significant increase in electromagnetic losses, leading to a less efficient motor. These losses could be taken into account in the electromagnetic – thermal coupling in order to compute the overall heating of the machine.

Another phenomenon caused by the mutual influence between the power converter and the machine is the overvoltage at machine terminals. This overvoltage is caused by the characteristic high dv/dt of PWM waveforms and it may appear when connection wires between power converter and machine are too long. In the worst case, this effect could irreversibly destroy the insulation of stator winding or even could damage some semiconductor of the converter.

Besides these two phenomena, a pulse width modulated voltage has another effect on machine performance. That kind of supplying voltage causes a non-sinusoidal current, which contributes to increase the ripple of the torque.

All these phenomena could be analyzed by multi-domain simulations concerning the electric and electromagnetic domains (see Fig 2-C). After identifying all these problems, several design modifications can be made in order to improve the performance of the electric drive.

2.3. Coupling between electromagnetic and thermal domains

It is widely known that the temperature is very critical parameter in the performance of electric machines. The overheating of the machine depends on various factors such as core losses, current amplitude and working cycles. Typically conservative thermal designs are applied which may lead to bigger and more expensive machines. In order to get more compact motors, it is necessary to optimize the thermal design according to the requirements of the application.

This electromagnetic-thermal coupling simulation (see Fig 2-B) permits the evaluation of the overheating of the machine depending on the working cycles of the application. Thus, the design can be adjusted to the application performance achieving a thermally optimized machine.

As a particular case, machines excited by means of Neodymium permanent magnets are very sensitive to the temperature. The main drawback of this kind of rare earth magnets is the low Curie temperature, about 300ºC, and the high temperature sensitivity. If the operation point of the magnet remains above the typical knee in its characteristic curve, the change due to temperature is fully reversible. However, if temperature causes a change in the operation point below the knee, the magnet gets partially de-magnetized. So, it is very important to evaluate the rotor temperature in order to prevent any problem related to partial de-magnetization of magnets.

3. FEM analysis of electrical machines

Software based on the Finite Element Method is known as very powerful tool for the design and analysis of electrical motors. FEM analysis makes possible an accurate evaluation of electromagnetic characteristics such as magnetic saturation level, torque or induced electromotive force, and it is also very useful for the computation of the lumped parameters of the equivalent electric circuit, i.e. inductances and iron losses resistance. In this section the basis of this numerical method are explained and some important details related to the FEM analysis of electrical machines are pointed out.

Simulations based on the FEM can be divided in three stages: Definition of the problem, resolution of the problem and results analysis. Next, these three stages are described briefly.

3.1. Definition of the electromagnetic problem

The definition of the problem is faced in three steps. First the geometry must be defined. After that the mesh may be created. And finally physical properties have to be defined. Undoubtedly meshing the problem is the most critical issue in the definition stage. It is widely known that the results accuracy depends strongly on the quality of the mesh. Therefore, normally the higher the amount of nodes is the more accurate the results are. On the other hand, the computational load also depends on the amount of nodes, so that the bigger the amount of nodes is the longer the required simulation time is. Thus, the meshing process is a trade of between accuracy and computational time.

3.2. Resolution of electromagnetic problems

Electric machines are non-linear problems so that they must be solved applying numerical iterative methods. In case of FLUX electromagnetic non-linear problems are solved by the Newton-Raphson method. Regarding to the equations, in the solving process of electromagnetic problems the Maxwell’s equations are considered.

$⇿X\stackrel{\text{⃘}}{E}=-\frac{d\stackrel{\text{⃘}}{B}}{\mathit{dt}}$$⇿\circ \stackrel{\text{⃘}}{B}=0$ $⇿X\stackrel{\text{⃘}}{H}=\stackrel{\text{⃘}}{J}+\frac{d\stackrel{\text{⃘}}{D}}{\mathit{dt}}$

$$⇿\circ \stackrel{\text{⃘}}{D}={ϟ}_v$$E1

(1)$E$is the electric field, $B$is the magnetic flux density, $H$is the magnetic field, $J$is the current density, ${ϟ}_v$is the volume charge density and $D$is the electric flux density. It is important to notice that in electromagnetic problems the fourth Maxwell’s law, that concerning to the electrostatic Gauss’s law, is not considered. In addition to the considered three Maxwell’s laws, the following constitutive laws are taken into account as well.$\stackrel{\text{⃘}}{B}=Ϛ\stackrel{\text{⃘}}{H}\mathrm{}\left[T\right]$

$$\stackrel{\text{⃘}}{J}=ϡ\stackrel{\text{⃘}}{E}\mathrm{}\left[A/m^2\right]$$E2

(2)$Ϛ$is the permeability of the medium, and $ϡ$is the conductivity of the medium. To solve these equations a vector model is used, in which the magnetic vector potential $\stackrel{\text{⃘}}{A}$ is considered as state variable. Applying this vector model the objective function (4) is solved at each node of the mesh, and the value of this state variable is obtained as result at each node.

The relationship between magnetic flux density and magnetic vector potential can be described as,

$$\stackrel{\text{⃘}}{B}=⇿\times \stackrel{\text{⃘}}{A}\mathit{}\left[T\right]$$E3

(3)

Figure 3.

Meshing of the electric machine

Figure 4.

Description of different stages in the FEM simulations

Figure 5.

Representation of the magnetic flux density and the EquiFlux lines

In transient magnetic model the problem can be considered as quasi-steady-state so that the variation in time of the electric flux density can be neglected ($\frac{d\stackrel{\text{⃘}}{D}}{\mathit{dt}}=0$) [1]. According to this assumption, the computed function at each node is simplified.

$$⇿\times \left(v_o\left[v_r\right]\left(⇿\times \stackrel{\text{⃘}}{A}\right)-{\stackrel{\text{⃘}}{H}}_c\right)+\left[ϡ\right]\left(\frac{\complement \stackrel{\text{⃘}}{A}}{\complement t}+\stackrel{\text{⃘}}{\mathit{grad}}\mathrm{}V\right)=0$$E4

(4)$v_o$is the reluctivity of the vacuum ($v_o=\frac{1}{{Ϛ}_o}$), $\left[v_r\right]$is the tensor of the relative reluctivity of the medium, $H_c$is the coercive magnetic field which is considered only in case of permanent magnets and $V$is the electric scalar potential.

3.3. Results analysis

During the resolution of the problem the magnetic vector potential is computed at every node of the mesh. After that, in the third stage different local and global results are analyzed which are derived from the magnetic vector potential. For example, the flux density is computed applying expression (3). Meanwhile, the magnetic flux crossing certain S surface can be computed by the following expression.

$$ϳ=\underset{S}{\cup }\stackrel{\text{⃘}}{B}\circ d\stackrel{\text{⃘}}{s}=\underset{S}{\cup }\left(⇿\times \stackrel{\text{⃘}}{A}\right)\circ d\stackrel{\text{⃘}}{s}=\underset{l}{\cup }\stackrel{\text{⃘}}{A}\circ d\stackrel{\text{⃘}}{l}\mathrm{}\left[\mathit{Webber}\right]$$E5

(5)

The torque produced by the motor is computed applying the well known virtual work method [2].

$$T={\left.\frac{\complement W\mathrm{\text{'}}\left(\mathrm{i},\varpi \right)}{\complement \varpi }\right|}_{\mathrm{i}=\mathit{cte}}=\frac{1}{2}\frac{\complement L}{\complement \varpi }{i}^2+\frac{1}{2}\frac{\complement {\mathfrak{R}}_g}{\complement \varpi }{ϳ}_{\mathit{pm}}+N\frac{{\complement ϳ}_{\mathit{pm}}}{\complement \varpi }i\mathrm{}\left[\mathit{Nm}\right]$$E6

(6)$\mathit{W\text{'}}$is the co-energy stored in the air-gap, $L$is the self-inductance of stator winding, $i$is the supplying current, ${\mathfrak{R}}_g$is the air-gap reluctance, ${ϳ}_{\mathit{pm}}$is the magnet flux, $N$is the number of turns per phase and $\varpi$is the mechanical angle of the rotor. As it is shown in expression (6) the torque produced by permanent magnet machines is the sum of three components. The first component is known as reluctance torque, the second term is known as cogging torque and the third term is known as electromagnetic torque. In the particular case of surface mounted motors, the first torque component is null. On the other hand the cogging torque is oscillatory so that it contributes to the torque ripple. Regarding to the electromagnetic torque it has a main average component but it has also certain torque ripple component. So the total torque ripple in permanent magnet motors is due to cogging torque and electromagnetic torque ripple. In the Fig 6 the cogging torque and the overall torque produced by the motor are shown during an electric period. In this particular case the order of the main cogging torque component is $n=6$meanwhile the order of the main total torque ripple component is$n=12$. In the case study of section 6 de influence of this torque ripple upon the ride comfort in a lift is evaluated.

3.4. Iron losses computation by the FEM

The total iron losses in permanent magnet synchronous machines (PMSM) can be separated in two components: The Eddy current losses in permanent magnets and the iron losses in electrical sheets. Hence, the computation of these losses components in FEM is carried out in two different ways. The Eddy current losses in permanent magnets are calculated as the active power dissipated in these regions considering the magnets as solid conductors. Meanwhile the losses in electrical sheets are computed by the tool so called Loss Surface Model (LSM) which is integrated in the FEM software FLUX for that purpose.

Figure 6.

Different torque components in permanent magnet motors

Loss Surface Model (LSM) for Iron Losses Computation in Electrical Sheets

Normally FEM software do not consider the magnetic hysteresis cycle of electrical sheets in the solving process, but an approximated so called normal magnetization curve is taken into account. As the iron losses depend on the area enclosed by the hysteresis loop, these losses are not computed in the solving process and their calculation is addressed typically as a post processing task.

One calculation option consists in the local computation of the iron losses using the loss separation method and considering local induction waveforms calculated in the solving process. As main disadvantages, the accuracy of this computation method depends on some coefficients, and additional phenomenon as skin effect and minor loops are not normally considered.

FLUXpresents an alternative tool so called Loss Surface Model for iron losses computation. This tool is supposed to be more accurate than the method based on the losses separation theory, because it does not depend on empirical coefficients and it considers additional phenomena such as the skin effect or the minor loops [3]. Due to this fact, in this case the LSM has been used as post processing tool for the calculation of the losses in electrical sheets.

In Fig 7-B the iron losses computation process using the LSM is explained. First of all the electromagnetic problem is solved in the FEM calculating the temporal waveform of induction at each node of the mesh. In the post processing analysis the temporal waveforms of the magnetic field strength are reconstructed from the induction waveforms obtaining in this way local hysteresis loops. This reconstruction of the magnetic field strength is performed using a model based on one $H\left(B,\mathit{dB}/\mathit{dt}\right)$surface curve which is characteristic to each material. This curve must be obtained by experimental characterization of the material. In Fig 7-A the surface curve of a given material is shown. The model that reconstructs the magnetic field strength can be described by the following expression.

$$H\left(B,\mathit{dB}/\mathit{dt}\right)=H_{\mathit{static}}\left(B\right)+H_{\mathit{dynamic}}\left(B,\mathit{dB}/\mathit{dt}\right)$$E7

(7)

In Fig 8 two different hysteresis loops reconstructed by the LSM tool at the same local point are shown. Resulting hysteresis loops with sinusoidal and PWM voltages are plotted. In this way the minor loops caused by the PWM voltages can be appreciated.

Figure 7.

Electric sheet losses computation methodology based on the LSM tool

Figure 8.

Local hysteresis loops reconstructed by the Loss Surface Model with sinusoidal and PWM voltages

Once the local hysteresis loops are obtained, the specific iron losses are computed integrating the area enclosed by each hysteresis loop.

$$P=\frac{1}{T}\underset{T}{\cup }\underset{B}{\cup }H\circ \mathit{dB}\mathrm{}\left[W/m^3\right]$$E8

(8)

Eddy Current Losses in Permanent Magnets

In case of non-fractional machines normally the Eddy current losses in permanent magnets can be neglected. While in fractional machines these losses are more significant due to the high harmonic content of the armature field, so they cannot be neglected [4]. In order to make the proposed method usable to both non-fractional and fractional machines, the calculation of the permanent magnet losses is also addressed. Eddy current losses are due to the asynchronous harmonics of the magnetic field induced by the armature winding in the air-gap. The induced Eddy current density in the magnets can be expressed as in the following equation.

$$J_{\mathit{eddy}}=-\frac{1}{ϟ\mathrm{}}\cup \frac{\complement B}{\complement t}\circ r\circ \mathit{dt}+C\mathrm{}\left[A/m^2\right]$$E9

(9)

Where $J_{\mathit{eddy}}$ is the induced Eddy current density in the magnets, $ϟ$is the electrical resistivity of the magnets, $B$is the armature magnetic field and $r$is the radius of the point in which the current density is represented. The Eddy current losses are computed using the following formula.

$$P_{\mathit{eddy}}=\frac{2p}{T}{\cup }_0^{Ϟ}{\cup }_{-}^{\frac{Ϗ}{2}}{\cup }_{R_1}^{R_2}\mathrm{}{J_{\mathit{eddy}}}^2\circ ϟ\circ r\circ \mathit{dr}\circ \mathit{d\varpi }\circ \mathit{dt}\mathrm{}\left[W\right]\mathrm{}$$E10

(10)

Where $R_1$ and $R_2$are the inner and the outer radius of the magnets respectively, $T$is the period of the fundamental component of the armature field and $Ϗ$is the magnets span in mechanical degrees [5]. It is necessary to consider in the formula the constant C in order to assure that the induced Eddy currents are limited to one pole piece. That means that the different magnet poles are supposed to be electrically isolated.

Regarding to the computation of the magnetic field induced by the armature winding in the air-gap, there are two options for that purpose depending on if the Eddy current reaction field is neglected or not. In case this reaction field is neglected, the Eddy currents are considered resistance limited and this way they can be calculated performing magneto static simulations. However, in case the Eddy current reaction field has to be considered, it is necessary to perform time step simulations. In the proposed method, as the machine is co-simulated with MATLAB-SIMULINK⁽, the simulations are time stepped so the Eddy current reaction field is taken into account.

As far as the FEM simulation is concerned, magnets have to be defined as solid conductors. As the permanent magnets are electrically isolated, each one has to be associated to one solid conductor. In Fig 9 the electrical circuit associated to the permanent magnets is represented. It can be seen how each solid conductor has connected in parallel one resistor. It is because in the FEM software all the electrical components have to be connected in close circuit. These parallel resistances are set to very high values in order not to have influence on the losses computation.

Figure 9.

Associated electrical circuit to the permanent magnets

4. Matlab-simulink coupling to FEM software

Being MATLAB numerical calculus software, FEM calculus can be carried out in MATLAB without any external software. However the tendency is to use specialized FEM software for each physical domain like Electromagnetic, Fluid dynamics and so on. For example Cedrat’s FLUX software is well known electromagnetic FEM analysis software. Using FLUX enables to design and study a wide range of electromagnetic devices such as different electrical machines, actuators or sensors.

However this software is usually generalist respect to the domain and an application approach of the software is usually desired for a more friendly use of the software by beginner users. That’s why companies like Cedrat are providing programmatic tools for external piloting of their calculus Kernel.

Using some libraries, user can pilot a FLUX file that has been prepared previously from MATLAB for instance. Although several parameters are required to prepare the flux file, using external piloting, it is possible to limit the parametric study to a short number of parameters that are required to be visible by application user. The use of FLUX calculus Kernel is now more application oriented and any user without FEM knowledge can use the calculus Kernel.

Another MATLAB FEM collaboration example is the co-simulation. Usually the electromagnetic device that is being designed is connected to other external devices. Complex mechanical models, electrical circuits or hydraulic systems are some examples of models that can be easily modeled by MATLAB-SIMULINK but cannot be included in the FEM software.

The simplest approximation to solve this multi-domain problem is the simplification of the electromagnetic device to an equivalent electric circuit based on lumped parameters such as resistance and inductances. This can be done with different levels of accuracy. For example, the inductance of a machine can be considered constant respect to the current of the machine. In a more precise approximation the saturation effect can be considered introducing a Look-up table to take into account the change of the inductance respect to the magnetic flux of the machine. In both cases, first FEM simulations are carried out and the values for the inductance are post processed to include them later in the MATLAB-SIMULINK simulation. However this approximation cannot take into account torque ripples, instantaneous saturations and so on. If more accuracy is required, it is possible to couple the FEM problem to the MATLAB-SIMULINK simulations. Cedrat provides a block to use in MATLAB-SIMULINK. In this interface block inputs and outputs should be defined". Then, the FEM problem is solved at each time step of MATLAB-SIMULINK simulation. The convergence of this type of simulations is difficult but the accuracy is rather high.

5. Case study 1: Iron losses computation in electrical machines

It is widely known that the performances of electric machines depend strongly on the losses. The efficiency, the thermal behavior or the compactness are some of the design constrains which are strongly influenced by the losses. So it is very important to predict these losses accurately if an optimum design of the electric machine is required. The losses in electric machines can be divided in three main components: Joule losses, iron losses and additional or stray losses (mechanical losses, induced Eddy current losses in frames, etc). The Joule losses are relatively easy to compute from the armature resistance and current values. However, the other two components are rather more complicated to estimate. This work is focused on the calculation of iron losses which can be relevant in some cases depending on the machine topology, speed, saturation levels, etc.

Nowadays the majority of the electric machines are integrated in variable speed drives. That means that the motors are supplied by modulated voltages such as Pulse Width Modulated (PWM) voltages. It is well known that the additional current harmonics introduced by PWM voltages increase the iron losses. So it is necessary to consider the voltage characteristic in the computation process in order to achieve accurate results in loss computation.

The most extended way to define the iron losses is by the widely known loss separation method in which they are divided in three main components: hysteresis losses, classical Eddy losses and excess losses [6,7]. In case of non sinusoidal voltages, the high order harmonics of currents can accentuate the presence of some phenomena such as the skin effect or the minor loops which are not directly taken into account by theses classical methods. Many efforts have been dedicated to the adaptation of the initial definition of the iron losses given by Bertotti in [6] to consider all these phenomena. For instance, in [8] Boglietti et al. completes the expression of the classical Eddy losses component with the aim of taking into account also the skin effect. In [9] Toda et al. proposes a modification of the hysteresis losses formula in order to consider the minor loops. It is also possible to find in the literature some papers that deals with the iron losses computation under arbitrary supply voltages, mainly focused on pulse width modulated (PWM) voltages [8,10,11].

In this section a novel method for iron losses computation under non-sinusoidal supply voltages is proposed. The method consists in coupling the FEM software FLUX from CEDRAT Company with the system simulator software MATLAB-SIMULINK. The iron losses are computed by a post-processing analysis separating the computation of the Eddy current losses in the permanent magnets (PM) and the iron losses in electrical sheets. The simulation results have been compared with experimental results in order to validate the proposed method.

5.1. Iron losses computation by co-simulations

In the proposed methodology the steady state iron losses are computed by a co-simulation between the FEM software and the system simulator software. The power system comprising a two-level converter is simulated in MATLAB-SIMULINK⁽, whereas the electrical machine is modeled in the FEM software. This way it is possible to evaluate the performance of the electrical machine under the functioning conditions imposed by the converter, in a relatively short and rather accurate way.

The FEM software permits the coupling of electric circuits for transient analysis. Generally these circuits comprise passive elements such as resistances and inductances along with sinusoidal voltage and current sources. Nevertheless, if more sophisticated circuits are required to consider, for example a power converter, it is necessary to couple the FEM software with the circuit simulator software.

When the machine is supplied by a converter, it is driven by voltage. In these cases there is an inevitable transient period which takes several electrical cycles before the motor reaches the steady state. It is considered that at least a resolution of 10 samples per switching period is necessary to achieve accurate results. The simulation of several electric periods with so small step size would lead to such an amount of samples and such long time consumption that current driven simulations are chosen instead of voltage driven simulations. In current driven simulations it is enough considering only one electric period to compute the iron losses, which enables to decrease significantly the time consumption of the simulation. For example, taking 10 samples per switching period, one co-simulation characterized by a switching frequency of fs=4KHz and a period of T=25.82 sec, is carried out with a time step of T=25sec and it takes approximately 30 minutes.

The proposed method for iron losses computation is structured in three stages as it is shown in Fig 10. First of all the lumped parameters of the equivalent electrical circuit of the machine are calculated by the FEM. From this previous simulation the machine can be modeled with lumped parameters such as d-q axis inductances Ld and Lq, and the permanent magnet flux pm. In the second stage the transient model of the machine in d-q axes is implemented in MATLAB-SIMULINK using the lumped parameters previously estimated and a two level inverter is implemented as power supply. In this simulation the current waveforms in steady state are obtained for PWM voltages. Finally the co-simulation between the FEM software and MATLAB-SIMULINK is carried out driving the FEM model of the electrical machine by current waveforms obtained in the previous stage. In this third stage the total iron losses are computed as the sum of the Eddy current losses in permanent magnets and the iron losses in electrical sheets.

In Fig 11 the total iron losses estimated by co-simulations for sinusoidal and PWM voltages are compared. It can be observed that the iron losses with PWM voltages are higher. Meanwhile in Fig 12 the resulting currents and electromagnetic torque are plotted for sinusoidal and PWM voltages. It can be noticed a significant difference in the ripples because of the PWM voltages.

Figure 10.

Proposed methodology for iron losses computation under non-sinusoidal voltages

Figure 11.

Iron losses computed with sinusoidal and PWM voltages with $f_s=4\mathit{KHz}$

Figure 12.

Co-simulation results obtained under sinusoidal and PWM voltages with $f_s=4\mathit{KHz}$

5.2. Experimental results

In this section some experimental results are shown in order to prove the accuracy of the proposed method. For that purpose 7.4KW brushless motor has been analyzed in a test bench and the iron losses have been measured.

Description of the experimental test

In Fig 14 the test bench layout is shown. It consists of two machines, the one under test and the load motor, one torque sensor, one commercial controller for each machine and one power analyzer. During the test the controlled variables are the speed in case of the load machine and the torque in case of the machine under test. The test has been carried out controlling the d-axis current to zero avoiding in this way the functioning of the machine in the flux weakening operation mode. To calculate the losses, the loss segregation method has been implemented.

$$P_{\mathit{in}}-P_{\mathit{cu}}-P_{\mathit{loss}}-P_{\mathit{mec}}=0\mathrm{}$$E11

(11)

Figure 13.

Description of the performed tests in order to validate the sinusoidal losses computation method

Figure 14.

Test bench configuration

Where $P_{\mathit{in}}$ is the active power consumed by the machine, $P_{\mathit{cu}}$are the Joule losses, $P_{\mathit{mec}}$is the mechanic power in the motor axis and $P_{\mathit{loss}}$ is the sum of the iron losses and the mechanical losses. The value of $P_{\mathit{in}}$ is calculated by the power analyzer which is connected to the input of the machine. To calculate the Joule losses, a high accuracy ohmmeter has been used due to the fact that the power balance is very sensitive to the errors in the resistance measurement. The mechanical power is calculated from the torque and speed measurements. As the test is carried out at constant speed, the mechanical power can be calculated as the following.

$$P_{\mathit{mec}}=T\circ {\mathrm{\chi }}_m\mathrm{}\left[W\right]\mathrm{}$$E12

(12)

Where $T$is measurement provided by the torque sensor and ${}_m$ is the mechanical speed measured by the encoder. In the power balance the only unknown variable is $P_{\mathit{loss}}$ which involves two loss components, the iron losses and the mechanical losses. As the tests are performed at relatively low speeds (290 rpm) and the machine has been heated previously, the mechanical losses can be considered negligible in comparison with the iron losses. In addition all the tests have been carried out at constant speed in order to maintain the mechanical losses constant.

The configuration of the test bench does not allow supplying the machine with sinusoidal voltages. The possible solution would be to add RLC filters between the inverter and the machine or to use sinusoidal voltage supplies instead of the inverters. However in this case it has not been possible to implement neither of these solutions so finally it has been necessary to implement an alternative solution. The power analyzer comprises low pass band filters for the measured currents and voltage, which has been set to 500Hz. This way only the fundamental components of voltages and currents are considered for the computation of the power. The calculated power in this way is considered to be the power for sinusoidal voltages. However, before implementing this method, it has been validated performing two different tests which are described in Fig 13.

First of all the so called generator test is carried out. In this test a resistive load is connected to the terminals of the machine under test, which operates as a generator, and it is driven by the load motor at constant speed. This way the currents are sinusoidal and so the computed iron losses correspond to the sinusoidal voltages. After that the second test is performed using the inverter and with the machine working exactly at the same operation conditions that in the previous test. In this case the input filters of the power analyzer are configured to 500Hz. In both cases the resulting P_loss are very similar, with a difference less than 5%. Hence it is considered that the losses calculated with the filter are equivalent to the sinusoidal losses.

Validation of the simulation results

The proposed method to calculate the iron losses has been validated with experimental results. As the machine iron losses are closely related to the stator current density, the simulations and the tests have been done at same phase current magnitude. The switching frequency of the PWM voltage is fixed to $f_s=4\mathit{KHz}$ during all the tests.

In Fig 15 the simulated torque and the experimental torque at different current values are shown for both cases, sinusoidal and PWM voltages. It can be noticed how the simulated and the experimentally measured torques are very similar.

Figure 15.

Torque with Sinusoidal and PWM voltages

In Fig 16 the computed iron losses with sinusoidal and PWM voltages are shown. It is remarkable that certain differences arises between the co-simulation and the experimental results either with sinusoidal or PWM voltages. These differences are considered to be because of the mechanical losses which are not taken into account by the co-simulations. This difference is rather constant respect to the current and very similar in both cases with sinusoidal and PWM voltages. This behavior of the difference between simulations and test fits perfectly with the supposition that the difference is due to the mechanical losses.

Although there is a slight error in the prediction of the iron losses, an important point of this work is that the proposed methodology enables to calculate rather accurately the increase of the iron losses due to PWM voltages. In Fig 17-A the comparison between sinusoidal and PWM iron losses is shown. A coefficient can be defined to represent the increase of the losses due to the PWM modulation. This coefficient is defined as PWM losses factor $K_{\mathit{pwm}}$ and it is calculated as following.

$$K_{\mathit{pwm}}=\frac{P_{\mathit{PWM}}}{P_{\mathit{sin}}}\mathrm{}$$E13

(13)$P_{\mathit{sin}}$are iron losses with sinusoidal voltages and $P_{\mathit{PWM}}$ are iron losses under PWM voltages. As it is shown in Fig 17-B the PWM loss factor is approximately constant for every load conditions. So it can be concluded that the losses increase due to PWM voltages hardly depends on the current value. In the procedure presented in this work, the loss terms which are not taken into account are independent from the modulation which enables us to calculate very accurately the loss increase in spite of making a slight error in the loss estimation.

Figure 16.

Iron Losses with Sinusoidal and PWM voltages

Figure 17.

Iron Losses and PWM losses factor $K_{\mathit{pwm}}$

6. Case study 2: Comfort analysis in elevators

The performance of the electric motor in an elevator is going to be studied by co-simulation. As the elevator is a very sensitive direct drive application, it is very important the smoothness of the torque. So in this case the influence of the torque ripple on the elevator performance is going to be evaluated.

In this study the overall elevator system is considered taken into account the control unit, the traction motor and the mechanical system. Thus, the control unit and the mechanical system are simulated in MATLAB-SIMULINK, while the electric motor is simulated in the FEM software FLUX.

6.1. Application description

An elevator is a complex system in sense that it may present resonances due to the elastic properties of ropes and the elasticity and damping effect between the frame and the cage of the lift [12]. Moreover, the resonant frequency of the system is variable, and depends upon the passenger load.

The elevator taken into analysis is a so-called 1:1 roping configuration which is shown in Fig 18-A. As main feature, it has a maximum speed of 1 m/s and a maximum acceleration of 0.5 m/s², with a capacity of 6 passengers. Considering a pulley diameter of 50 mm, rated value of 20 rad/sec and 200 Nm are required for the drive.

In Fig 18-B the frequency response of the acceleration of the cage to pulley torque variations for different load states is shown.

Figure 18.

Description of the analyzed elevator system

6.2. Co-simulation scheme

In Fig 19 the co-simulation block diagram is shown. The control strategy, the speed reference generator and the elevator model are implemented in MATLAB-SIMULINK, while a precise machine model is simulated in FLUX.

“Coupling with flux2d” function of MATLAB-SIMULINK library allows data transfer between both simulation tools. In this function the following parameters have to be specified: the name of the FLUX file, the input and output parameters, the sample time and the memory size to allocate. The inputs are the currents of phase A and B, and the mechanical speed. With these inputs, FLUX model returns as outputs the electromagnetic torque and the position of the rotor. The torque is the input for the elevator model which computes the linear acceleration of the cage and the angular speed of the pulley. The speed is fed back along with the position which is aligned with the flux vector.

Figure 19.

Co-simulation block diagram

Regarding to the speed of the elevator, it has to follow a specific track in order to get a good ride quality. In this case, this track is obtained by means of a speed track generator. The speed control is a vector control implemented in a d-q frame oriented to the magnet flux vector. The block diagram of that control is shown in Fig 20.

As it is shown in Fig 20, the electromagnetic torque is controlled acting upon Iq component. In this application flux weakening is not necessary, so I_d component is null. Notice that current control is not implemented. The dynamic of current loops is so quick in comparison with the dynamic of the mechanical system that they can be omitted. So, control consists in a PI for speed regulation, a frame transformation from d-q to - and an inverse CLARKE transformation from - to abc.

Figure 20.

Speed Control Diagram

On the other hand, a proper choice of sample time is very important. It must not be too small in order to avoid too long simulation times. Otherwise, it has to have enough precision if accurate results are wanted. In this case, the aim of co-simulation is to evaluate the effect of torque ripple on vertical acceleration of the lift, so the sample time has to be small enough to get a proper resolution of that ripple. As it has been already explained, the frequency of the main torque ripple component is six times the fundamental frequency which is equal to 19 Hertz. So the torque ripple has a period of 8.77ms. Thus, a sample time of 500s is chosen so as to use 17 samples for the definition of a period of the torque ripple.

Another critical aspect that affects final results is the mesh resolution of the machine model in FLUX. A too accurate resolution may slow down the simulation. As the sample time, the mesh resolution has to be carefully set in order to get enough accuracy in the results in a reasonable simulation time.

6.3. Results

In this subsection the co-simulation results are presented. The simulation has been carried out under 50% of rated conditions, with a load of 3 passengers. Two machine models have been taken into account: an ideal one which consist in a constant gain equal to the torque constant$K_T$, and the FEM model (see Fig 19). In Fig 21-A it could be seen the track of the lift speed in which a little overshoot can be noticed. That is because the dynamic of the regulator has been limited in order to avoid the excitation of the system resonance by the control. In Fig 21-B the torque produced by the FLUX model is presented.

Figure 21.

Torque and speed obtained by co-simulation

In Fig 22 the resultant lift accelerations for an ideal machine model and for the FEM machine model are shown. Unlikely to the ideal case, the torque of the FLUX model is not as smooth as application requires. Consequently, that torque ripple excites the mechanical resonances of the system, causing vibrations in the vertical acceleration. This phenomena lead to a notary worsening in the ride quality of the elevator.

Figure 22.

Vertical acceleration of the cage obtained by simulation

7. Conclusions

In this chapter a co-simulation tool for analysis and design of electrical machines is proposed.

Two different case studies have been reported. In the first one co-simulation between the system simulator software and the FEM software are presented for the iron losses calculation under arbitrary voltages. A particular procedure is described to carry out current driven co-simulations in order to reduce as much as possible the computation load and the time consumption of the calculus.

The proposed method is validated by experimental tests. The estimated iron losses under sinusoidal and PWM voltages have a good agreement with the iron losses measured experimentally. So it can be stated that using this method it is possible to estimate the iron losses under arbitrary voltages without the need of prototyping, which could lead to less expensive and more flexible design processes.

In the second case study the performance of the electric machine in an elevator is evaluated. Thus, the vibrations caused by the motor in the elevator cage are studied applying the proposed co-simulation. It has been demonstrated that the torque ripple of the motor can affect negatively on the elevator comfort.

Probably the main drawback of these co-simulations is the required computation time. They ought to be carefully configured, on one hand to avoid too large computation times, and on the other hand to achieve an adequate accuracy.

As main conclusion it could be remarked that these co-simulations are very helpful tools which may allow accurate and relatively easy analysis of multi-domain systems in a relatively short development time.