**Continuous Control Set Model Predictive Control of a Switch Reluctance Drive Using Lookup Tables**

**Alecksey Anuchin 1,2,\*, Galina L. Demidova 2, Chen Hao 3, Alexandr Zharkov 1, Andrei Bogdanov <sup>2</sup> and Václav Šmídl <sup>4</sup>**


Received: 18 May 2020; Accepted: 22 June 2020; Published: 29 June 2020

**Abstract:** A problem of the switched reluctance drive is its natural torque pulsations, which are partially solved with finite control set model predictive control strategies. However, the continuous control set model predictive control, required for precise torque stabilization and predictable power converter behavior, needs sufficient computation resources, thus limiting its practical implementation. The proposed model predictive control strategy utilizes offline processing of the magnetization surface of the switched reluctance motor. This helps to obtain precalculated current references for each torque command and rotor angular position in the offline mode. In online mode, the model predictive control strategy implements the current commands using the magnetization surface for fast evaluation of the required voltage command for the power converter. The proposed strategy needs only two lookup table operations requiring very small computation time, making instant execution of the whole control system possible and thereby minimizing the control delay. The proposed solution was examined using a simulation model, which showed precise and rapid torque stabilization below rated speed.

**Keywords:** switched reluctance motor drive; model predictive control; continuous control set; pulse-width modulation; magnetization surface; electrical drive

### **1. Introduction**

Switched reluctance drives (SRDs) have attracted a significant amount of attention in recent decades as the most promising type of electric drive. With the spread of 3D-printing technologies, this drive is one of the best candidates to become the first commercial 3D-printed machine [1] and modular designed machine [2]. Being very simple, the machine requires a sophisticated control strategy, which currently has no general approach, unlike the field-oriented control (FOC) or direct torque control (DTC) used in AC electric drives. Many papers have considered the mitigation of torque ripple by adjustment of the commutation angles [3–5] and current profiles [6,7], or optimization of the motor magnetic geometry [8]. However, solutions suggested in [3–5] help in the limited range of speeds as the current slope varies with the speed also affecting the produced torque. Current profiling suggested in [6,7] cannot be applied to each particular motor without hand tuning in all operation modes. In addition, optimization of the magnetic geometry for torque ripple minimization [8] contradicts the goal of efficiency optimization of the electrical machine.

The approach with model predictive control using a finite control set (FCS) was suggested in [9] to help decrease torque ripple and make the control strategy applicable for any motor. It was improved in [10] by decreasing the losses in the motor and in [11] by achieving equal distribution of the temperatures of the power converter's IGBT-modules. These approaches utilize the magnetization surface of the switched reluctance motor for torque estimation while selecting the best control command. These systems stabilize the output torque by profiling the phase current shape, not only with respect to the commanded torque, but also minimizing ohmic losses. The disadvantage of the finite control set model predictive control is the random switching rate and the error in the commanded torque or speed if a short prediction horizon is used. Continuous control set (CCS) model predictive control solves the voltage equations to obtain the references to the inverter, but this approach is problematic for a switched reluctance machine (SRM) due to the high nonlinearity of the controlled plant as the motor inductance varies with current and the power converter frequently saturates.

During analysis of the model predictive control (MPC) operation it was noticed that, for each electrical revolution for the same torque command, the current references repeat. The cost function, which considers the difference between the commanded and predicted torque, and minimizes ohmic losses, gives the predefined shape of the current references with respect to the rotor position. Together with the magnetization surface of each phase, this can be used to evaluate the desired voltage command, which adjusts the state variables of the drive to follow the references. Compared to the existing methods, the proposed solution uses information about the magnetization map of the machine to obtain the current references for the demanded torque, and then evaluates voltage references using the magnetization surface. The control strategy uses only lookup table operations, which allows evaluation of the voltage command nearly instantly, thus requiring very small computation efforts.

### **2. Switched Reluctance Drive Model**

### *2.1. Drive Topology*

The switched reluctance drive can have a varying number of phases and its power converter topology can vary. The most common configuration of the drive is shown in Figure 1. The stator of the motor contains three phases with concentrated windings located at six stator teeth. The rotor has four teeth, and the electromechanical reduction ratio is 4. The power converter is represented by three asymmetrical H-bridges, which feed the motor phases with unipolar current.

**Figure 1.** Switched reluctance drive topology.

### *2.2. Motor Equations*

Each phase has its voltage balance equation [10]:

$$\frac{d\psi}{dt} = \upsilon - i\mathcal{R}\_\prime \tag{1}$$

where *v* is the applied voltage, *i* is the flowing current, *R* is the phase resistance, and ψ is the flux linkage of the winding.

The magnetic system of an SRM is highly non-linear, which is usually defined by means of a magnetization surface/map or flux linkage map. For AC machines these maps can be identified online using a testbench that consists of the motor under test and a prime mover motor [12]. This testbench can realize all possible operation conditions and uses high-frequency injection to estimate differential inductance maps in order to build flux linkage maps. Such a method is applicable for SRMs; however, other methods are more popular due to the specific nature of this machine.

The motor phases of an SRM are usually considered to be independent without, or with very little, magnetic coupling between them [11]. Therefore, the magnetization surface of each phase depends on its flowing current and current rotor angular position as shown in Figure 2. The magnetization map of an SRM can be obtained experimentally [12,13]. Paper [12] provides a thorough description of experimental setups for making flux linkage maps. In [13], authors represent two methods of flux linkage map measurement based on the example of a four-phase SRM. The first method utilized a testbench with a mechanical rotor locking device and measured current and voltage. The flux linkage was calculated offline using experimental data. The second method was an online method. It did not require a rotor locking device, and online data from current and voltages sensors were used to evaluate the magnetization map when the SRM was running in a normal operation mode. For the considered solution, the offline method is preferable, as the online method provides an incomplete data set limited by the operation conditions of the drive, which makes estimation of the torque surface impossible.

**Figure 2.** Magnetization surface of the switched reluctance motor (SRM).

The torque of the phase at the electrical speed can be defined as a function of the phase energy:

$$T = -\frac{\partial \int\_0^{\psi} \mathbf{i}(\psi, \theta) \mathbf{d}\psi}{\partial \theta},\tag{2}$$

or co-energy:

$$T = \frac{\partial \oint\_0^i \psi(i, \theta) di}{\partial \theta} \tag{3}$$

where θ is the electrical rotor angular position.

The total torque of the machine is the sum of the torques produced by each phase. For a considered motor with number of phases *N* equal to 3, the torque equation can be written as follows:

$$T = \sum\_{n=1}^{N} T\_n = T\_1 + T\_2 + T\_3. \tag{4}$$

### *2.3. Power Converter*

Each phase of the motor is connected to the asymmetrical bridge as shown in Figure 1. The phase current can only flow in one direction, which is considered to be the positive direction. It is possible to apply a positive voltage of the DC link by switching on both switches in the bridge, zero voltage by switching on only one of the switches, and negative voltage by switching off both switches. In the last case the phase current flows through freewheeling diodes, and the negative voltage is applied only while the flowing current is bigger than zero.

### **3. Continuous Control Set Model Predictive Control for SRD**

The FCS MPC, which is taken as the reference, was proposed in [10]. The MPC algorithm checked 27 possible states of the power converter and estimated the behaviors of the phase currents and torque for each case. The cost function was developed to keep total torque closer to the reference with minimum ohmic losses in the phase windings.

As a result, for the infinite switching rate, the current shape was profiled to achieve nearly constant torque as shown in Figure 3a. The waveform of the phase current contains the period of growing current, when the phase operates together with the previously on phase, which is reaching the point of maximum inductance (aligned position). When the previous phase reaches the aligned position, it is not desirable to keep current in it anymore as it does not produce positive torque. So, the previously on phase turns off, and all the torque, produced by the motor, is delivered to the shaft by a single phase at its second time interval. The current is regulated according to the magnetization surface of the motor to produce the commanded torque. Then, the inductance of this phase reaches its maximum value. Its current is reduced, and the next phase is switched on.

This current waveform is repeated each revolution; therefore, it is possible to evaluate it in advance (offline) for each torque reference.

Having variable rotor angular positions and torque references, the solution for each phase current can be represented as a function of these two variables. This function is easy to implement using a lookup table. An example of the current reference surface for phase A is shown in Figure 3b.

The current regulation can be performed using the magnetization surface of the motor phase. With the rotation of the machine and change of the current reference, the system should move from one point of the surface to another.

At the end of the pulse-width modulation (PWM) cycle, the current is to be measured and the control system should be executed. The control system has information about the measured current *i*[*k*] (see Figure 4) and current rotor electrical position θ[*k*]. Knowing the motor speed ω, the rotor position in the end of the next PWM cycle can be evaluated by:

$$
\Theta[k+1] = \Theta[k] + a \cdot T\_{\text{PWM}} \tag{5}
$$

where *TPWM* is the duration of the PWM cycle.

Knowing the torque reference *Tref*, the current reference for the next PWM cycle can be evaluated from the lookup table as:

$$i\_{ref}[k+1] = \mathbf{f}\_{\text{current reference}}\{\hat{\Theta}[k+1], T\_{ref}\}.\tag{6}$$

Thus, the control system knows the phase currents and current rotor position, and current references and predicted rotor position for the next PWM cycle. The system moves along the magnetization surface of the motor phase from point *k* to point *k* + 1 as shown in Figure 4.

**Figure 3.** Current waveform for constant torque operation (**a**) and the phase A current reference surface example (**b**).

**Figure 4.** Transient during a single pulse-width modulation (PWM) cycle.

For both points the flux linkage can be estimated as:

$$\begin{aligned} \label{eq:SDAR} \Phi[k] &= \text{f}\_{\text{magnetization surface}}(\Theta[k], i[k]);\\ \Phi[k+1] &= \text{f}\_{\text{magnetization surface}}(\Theta[k+1], i\_{ref}). \end{aligned} \tag{7}$$

These data can now be used to solve Equation (1) in order to find the voltage command:

$$v\_{ref} = \frac{\hat{\Psi}[k+1] - \hat{\Psi}[k] - \frac{i[k] + i\_{ref}}{2}R \cdot T\_{PWM}}{T\_{PWM}} \,\tag{8}$$

which depends on the difference between the next and current flux linkage estimations taking into account the voltage drop across the phase resistance by the average current that will flow in the winding during the next PWM cycle.

The duty cycle for the asymmetrical bridge control can be evaluated using the voltage command and the current DC link voltage *vDC*:

$$\gamma = \frac{v\_{ref}}{v\_{DC}} \tag{9}$$

which should be limited in the range of [−1; +1], where +1 corresponds to fully on switches and −1 corresponds to fully off.

### **4. Control System Implementation**

### *4.1. Evaluation of the Magnetization Surface and Torque Surface as the Functions of Phase Current and Rotor Angular Position*

The surface from Figure 2 should be represented by an array of *N*·*N* points. As the resulting array lies in the memory of the control system (microcontroller or field-programmable gate array (FPGA), the number of points *N* should be selected taking into account this constraint as well as the accuracy of the surface representation. If the number of points is sufficiently large, it is possible to use simple bilinear interpolation or fetch the nearest value for fast evaluation of the control Equation (7). The array of 100 points in both dimensions requires 10,000 words of memory for surface representation, which is suitable for most modern microcontroller devices. This representation is convenient if the online estimation of the flux linkages is running in parallel with the control system adjusting the reference points in the magnetization map [14].

Having a magnetization surface representation, it is possible to obtain the torque surface, which can be also represented as an array of *N*·*N* points. The torque can be evaluated as the derivative of the co-energy with respect to the angle. The co-energy integral can be expressed with series for numerical evaluation using Tustin's method as:

$$
\Delta \mathsf{W}'(\dot{\imath}\_{N\_i}, \boldsymbol{\Theta}\_{N\_\Theta}) = \int\_0^{\dot{\imath}\_{N\_i}} \boldsymbol{\upmu}(\dot{\imath}, \boldsymbol{\upTheta}\_{N\_\Theta}) \mathrm{d}\boldsymbol{\upepsilon} \approx \sum\_{n\_i=2}^{N\_i} \left( \frac{\boldsymbol{\upmu}[n\_i - 1, N\_\Theta] + \boldsymbol{\upmu}[n\_i, N\_\Theta]}{2} \right) . \tag{10}
$$

where *iNi* and θ*N*<sup>θ</sup> are the current and angle of *Ni*- and *N*θ-points in the magnetization surface table, respectively. Each point on the surface can be evaluated as a difference between co-energies in a one-step clockwise direction from the current rotor angular position and a one-step counterclockwise direction divided by the angle difference using the following equation:

$$\begin{array}{rcl} T[N\_i, N\_\Theta] &=& \frac{\Lambda W'(i, \Theta)}{\Lambda \Theta} \\ &=& \frac{\sum\_{i\_j}^{N\_j} \left( \left( \frac{\Phi\left[n\_i - 1, N\_\Theta + 1\right] + \Phi\left[n\_i N\_\Theta + 1\right]}{2} \right) \frac{l\_{\text{max}}}{N} \right) - \sum\_{i\_j = 1}^{N\_j} \left( \left( \frac{\Phi\left[n\_i - 1, N\_\Theta - 1\right] + \Phi\left[n\_i N\_\Theta + 1\right]}{2} \right) \frac{l\_{\text{max}}}{N} \right)}{\left(2 \cdot \frac{2\pi}{N}\right)} \end{array} \tag{11}$$

which can be simplified as:

$$T[N\_i, N\_\Theta] = \frac{I\_{\max}}{8\pi} \cdot \sum\_{n\_i=2}^{N\_i} (\psi[n\_i-1, N\_\Theta+1] + \psi[n\_i, N\_\Theta+1] - \psi[n\_i-1, N\_\Theta-1] - \psi[n\_i, N\_\Theta-1]). \tag{12}$$

### *4.2. Evaluation of the Current Reference Surface*

The referenced torque can be produced by an infinite number of combinations of the phase currents. If minimum ohmic losses are desired, then there is a single solution for a given torque reference and angular position. Thus, the goal is to evaluate the current reference for each phase as a function of the torque reference and the rotor angular position.

For the three-phase SRM, there are angular positions when only one phase can produce a positive torque, or positions when two phases can produce it. Consider the situation for two phases generating positive torque.

For the used linear model of the motor, it is possible to find an analytical solution for the optimal current references, but this model is used only to simplify the simulation model. In the general case, the magnetization surface has high nonlinearity, and cannot be represented as a Fourier or Taylor series with a reasonable number of coefficients for desired accuracy [15]. Thus, the best option is to use a lookup table for the current reference surface as well as for the torque and flux linkage surfaces.

Consider some small positive angle (point O in Figure 3a) where the positive torque can be produced by phases A and C. Set the initial phase A and C currents to zero. The phase C current can be adjusted by integral torque controller, which has some gain *k* as shown in Figure 5a. The integral torque controller regulates the phase C current until the error in the torque reaches the permitted tolerance. Thereafter, the phase A current starts to increase (see Figure 5b). The torque controller maintains its operation, continuously time varying the phase C current and keeping the total torque from both phases A and C close to the reference.

From the moment when the phase A current starts to grow, the system begins to track the minimum total current as it determines the ohmic losses in the motor windings. In the end of the transient from Figure 5b, the outputs hold the most efficient references of the phase currents *ia ref* and *ic ref*, which can be used as points on the current reference surfaces for each phase as shown in Figure 3b.

A similar approach can be used for the angles when the positive torque is generated by one phase only. The only difference is that the search for minimal loss points is no longer needed due to the fact that total torque is produced by a single phase, and there is no extra degree of freedom.

**Figure 5.** Evaluation of the current references for the given torque command: (**a**) block diagram, (**b**) transient during optimal reference evaluation.

### *4.3. Control System Flowchart*

The control system is divided into an offline time-consuming algorithm placed into the initialization stage, as shown in Figure 6a, and a real-time control interrupt (see Figure 6b), which implements the proposed model predictive control strategy.

First, for the given representation of the magnetization surface, the initialization procedure evaluates the flux linkage surface for the desired resolution of the two-dimensional array, for example, 100·100 points. Next, these data are used to evaluate the torque surface using Equation (12). Finally, the most time-consuming algorithm is executed, which evaluates current references with respect to the range of all possible torque commands and angles using the iterative algorithm shown in Figure 5. Thereafter, the control system can start operation and rotate the machine.

The interrupt routine, which is shown in Figure 6b, should be executed at the end of the PWM cycle. The control system reads the feedback from the current, voltage, and position sensors inside this interrupt. At first, it uses the current reference surface to obtain current references for the next computed rotor angular position using Equation (5). Then, it obtains flux linkages for the current and next referenced steps from the magnetization surface. Finally, the voltage command and the duty cycles are evaluated and applied as the references to the PWM generator just before the new PWM cycle starts.

**Figure 6.** Flowchart of the control system: (**a**) initialization procedure and background loop, (**b**) interrupt with MPC strategy.

### **5. Simulation Results**

There are various approaches to building a model of SRM. One of the most convenient options is to use a linearized magnetization profile [16], as shown in Figure 7, which allows fast simulation with simple equations for torque estimation to be undertaken. Phase inductance below the saturation knee was found from the following equation:

$$L = L\_{\text{av}} - \Delta L \cos \Theta,\tag{13}$$

where *Lav* is the average inductance, Δ*L* is the half-sum of maximum and minimum inductances, and θ is the rotor electrical angular position. When operating above saturation current *isat* the phase differential inductance becomes equal to the minimal inductance of the motor for a misaligned position as assumed in Figure 7.

**Figure 7.** Linearized magnetization profile of SRM model.

The Euler integration method was used to obtain flux linkage in each phase and the resulting equation is given below:

$$
\Psi[m] = \Psi[m-1] + (\upsilon - i\mathbb{R})h,\tag{14}
$$

where ψ[*m*] and ψ[*m* − 1] are the new and previous flux linkages of the motor model, respectively, and *h* is the integration step size.

The power converter, which contains an asymmetrical H-bridge, is depicted in Figure 1. It can only produce a positive current in any phase, which should be taken into account in the model when zero or negative voltage is applied to the winding. If the numerically integrated flux linkage value in the model becomes negative, then it should be set to zero:

$$
\Psi[m] = \begin{cases}
\,\,\,\Psi[m] \,\,\,\Psi[m] \ge 0; \\
0, \,\,\,\Psi[m] < 0.
\end{cases} \tag{15}
$$

The value of the current in the phase winding was calculated in accordance with the estimated flux linkage and the inductance value found from Equation (13):

$$i = \frac{\Psi}{L}.\tag{16}$$

If the evaluated current lies above the saturation knee, which can be checked by comparing it with the saturation current, then Equation (13) is not applicable and the phase current should be evaluated taking into account the differential inductance of the phase using the following equation:

$$i = i\_{\rm sat} + \frac{\Psi - L \cdot i\_{\rm sat}}{L\_{\rm min}} \tag{17}$$

instead of Equation (16), where *Lmin* is the minimum or differential inductance of the winding. If the actual current is smaller than the saturation current value, then the torque of a single phase was calculated using:

$$T = \frac{i^2}{2} \frac{\mathrm{d}L(\theta)}{\mathrm{d}\theta} = \frac{i^2}{2} \frac{\mathrm{d}(L\_{av} - \Delta L \cos\theta)}{\mathrm{d}\theta} = \frac{\Delta L}{2} i^2 \sin\theta. \tag{18}$$

*Energies* **2020**, *13*, 3317

The torque in the operation point above the saturation knee, in addition, was evaluated using co-energy, which can be expressed as:

$$\mathcal{W}'(\dot{i}, \boldsymbol{\Theta}) = \left(L\_{\text{av}} - \Delta L \cos \theta\right) \left(\frac{\dot{i}\_{\text{sat}}^2}{2} + i\_{\text{sat}}\left(\dot{i} - i\_{\text{sat}}\right)\right) + L\_{\text{min}}i\_{\text{sat}}\left(\dot{i} - i\_{\text{sat}}\right). \tag{19}$$

Then torque was found using:

$$T = \left. \frac{\partial \mathcal{W}'(i\_r \theta)}{\partial \theta} \right|\_i = \left( i\_{\rm sat} i - \frac{i\_{\rm sat}^2}{2} \right) \Delta L \sin \theta. \tag{20}$$

In this research, simulation was performed using C++ Builder. This tool was utilized in order to produce code for a control system suitable for further implementation using a microcontroller. The simulation was performed with the integration step size equal to 0.1 us. The parameters of the drive are listed in Table 1. As the control strategy contains only lookup table operation, it was supposed that the computation time is close to zero and that it is possible to evaluate the voltage commands for the next PWM cycle immediately at the end of the current PWM period.

**Table 1.** Switched reluctance drive (SRD) parameters.


The motor was running at 80 rad/s. The torque reference was changed during simulation as shown in Figure 8a. Initially it was set to 30 Nm, then after 30 PWM cycles it was set to 10 Nm, and finally to 45 Nm after 60 PWM cycles. The torque response transient time was strongly dependent on the current inductance of the winding. For example, at the beginning of the simulation, the torque slope was limited by the high unsaturated inductance of phase C. As the current reached the saturation knee, the torque on the shaft started to grow more quickly. Then, inside the saturated region, the phase currents and output torque were regulated more rapidly. At the end of commutation cycle of the phase, the current reference varied faster than the minimum possible current slope due to the voltage limit.

The torque stabilization is sufficiently precise over the entire range of the loads. It has visible deviations due to the PWM nature of the inverter, which can be minimized by increasing the PWM frequency. The same load profile was applied to the motor operating at a higher speed equal to 480 rad/s (see Figure 8b). With the growth of the electrical speed, the voltage limit made stabilization of the commanded current impossible. Therefore, torque pulsation appears. This behavior is natural for switched reluctance drives [17,18] and can be slightly improved by adjusting the current reference profile with respect to the voltage limit.

**Figure 8.** Operation at constant electrical speed: (**a**) 80 rad/s—mostly absent of voltage constraints; (**b**) 480 rad/s—operation with regular voltage constraints.

### **6. Conclusions**

Although the problem of precise torque control in switched reluctance motor drives was partially solved by means of finite control set model predictive control, the drive still suffers from an unpredictable switching rate and acoustic noise. The proposed solution implements a continuous control set MPC controlling power converter using pulse-width modulation and requires small computation efforts. It operates utilizing the assumption that optimal current reference profiles for each torque reference and angular position can be evaluated offline from the magnetization surface of the electrical machine. By knowing the current reference and magnetization surface, the voltage commands for the PWM-driven inverter can be evaluated using simple lookup tables, which has the advantage of practically zero computation delay and allows the calculation of the control law immediately at the end of each PWM cycle before applying the new duty cycles for the next PWM period.

The proposed CCS MPC can be implemented using modern microcontrollers. It was verified by a simulation model where accurate torque stabilization was achieved from zero to the rated speed. The settling time was limited by the supply voltage and the phase inductance.

Future work will be devoted to the practical implementation and verification of the proposed control strategy. In addition, it will address several questions that arose during this research.

Due to the current and torque slope limit when operating below the saturation knee, the proposed current reference profile is not optimal in terms of settling time. Nonetheless, it provides good efficiency of the motor and drive. Similar optimization of the current reference profiles is possible for operation at high speeds considering the voltage limit.

The sensitivity of the method to the inaccuracies of the model should be investigated, as well as the ability to use the current tracking error for identification of the rotor angular position for encoderless control.

**Author Contributions:** General idea, A.A.; Simulation software, A.Z. and A.A.; Simulation model verification, C.H.; Methodology, V.Š.; Simulation and data analysis, G.L.D. and A.B.; Validation, A.B.; Writing original draft, A.A. and G.L.D.; Writing review & editing, G.L.D., C.H. and V.Š. All authors have read and agreed to the published version of the manuscript.

**Funding:** This work was financially supported by Government of Russian Federation, Grant 08-08.

**Conflicts of Interest:** The authors declare no conflict of interest.

### **Nomenclature**


### **References**


© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

## *Article* **Additive Manufacturing and Performance of E-Type Transformer Core**

**Hans Tiismus 1,\*, Ants Kallaste 1, Anouar Belahcen 2, Anton Rassolkin 1, Toomas Vaimann <sup>1</sup> and Payam Shams Ghahfarokhi 1,3**


**Abstract:** Additive manufacturing of ferromagnetic materials for electrical machine applications is maturing. In this work, a full E-type transformer core is printed, characterized, and compared in terms of performance with a conventional Goss textured core. For facilitating a modular winding and eddy current loss reduction, the 3D printed core is assembled from four novel interlocking components, which structurally imitate the E-type core laminations. Both cores are compared at approximately their respective optimal working conditions, at identical magnetizing currents. Due to the superior magnetic properties of the Goss sheet conventional transformer core, 10% reduced efficiency (from 80.5% to 70.1%) and 34% lower power density (from 59 VA/kg to 39 VA/kg) of the printed transformer are identified at operating temperature. The first prototype transformer core demonstrates the state of the art and initial optimization step for further development of additively manufactured soft ferromagnetic components. Further optimization of both the 3D printed material and core design are proposed for obtaining higher electrical performance for AC applications.

**Keywords:** additive manufacturing; soft magnetic materials; selective laser melting; iron losses; magnetic properties; transformer

### **1. Introduction**

Metal additive manufacturing (AM) is maturing, enabling previously unavailable production possibilities in terms of feasible product complexity and personalization. As currently, the cost per part of AM is still relatively high, it has been most applicable for parts for high tech industries: producing specialized parts benefiting the most from the topology optimization possibilities of AM. For example, 3D printing has been utilized for the production of more efficient and long-lasting inductor coils [1], stronger, cheaper and lighter aircraft fuel nozzles [2], and high performance heat exchangers [3].

In parallel to the printing of structural, thermal, and electrical components, research interest in printed soft magnetic materials and topology optimized electromechanical components has spiked drastically over recent years. It has been proposed that with the easily available computational power and free-form printing capabilities of AM systems, next generation electrical machine designs could be modelled and constructed by the research community. These topology optimized designs (with reduced weight, integrated cooling channels, reduced inertia, increased heat exchange etc.) could be prototyped in-house, significantly reducing the lead time, cost, and machinery involved [4].

State of the art additive manufacturing of electromagnetic devices involves selective laser melting (SLM) printing of conductive and soft magnetic materials with air gaps

**Citation:** Tiismus, H.; Kallaste, A.; Belahcen, A.; Rassolkin, A.; Vaimann, T.; Shams Ghahfarokhi, P. Additive Manufacturing and Performance of E-Type Transformer Core. *Energies* **2021**, *14*, 3278. https://doi.org/ 10.3390/en14113278

Academic Editor: Salvatore Musumeci

Received: 5 May 2021 Accepted: 31 May 2021 Published: 3 June 2021

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2021 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

partitioning the material structure for separating individual turns in coils and reducing the induced eddy currents in soft magnetic cores [5,6]. The air gaps are printed due to the current lack of multi-material printing capacity of SLM systems, limiting the parallel printing of conductive, core, and insulation materials. The introduction of airgaps considerably reduces the power density of the components, however, as gapped printed component fill factor is typically relatively low (in the range of 60%) [6,7].

Despite extensive material optimization of different soft magnetic alloys, relatively few functional components or devices have actually been printed and characterized. For this reason, in this work, a full small-scale transformer core is printed, characterized, and compared with a commercial transformer. The simplistic design of an E-type transformer makes it ideal for the next step of testing additively manufactured magnetic material capacity and performance for electrical machine applications (succeeding the characterization of small-scale toroidal samples). In this paper, a novel interlocking core design is employed for eddy current reduction, which exhibits a competitive component fill factor. The paper is divided into two larger sections. The first part describes the 3D printed core design and its fabrication process, and the second the characterization and comparison of the printed core with conventional cores.

### **2. Transformer Core Design**

### *2.1. Commercial Transformer*

The 3D printed core design investigated in this paper was based on the commercially available 30 VA single phase isolation transformer provided by MS Balti Transformers Ltd (Tallinn, Estonia). The transformer was chosen based on its suitable size, type, and availability. Its shell-type transformer core is constructed from E-type stampings of grainoriented M 165-35S silicon steel. The conventional transformer design with its dimensions are detailed on Figure 1.

**Figure 1.** Investigated conventional transformer: (**a**) Core dimensions, (**b**) E-type stamping dimensions in detail.

The fully encapsulated modular windings of the transformer are utilized in both the conventional and 3D printed core designs. The modular windings are incorporated in both designs in order to improve the comparability of the transformer core performance and to demonstrate the compatibility of 3D printed and conventional parts. The nominal parameters of the windings are characterized in Table 1.

### *2.2. 3D Printed Design*

Next, an SLM printing system was utilized for the 3D printing of the full transformer core. The 3D printed core design was required to exhibit compatibility with the modular windings, incorporate the segregated structure for classical eddy current loss reduction with high filling factor, and adhere to the printing system requirements. Lamination thickness of 0.95 mm was chosen to obtain high fill factor and mechanical strength of the first prototype. For all segregated designs considered, it was critical to achieve continuous geometries (with minimal air gaps dividing the flux paths) with maximal flux path cross sectional area (high fill factor). Furthermore, since the printed transformer must be comprised of at least two parts (to accommodate the modular winding), optimization of the inter-part air gap must be considered. In conventional transformers, the influence of the inter-stamping airgaps is typically reduced by overlapping stamping layers: which facilitates the flux paths through the adjoining stampings. Similar overlap between the flux-guides can be realized in printed designs.


**Table 1.** Nominal parameters of the modular transformer coil.

For simplicity, in this paper, only conventional stamping inspired designs were considered for 3D printing. In Figure 2, three considered transformer core designs are illustrated: (a) a laterally laminated interlocking design from four parts, (b) an axially laminated gapped design from two parts, and (c) an axially laminated interlocking design from four parts. The axially laminated interlocking design was chosen for printing due to its simplicity and similarity to the conventional design, its high achievable fill factor and its post-processing possibilities: all of the unmelted powder can be removed between the laminations post-printing and, if needed, all of the surfaces can be cleaned and oxidized or varnished for enhanced inter-lamination electrical resistance.

**Figure 2.** Considered lamination strategies: (**a**) Laterally laminated interlocking design, (**b**) Axially laminated design with air-gapped core structure, (**c**) Exploded view of the 3D printed transformer core design with interlocking axial laminations comprising four individual components.

### **3. Methods**

### *3.1. Powder Characteristics*

Transformer parts were printed with identical powder, processing, and annealing parameters to the previous study characterizing the AC and DC losses of the printed material [6]. Pre-alloyed, gas-atomized Fe-Si provided by Sandvik group was utilized for printing. The powder exhibited roughly spherical particle shape with a median diameter of 38 μm, and its chemical composition is described in Table 2. The powder size, shape, and chemical composition were verified to verify the manufacturer declaration.

**Table 2.** Chemical composition of the employed Fe-Si powder.


### *3.2. SLM Printing of the Transformer Core*

Transformer core parts were printed on the SLM Solutions GmbH Realizer SLM-280. The printing system provides a 280 × 280 × 350 maximum build envelope and a single 1070 nm yttrium scanning laser (1 × 700 W). Custom smaller build platform (D100 mm) and re-coater were used for printing of the transformer core, designed for streamlining the powder substitution between projects for different raw powders.

Laser re-melting strategy was used to prevent the powder balling related uneven growth of the relatively large transformer parts during printing, which can result in rough porous material structure or the termination of the print job due to re-coater jamming. The phenomenon is related to an oxide film on the preceding layer impeding interlayer bonding and leading to balling, due to insufficient wetting of the molten metal on the oxide layer [8]. The balling phenomenon can be reduced in a higher purity environment (oxygen level below 0.1%), applying a combination of high laser powder and low scanning rate or applying re-melting scanning on the part [9].

Stripe (10mm wide) scan pattern was utilized with 30◦ rotation between layers. All of the printing was conducted in a nitrogen inert gas environment because of its relatively low cost. Platform pre-heating was not utilized as the custom reduced platform is not equipped for it. A summary of the main laser printing parameters is presented in Table 3.


**Table 3.** Summary of the printing parameters.

Transformer printing was completed in three parts in a total of 16 h: interlocking E-profiles separately (2 × 6 h) and the I-profiles in the same build (1 × 4 h). The printed components are illustrated in Figure 3: showing the surface finish, support structure, and the powder bed post-printing. Some concave warpage of the E-profiles was observed after separation from the build platform due to internal part stress, which obstructed the transformer assembly. Its causality can be traced to the relatively high internal stresses induced in part by the micro-welding process of SLM, and it can be resolved through the annealing of the printed parts at moderate temperature, pre-cutting from the platform for stress relief. Next, the support surfaces were polished and the inter-lamination air-gaps were lightly sanded for improved surface finish and fitting of the components.

**Figure 3.** SLM printed transformer components: (**a**) E-profile component post-printing, (**b**) I-profile component welded on the baseplate.

### *3.3. Annealing*

After mechanical post-processing, the printed transformer parts were annealed at 1200 ◦C in a low vacuum environment (~0.1 mBar) with a heating rate of 300 K/h, maintained at the target temperature for 1 h and then slowly furnace-cooled to room temperature.

### *3.4. Material Properties*

The additively manufactured 3.7% silicon steel shows comparable magnetic performance to non-oriented conventional silicon steels after thermal treatment. Magnetization of 1.5 T is achieved at 1800 A/m, exhibiting electrical resistivity of 56.9 μΩ·cm and hysteresis losses of 0.61 (W10,50) and 1.7 (W15,50) W/kg [6]. In comparison, a typical non-oriented steel M235-35A used for electrical machine fabrication exhibits total core losses of 0.92 (W10,50) and 2.35 W/kg (W16,50), resistivity of 59 μΩ·cm, and magnetization of 1.53 T at 2500 A/m. In this paper, we are comparing the additively manufactured core with a conventional Goss textured silicon steel M165-35S (equivalent to M111-35N) core, which shows superior magnetic properties to the non-oriented materials for transformer applications, as presented on Figure 4. The grain-oriented transformer steel shows approximately 0.3 T greater saturation magnetization than both of the non-oriented steels.

**Figure 4.** Magnetization curves of the studied materials: grain-oriented silicon steel M165-35S [10,11], non-oriented silicon steel M235-35A [12], printed annealed 3.7% silicon steel, and printed unannealed 3.7% silicon steel [6].

The grade designation of M165-35S of the Goss textured steel specifies 1.65 W/kg losses at 1.7 T (W17,50), and a lamination thickness of 0.35 mm. The materials' exact silicon content, resistivity, and other typical properties are unspecified and depend on the manufacturer (manufacturing freedom in the range of grade specifications).

### *3.5. Transformer Characterization*

The nominal performances of both the 3D printed and the conventional magnetic core transformers were characterized through open circuit and full load testing. The nominal parameters of the conventional transformer were obtained from the manufacturer's declaration. A drop in the nominal voltage is expected for the printed transformer due to its reduced fill factor, possible fitting defects (air-gaps between laminations), and lower saturation magnetization of the printed material. Its nominal voltage and iron losses were determined from the open circuit tests of the conventional transformer. To determine the transformer efficiencies, a load test was performed, where the transformer was energized up to nominal power. For thermal performance assessment, steady-state thermal images of the fully loaded transformers were captured with a Fluke Ti10 Thermal Camera.

The open circuit test setup is described in Figure 5, consisting of an autotransformer for variable voltage input and digital multimeters for measuring the voltage, current, and active power consumed in the transformer coil. In the open circuit test, the current drawn by the transformer establishes the magnetic field in the core. The active power consumed by the transformer signifies its total power loss, consisting mainly of magnetizing, and some ohmic, losses. The magnetizing losses summarize the energy lost from each magnetizing cycle, which are classically segregated into the hysteresis, classical, and excess eddy current loss.

**Figure 5.** Open circuit transformer: test setup (**a**) and its schematic (**b**).

The ohmic losses are induced from joule heating of the coils due to the magnetizing current drawn. The total specific transformer core losses can be calculated from (1), where *W* is the active power loss measured in the open circuit test, *I* is the magnetizing current, *R* is the magnetizing coil resistance, and *m* is the weight of the core.

$$P = \frac{\left(W - IR^2\right)}{m} \tag{1}$$

Traditionally, the magnetic material loss behaviour is discussed in terms of cycle peak polarization (*Bmax*) of the core. Unlike in the toroidal cores for magnetic material characterization [6,13], however, the flux density in the investigated transformer core can only be evaluated as an approximation, due to its uneven flux distribution. The analytical expression for calculating the peak polarization in a transformer can be derived from the differential form of Faraday's law (2), where *E* is the induced electromotive force by the switching magnetic field, *N* is the number of turns on the primary coil (1370), *f* is the excitation frequency of the magnetic field (50 Hz), *Bmax* is the peak material polarization, *S* is the core cross sectional area, *F* is the core filling factor, *U* is the applied voltage on the primary coil, and *Ur* is the voltage drop over the primary coil.

$$E = N \frac{d\Phi}{dt} \to E\_{\text{max}} = N 2\pi fSFB\_{\text{max}} \to \quad B\_{\text{max}} = \frac{E\_{\text{max}}}{N 2\pi fSF} = \frac{\mathcal{U} - \mathcal{U}\_r}{N 2\pi fSF} \tag{2}$$

Alternatively, the approximate material polarization can be evaluated from the material B-H curve (as presented on Figure 4) or by the finite element method (FEM) simulation. In both methods, the actual B-H curve of the transformer core can differ from the previously characterized material, most prominently due to air-gap related curve shearing. For *Bmax* evaluation, the magnetic field strength in the transformer is calculated from (3), where *N* is the number of turns on the primary coil, *i* is the peak magnetizing current and *l* is the length of the mean magnetic flux path of the core. All FEM simulations are performed in open source finite element analysis software package Finite Element Method Magnetics (FEMM). The model accounts for the transformer cross sectional geometry, magnetized up to the peak magnetizing current measured from the open circuit test, including the material magnetization curve and fill factor, but excluding any gaps in the core internal structure.

$$H = \frac{Ni}{l} \tag{3}$$

### **4. Results**

### *4.1. Assembled Transformer*

The conventional and finished assembled printed transformer cores are presented in Figure 6. The overall transformer core dimensions correlated well, with the printed transformer exhibiting a slightly thinner and lighter core. The fill factor of the 3D printed core was measured from the axial centerline of the interlocking E-cores. For the conventional transformer, the fill factor was adopted from the stamping datasheets. The physical comparison of the transformer cores is presented in Table 4. No additional oxidation, treatment, or varnishing was applied to the surfaces of the 3D printed transformer core for increased eddy current reduction—the insulation is provided by the high natural surface roughness of the printed parts.

**Figure 6.** Printed (**a**) and conventional (**b**) transformer cores.

**Table 4.** Physical comparison of the transformer cores.


### *4.2. Performance*

Open circuit tests of the transformers confirmed the flux drop in the core and the reduction of the sustainable operating voltage of the printed transformer. In Figure 7, both the magnetizing current drawn from the supply for generating the desired voltage and the iron loss behavior calculated from (1) are presented. At 40 mA magnetizing current, the conventional transformer is energized up to 230 V, while the printed transformer is energized to a 30% lower voltage of 160 V. This is due to the lower flux density sustained by the printed material. For energizing the printed transformer up to 230 V, a magnetizing current of 220 mA is required. This is inefficient, however, due to deep core oversaturation, requiring 450% more current than for magnetizing the conventional core and 30% more current than the rated full load current of the winding.

**Figure 7.** Magnetizing current drawn (**a**) and the specific core loss (**b**) of the tested transformers.

At 40 mA RMS excitation current (60 mA peak current), the analytically calculated (3) average H field generated in the magnetic core is 668 A/m, which corresponds to the magnetization of 1.72 T for M165-35S and 1.42 T for the annealed 3D printed material as determined from the magnetization curves in Figure 8. At 40 mA RMS excitation current, analogous excitation of both cores is achieved. Both are magnetized slightly above the approximate material knee-point and exhibit identical copper losses. Excitation of the conventional core to 160 V or the 3D printed core to 230 V would be impractical comparison-wise, as both states exhibit significantly differing magnetic behavior. At 160 V, the conventional core is still at the linear magnetic behavior: drawing only 7.6 mA magnetizing current and exhibiting 0.005 W of copper losses and 0.35 W of iron losses. At 230 V, the printed transformer shows deep saturation behavior, drawing 220 mA of magnetizing current, resulting in a significant voltage drop of 21.6 V, copper losses of 4.7 W, and iron losses of 2.6 W.

**Figure 8.** Core material magnetization curves correlated with the no load measurements of the investigated transformers.

FEM simulation of the transformer cores shows similar values of material magnetization: reaching 1.68 T for the conventional and 1.39 for the 3D printed core (Figure 9). Additionally, the simulation illustrates the uneven flux distribution in the core due to variations in transformer limb width. Analytical calculations with (2) show lower core flux density required for inducing a specific voltage in the core. For energizing the transformer up to 230 V, a flux density of 1.65 T is required, while for 160 V, a flux density of 1.26 T is required. The higher magnetization calculated from the experimental excitation current and FEM simulation is most likely the result of intra-lamination air-gaps, which shears the material magnetization curve and requires more current for achieving the same material polarization.

**Figure 9.** Flux distribution in the (**a**) conventional and (**b**) 3D printed transformer core.

Iron losses were identified as 1.82 W/kg for the conventional core at 230 V (at approximately 1.7 T, 50 Hz) and 3.05 W/kg for the 3D printed core at 160 V (in the range of 1.26–1.4 T, 50 Hz). Efficiency of the transformers was calculated from the load test measurements at both the ambient core temperature and the steady state temperature at full load conditions. The transformers reached steady state temperature after four hours of loading. The thermal images of the transformers are shown in Figure 10, with slightly higher heating observed for the 3D printed transformer core. The measured coil hotspot temperature was measured at 91.1 ◦C for the conventional core and at 95.1 ◦C for the 3D printed core. The core hotspots were measured with a thermocouple sensor due to the high reflectivity of the printed core, exhibiting temperatures of 71 ◦C (conventional) and 75 ◦C (3D printed).

At full load, the measured efficiency of the transformers ranged from 83.8% (21 ◦C) to 80.4% (71 ◦C) for the conventional transformer and 74.7% (21 ◦C) to 70.1% (75 ◦C) for the 3D printed transformer. The efficiency-load characteristic is presented in Figure 11. The highest efficiencies were measured at 41% load at ambient core temperature, reaching an efficiency of 88.7% for the conventional transformer and 80.5% for the 3D printed core. The efficiency of the 3D printed core was approximately 10% lower over the full measurement range. Due to the material saturation and inter-lamination air-gap related reduction of nominal voltage, the printed transformer core sustained reduced power density when compared to the conventional core. The transformer power density dropped 34% from 59 W/kg to 39 W/kg. The results of the transformer performance characterization are summarized in Table 5.

**Figure 10.** Steady state temperature of the studied transformers in the (**a**) conventional core and (**b**) 3D printed transformer core.

**Figure 11.** Efficiency-load characteristics of the studied transformers.

**Table 5.** Comparison of transformer performance.


### **5. Discussion**

The characterized transformers show typical performance values for small 20–30 VA power rating single-phase transformers. From manufacturer datasheets, the typical efficiency for a 30 VA rated power transformer is in the range of 83 [14]–81% [15], which decreases to 77% [14] at 22 VA and to 65% [14] at 4.5 VA. The rated power densities vary

significantly depending on the design (some designs are fully encased), and are typically in the range of 56 [15]–39 VA/kg [14] for 30 VA rated transformers and slightly lower (50 [16]–39 [14] VA/kg) for 20 VA rated transformers. In this study, we obtained an efficiency of 80.5% for the conventional transformer and 70.1% for the 3D printed transformer core at steady state temperature. The 10% reduced overall transformer efficiency can most prominently be attributed to the eddy currents generated in the 170% thicker laminations of the printed design. The reduced power density of the printed design can be attributed to both a larger degree of assembly defect related air-gaps within the core and the overall lower magnetic saturation of the printed material compared to the Goss textured conventional steel. Both designs are within the range of typical power density values for low power transformers.

The 3D printed core exhibited iron losses of 3.05 W/kg at 160 V transformer energization. Analytical calculations identify an average *Bmax* of 1.26 T at this transformer voltage level. Comparing the magnetizing values with previously measured 3D printed material magnetization curves, its shearing is proposed. Due to the air-gaps in the assembled printed design, more magnetizing current is required for the same material polarization and voltage generated by the transformer. Similar iron loss values have been measured by Plotkowski et al. for a 3D printed E-type transformer core [17]. In their work, they achieved a core loss of 3.5 W/kg (W10,60) at 1.0 T, 60 Hz magnetization for a printed 3% silicon steel lamination inspired core. They achieved considerably improved losses with more complex geometry, reaching approximately 1.5 W/kg (W10,60) to 3.2 W/kg (W15,60) with 'Hilbert pattern' 6% silicon steel. It is important to note, however, that in their work approximately 56% core fill factor was achieved, resulting in low power density and voltage generation of the transformer.

Further optimization of both the component topology and its material properties are unavoidable for achieving high performance 3D printed transformer cores. To obtain high magnetic polarization (high power density) of the printed material with minimal magnetomotive force, a higher degree of control of the printed material grain structure must be achieved. The effect of the grain structure orientation in relation to the magnetic field is significant as illustrated by Figure 12 [18]. In conventional stampings, the grainoriented pronounced Goss texture can be achieved with various hot and cold rolling stages of the steel sheets. In printed material, the optimization of the material grain structure is largely immature, with some grain structure evolution observed in [13], in heat treated laser-remelted printed silicon steel samples.

Several topological improvements can be applied to the transformer for enhanced performance. The printed transformer topology can be improved by increasing the fill factor of the assembled components, optimizing the lamination thickness for reduced eddy current loss, and increasing its power density through shape optimization for achieving uniform magnetization. Due to the limited multi-material printing capacity of current SLM systems, two methods are proposed for eddy current reduction: the interlocking and the gapped core designs. With next-generation powder deposition methods [19], multiple metal or intermetallic materials can be utilized in parallel, allowing for more options and more advanced core topologies.

First, for increasing the fill factor, higher accuracy of the printing system must be achieved. With the current settings, the printed parts still suffer from low surface roughnessrelated reduced fill factor for interlocking designs or inter-lamination short-circuits and sintered unremovable powder for the gapped designs. Secondly, the lamination thickness can be optimized to provide minimal core losses with maximum part fill factor, i.e., to achieve the optimal ratio of air gap to lamination width. Thirdly, the shape of the core can be optimized for achieving uniform magnetization, weight reduction and improved thermal capacity. Several methods for improving ferromagnetic part performance through topology optimization are discussed in further detail in [20,21]. For improved heat exchange of the printed transformer, enhanced convective heat transfer can easily be obtained by increasing its outer surface area with different surface relief structures [22].

### **6. Conclusions**

In this paper, a fully functioning, additively manufactured soft magnetic transformer core was fabricated and tested. For the first time in literature, an electromagnetic device with a fully 3D printed magnetic core was evaluated in terms of efficiency and performance. The prototype core showed uncompetitive performance when compared to modern conventional transformer cores. Although the printed material is not currently suitable for the production of commercial transformer cores, the analysis of the prototype core did allow us to demonstrate the current state of the art, identify the technical challenges involved, and propose next steps for realizing topology optimized 3D printing soft ferromagnetic components.

A novel, interlocking core design was developed and utilized successfully for achieving a relatively high fill factor of 89% (compared to other 3D printed cores) and eddy current reduction of the additively manufactured transformer core. For obtaining higher fill factor with this method, lower surface roughness of the printed parts must be obtained for more precise fitting of the components. Furthermore, the interlocking core design enabled the integration of modular winding to the transformer design, simplifying its assembly process.

The first prototype transformer core showed both lower efficiency (10% reduced) and power density (34% reduced), when compared to the conventional modern transformer at their respective optimal working conditions. These preliminary performance results of the first prototype core are likely to improve with more refined core designs and materials as part of future research. Currently, the main challenge in realizing high-performance 3D printed soft magnetic components is achieving a higher degree of control over the printed material grain texture, since the conventional post-processing methods for Goss textured silicon steel sheets are not suitable for processing geometrically complex 3D printed magnetic components. Even so, for non-grain-oriented applications (such as rotating electrical machines), the current material properties appear suitable, especially with the unprecedented prototyping freedom of 3D printing systems—which could enable the emergence of entirely new types of machines. Although the current 3D printed cores for AC applications suffer either from high eddy current losses or low filling factor, nextgeneration emerging multi-metal SLM printers can potentially improve the additively manufactured core performance considerably. Future work on this project will include further optimization of both the printed material and component topology for designing and constructing AM topology optimized electrical machines.

**Author Contributions:** Conceptualization: A.K. and H.T; methodology: A.R.; validation, A.K.; investigation, H.T. and P.S.G.; resources, P.S.G.; writing—original draft preparation, H.T.; writing review and editing, H.T. and A.K.; supervision, A.K. and A.B.; project administration, T.V.; funding acquisition, A.K. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research work has been supported by the Estonian Ministry of Education and Research (Project PSG-137).

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** Data is contained within the article.

**Acknowledgments:** The authors would like to thank Balti Transformers Ltd. for cooperation.

**Conflicts of Interest:** The authors declare no conflict of interests.

### **References**


## *Article* **Sliding Mean Value Subtraction-Based DC Drift Correction of B-H Curve for 3D-Printed Magnetic Materials**

**Bilal Asad 1,2,\*, Hans Tiismus 1, Toomas Vaimann 1, Anouar Belahcen 1,2, Ants Kallaste 1, Anton Rassõlkin <sup>1</sup> and Payam Shams Ghafarokhi <sup>3</sup>**


**Abstract:** This paper presents an algorithm to remove the DC drift from the *B-H* curve of an additively manufactured soft ferromagnetic material. The removal of DC drift from the magnetization curve is crucial for the accurate estimation of iron losses. The algorithm is based on the sliding mean value subtraction from each cycle of calculated magnetic flux density (*B*) signal. The sliding mean values (SMVs) are calculated using the convolution theorem, where a DC kernel with a length equal to the size of one cycle is convolved with B to recover the drifting signal. The results are based on the toroid measurements made by selective laser melting (SLM)-based 3D printing mechanism. The measurements taken at different flux density values show the effectiveness of the method.

**Keywords:** additive manufacturing; convolution; infinite impulse response (IIR) filters; additive white noise; DC drift; magnetic flux density; magnetic hysteresis; kernel; magnetic materials

### **1. Introduction**

Unlike subtractive and injection molding-based manufacturing techniques, additive manufacturing (AM) is gaining heightened popularity. It is also known as 3D printing, where any object can be built layer upon layer using any AM technique. The most common AM techniques include powder bed fusion (PBF), bed jetting (BJ), direct energy deposition (DED), material extrusion (ME), material jetting (MJ), sheet lamination (SL), vat polymerization (VP), etc. However, in electrical machines, direct laser melting (DLM) (a subcategory of PBF) has proven its effectiveness. Various materials such as metals, thermoplastics, ceramics, and biochemicals can be handled through AM techniques. The capability to design complex geometries, to create relatively lightweight products, to perform rapid prototyping, and to save time, as well as the lack of need for fixtures or dies, are the prominent advantages of this technology.

The growing AM-related technological advancements and complex geometries of electrical machines are persuading researchers to test 3D printing in this field. Electrical machines contain various materials such as copper, silicon steel (SiFe), insulation materials, etc. Although the printing of a complete machine in one set has not been achieved so far, the production of separate portions and their assembly afterward is used for various machines. The printing and assembly of motor parts, such as soft magnetic rotors [1,2], stators [3,4], electrical windings [5,6], bearings [7], heat exchangers [8,9], and insulations [5], can lead toward the production of any electrical machine.

Electrical machines are a combination of non-magnetic, with high electric conductivity, magnetic, with low electric conductivity, and dielectric materials, with no electric conductivity. These materials make the coils, cores, and insulation parts of electrical machines. The selection of any appropriate composite material with fixed proportions of

**Citation:** Asad, B.; Tiismus, H.; Vaimann, T.; Belahcen, A.; Kallaste, A.; Rassõlkin, A.; Ghafarokhi, P.S. Sliding Mean Value Subtraction-Based DC Drift Correction of B-H Curve for 3D-Printed Magnetic Materials. *Energies* **2021**, *14*, 284. https:// doi.org/10.3390/en14020284

Received: 3 December 2020 Accepted: 4 January 2021 Published: 6 January 2021

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2021 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

different elements is essential for better electrical machine performance. The composition of different materials can have a prominent effect on the core's magnetic properties, which can change the iron losses. The printing of several different materials using selective laser melting (SLM) is available in the literature. The SLM fabrication of conductive (Cu [10], AISi10MG [11]) and soft ferromagnetic (FeSi6.7 [12], FeSi6.9 [13], FeSi3 [14], Fe-Co-V [15], Fe-Ni-Si [16]) materials are well explained in the literature.

A slight variation in the proportion of the materials can produce dramatic changes in their magnetic characteristics. This makes the investigation of the *B*-*H* curve of the design material very important. The most common way to characterize any magnetic material is through its hysteresis loop measurement. Commercial hysteresis loop tracers are designed for traditional hard magnetic materials. For soft magnetic materials, custom measurement setups, designed according to the sample's particularities and the measurement objectives, are used. Unlike planar samples, closed circuits (e.g., toroid) are good choices to determine the material's magnetic characteristics. The toroidal structure is preferred because it has the least demagnetizing effects. Furthermore, they can be used to obtain magnetic anisotropy by torque magnetometry.

The biggest challenge in the measurement of an accurate hysteresis loop is the integrator DC drift. The leading causes of this drift include data acquisition setup offset values, thermal shift of parameters, impedance mismatch among various components of the measurement setup, and additive white gaussian noise (AWGN). This drift becomes very significant when many cycles need to be measured, such as toroidal samples where the maximum magnetic flux density is much larger in value than the coercivity [17]. The coercivity and remanence measurement depends upon the points where the hysteresis loop cut the *B* and *H* axis. Hence those points should be in a narrow range. Commercial integrators give a correction range up to 2 μWb/min [18], which is not enough for signals with a long measurement time. Some solutions to remove the DC drift are proposed in the literature, but they possess drawbacks, such as the following:


This work presents an algorithm to effectively remove the drift due to the DC offset of the data acquisition setup. This technique has not been previously presented in literature to the best of the authors' knowledge. Unlike constant and approximate polynomial data fitting functions, the proposed algorithm detects the drift at each sample, making it independent of signal length. For this purpose, the measured signals are convolved with a kernel function to get a sliding mean value function (SMVF). The subtraction of SMVF from corresponding measured signals reduces the DC drift significantly. Additionally, the high-frequency switching noise is reduced using infinite impulse response (IIR) low pass filters [22,23]. For the validation, a 3D-printed toroid is tested under different values of maximum flux density.

### **2. The Proposed Algorithm**

A flowchart diagram of the proposed algorithm is presented in Figure 1. A detailed explanation of all steps is provided below.

**Figure 1.** The algorithm flowchart.


**Figure 2.** The impact of sampling frequency on zero-crossing detection of flux density *B* measured at 50 Hz: (**a**) at a low sampling frequency (2 kHz), (**b**) at an improved sampling frequency (10 kHz).


$$H(t) = \frac{N\_1 i\_1(t)}{l} \tag{1}$$

$$B(t) = \frac{1}{N\_2 A\_c} \int\_{t\_i}^{t\_f} v\_2(t)d(t) \tag{2}$$

where *N*<sup>1</sup> and *N*<sup>2</sup> are the number of primary and secondary turns, respectively, *i*1(*t*) is the current in the primary winding, *l* is the average length of the toroid taken at the center, *Ac* is the core cross-sectional area, *v*2(*t*) is the secondary induced voltage, and *ti* and *tf* represent the measurement interval and signal length, respectively. It is important to know that the measured induced voltage contains the drift signal acting as a constant of integration.

The number of samples per cycle (SPC) is vital to define the length of the kernel function. For the known values of sampling (*Fs*) and supply frequency (*fs*), the size of the one cycle of *B* or *H* can be calculated as follows:

$$SPC = \frac{F\_s}{f\_s} \tag{3}$$

### *2.2. Defining Kernel Function and the Removal of the DC Drift*

Convolution is a potent tool in signal processing, where two functions generate a third function that describes how the shape of one signal is modified by the other. In digital signal processing, the convolution of two functions *f* and *g* can be defined as follows:

$$f(f\*g)[n] = \sum\_{m=-M}^{M} f[m]g[n-m] \tag{4}$$

where *f* is the input vector, which is *B*[*n*] in our case, and g is the kernel vector defined by (5), which acts as a filter. The size and shape of the kernel is the most significant part in signal processing. Since our goal is to find out the signal's envelope pattern, a constant kernel vector of size equal to one cycle of the input signal can be the right choice. Furthermore, to eliminate the effect of the kernel itself, the sum of its elements should be equal to unity, as shown below.

$$k[n] = \frac{ones(SPC, 1)}{SPC} \tag{5}$$

The convolution will give the moving mean value function across the entire input signal by making the kernel size equal to one cycle's length. A detailed description of the proposed convolution is presented in Figure 3. Figure 4 shows the rising trend of the envelope of *B*, which is recovered by convolving the flux density function (*B*[*n*]) with the proposed kernel function *k*[*n*]. The trend curve is brought to the signal's surface for ease of understanding, although it remains in the center of the signal as a mean value function. Subtracting the recovered mean value function from the original signal (*B*) removes the DC drift.

$$CB\_{l} = \frac{n\_{11}}{\nu \nu \nu\_{1}} + \frac{n\_{12}}{\nu \nu \nu\_{2}} + \dots + \frac{n\_{1n}}{\nu \nu \nu\_{n}}, \qquad CB\_{l+1} = \frac{n\_{12}}{\nu \nu \nu\_{1}} + \frac{n\_{13}}{\nu \nu \nu\_{2}} + \dots + \frac{n\_{21}}{\nu \nu \nu\_{n}}$$

$$CB\_{n} = \frac{n\_{n1}}{\nu \nu \nu\_{1}} + \frac{n\_{n2}}{\nu \nu \nu\_{2}} + \dots + \frac{n\_{nn}}{\nu \nu \nu\_{n}}. \qquad \qquad n \text{ ((cycle) ) ((sample)}$$

**Figure 3.** Description of the moving mean value technique.

**Figure 4.** The envelope of the flux density *B*(*t*) (**a**) with DC drift, shown by brown line recovered by the proposed algorithm; (**b**) the corrected signal after removing DC drift.

It is essential to know that when the kernel's overlapping size is not equal to the signal's cycle length, the non-overlapping elements of kernel functions are considered zero during the starting and ending intervals. This leads to the transient interval during the resultant signal's starting and ending points, represented by red dots in Figure 3. Since the length of those intervals is equal to one cycle's size, the mean value signal should be saved from *CBi* to *CBn*. The initial and last cycles can be discarded using zero-crossing detection.

### *2.3. Zero-Crossing Detection*

Almost all filters have transient and steady-state intervals. The transient interval appears at the beginning and the end of the signal. The duration of the transient interval depends upon the type of the filter and the windowing function. Since in the proposed algorithm, two filters, lowpass IIR and convolution, are used, the transient interval is very narrow. The resultant signal can be saved from the second zero crossing to the third last zero crossings to achieve the steady-state interval. The resultant acquired signal in the steady-state interval can have an integral or fractional number of cycles. Although the number of cycles does not significantly impact the *B-H* curve, they can affect the frequency resolution of fast Fourier transform (FFT) or other spectrum analysis-based signal processing techniques. Since each sinusoidal signal has three zero crossings, the integral number of cycles would have an odd number of zero crossings. This fact is used to save the signal in such a way that if the number of zero crossings is odd, then the signal is ready for further analysis; otherwise, the signal can be saved from first zero crossings until second last. Due to the ease in detecting the starting and ending fractional cycles, the method of zero-crossing detection for counting the integral number of cycles is preferred over the conventional approach where sampling frequency and the measurement length of the signal can be used for the same purpose.

### **3. Measurement Setup**

The experimental setup consists of two parts: printing and measurement. The sample toroid is prepared from Fe-Si powder with powder bed fusion printer Realizer SLM-50 (Germany). Laser re-melting strategy (each printed layer is scanned twice before applying the next layer of powder) is utilized for reducing the irregularities in the solidified part. The laser beam power was chosen as 50 W (1 m/s) for the primary and 75 W (0.75 m/s) for the secondary scan. The printed toroidal sample exhibited a 5 mm × 5 mm rectangular cross-section and a 60 mm outer diameter. The toroid was post-processed in a vacuum furnace for internal stress normalization and structure recrystallization. The toroid was heated 300 K/h up to the temperature of 1150 ◦C and was annealed in a vacuum chamber at 1150 ◦C for one hour and then furnace cooled. Q detailed description of the printing procedure is presented in [24,25].

The printed toroid was wound with primary (inner, magnetizing) and secondary coils (outer, measuring) with 150 (*N*1) and 50 (*N*2) turns. All the measurements were conducted per European standards EN 60404-4 and EN 60404-6 [26,27]. The primary coil was supplied with a sinusoidal current with different peak values to maintain flux density in the range of 0.5–1.6 T. The sinusoidal current was generated using an Omicron power amplifier CMS 356 (Austria), which was fed with a sinusoidal reference signal coming from a digital function generator, as shown in Figure 5. The current was measured using a precision resistor 75 mV/15 A connected in series with the primary winding. This current calculates magnetic field strength (*H*(*t*)) as in (1), while magnetic flux density (*B*(*t*)) is calculated from the output voltage across the secondary winding using (2). All these measurements were taken with a considerably high sampling frequency of 10 kHz using data acquisition setup Dewetron DEWE2-M (Austria). However, the measurements can be artificially improved using data interpolation, as discussed earlier.

**Figure 5.** Experimental setup of (**a**) the printing of the sample and (**b**) the measurement schematic diagram.

### **4. Results and Discussion**

The measured hysteresis loop with different maximum flux densities (*Bmax*) is shown in Figure 6. Figure 6a shows the results without any significant signal processing. As discussed earlier, the flux density and field strength vectors are calculated from the measured voltage and currents. Only the *B* vector is normalized across the zero line by subtracting it with its mean value. This is important; otherwise, different loops will have other locations due to significant DC shifts. The *B* and *H* intercept drift is visible by the lines' thickness as the drift due to DC offset is still there. Figure 6b shows the recovered loops after removing DC drift using 3rd-degree polynomial function. The polynomial function creates an approximate fitting function at the center of the vector. Hence, the drift can be removed by subtracting the *B*(*t*) from the polynomial. The drift is somewhat reduced, but still, the line widths are considerable. This is because of the approximation considered by the polynomial function. This problem worsens if the number of measurement cycles is less in number. Figure 6c shows the recovered hysteresis loops using the proposed algorithm. As compared to Figure 6a,b, the drift is considerably reduced in Figure 6c. Moreover, the proposed algorithm does not depend upon the number of cycles under consideration because it considers every cycle independently.

Figure 7 compares recovered loops using the conventional polynomial fitting function with ones obtained using the proposed algorithm. The zoomed windows in Figure 7a,b depict the *B-H* curve's sharpening using the proposed algorithm.

The approximate errors using different techniques are given in Table 1. This error is calculated by subtracting the higher (*B+h*) and lower intercept (*B+l*) values on the *B* axis at zero *H*. The error is considerably reduced from 0.05 T (intercept width) to 0.0002 T using the proposed algorithm.

**Figure 6.** The measured *B-H* curves at different maximum flux densities *Bmax* (**a**) without DC drift correction; (**b**) DC drift correction using the polynomial fitting function as in [17]; (**c**) the corrected *B-H* curves using the proposed algorithm.

**Figure 7.** A comparison of the different approaches, (**a**) with zoomed window near zero crossing of the *B* axis, (**b**) with zoomed window at the peak value.

**Table 1.** The comparison of different drift removal techniques in terms of *B* intercept width.


### **5. Conclusions**

Advancements in SLM-based AM techniques have opened a broad domain of intricate electrical machine designs. This complexity can be in the form of optimized complex mechanical geometry or the selection of different types of composite material for efficiency improvement. The selection of composite material to be used for fabrication purposes depends upon its magnetic properties, which directly influence efficiency. Hence, before the final selection of machine-fabricated material, an evaluation of its characteristics is mandatory. The percentage content of different elements in the fabricating power can be optimized based on evaluations of those characteristic. In electrical machines, the material's magnetic characteristics are the most crucial parameter among several others. For this purpose, an algorithm was proposed in this paper, which was shown to produce the *B-H* curve of the 3D printed sample with less error. The proposed model was shown to have the following benefits:


**Author Contributions:** Conceptualization: B.A., H.T., T.V., and A.B.; methodology: B.A., H.T., T.V., and A.B.; validation: A.K., A.R., and A.B.; data curation: B.A. and H.T.; writing—original draft preparation: B.A.; writing—review and editing: B.A., H.T., and A.R.; visualization: T.V. and P.S.G.; supervision: T.V., A.K., and A.B. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research work has been supported by Estonian Ministry of Education and Research (Project PSG137).

**Conflicts of Interest:** The authors declare no conflict of interest.

### **References**

