Article navigation
Purpose

This paper aims to address the limitations of existing transformer core models, which often lack applicability in real-time control tasks. The paper investigates properties of a control-oriented transformer core model that is shown to accurately represent the effects of hysteresis, saturation, eddy current losses and leakage flux while at the same time minimizing the computational burden.

Design/methodology/approach

The investigated control-oriented model is mathematically expressed in a state-space form, making it inherently suitable for model-based controller design. This model can be interpreted as an equivalent electric diagram, providing a clear and intuitive representation. The model is experimentally evaluated across various transformers which are selected based on an extensive data analysis aimed at identifying practically relevant representatives. The k-means clustering algorithm is used to ensure that the transformers exhibit fundamentally different characteristics.

Findings

The investigated control-oriented model successfully fulfils its purpose. It effectively replicates the transient current response of various current transformers and generalizes well across a wide range of input signals with minimal computational effort.

Originality/value

Comparing the proposed control-oriented model with the Jiles–Atherton model demonstrates its effectiveness in terms of simulation accuracy and computational efficiency. Therefore, it offers a practical solution for problems commonly found in real-time control applications.

A major challenge of modelling inductive transformers is capturing the dynamics of the ferromagnetic core, which includes material-specific effects like saturation, hysteresis, and eddy current losses. According to Mörée and Leijon (2023), the two primary approaches for modelling ferromagnetic hysteresis are Duhem-type and Preisach-based models.

Duhem-type models use an ordinary differential equation (ODE), with the Jiles–Atherton model being a prominent implementation (Jiles and Atherton, 1986). However, evaluating the Jiles–Atherton equations in the presence of measurement noise poses significant numerical challenges.

Preisach-based hysteresis models (Preisach, 1935) approximate the behaviour of ferromagnetic materials through the superposition of numerous independent hysteresis elements. In general, Preisach-based models need an extensive amount of parameters, making their identification challenging. Thus, it is often regarded as a data-driven approach.

A different approach modelling the hysteresis of ferromagnetic materials is the so-called dynamic hysteresis model. It is introduced in Zirka et al. (2014) and is based on the reconstruction of first-order reversal curves. As it describes the behaviour of the hysteresis with the magnetic flux density as the input and the magnetic field strength as the output, the model is particularly useful implementation in circuit simulators and transient simulation (Zirka et al., 2015). The model has been successfully applied to simulate the ferromagnetic losses of a laminated core under PWM excitation in Simulink (Rasilo et al., 2019). This model is considered to be an “equation-free” model. It is implemented as an algorithm instead of equations – which makes the model inconvenient to be applied for control tasks.

Model-based controller design methods often benefit from using models in state-space representation (Khalil, 2014). However, neither the Duhem-type model nor the Preisach model can be expressed in this format, highlighting the need for control-oriented models.

A hysteresis model which can be formulated in a state-space form is the Chua–Stromsmoe model (Chua and Stromsmoe, 1970). An implementation of this approach is described in Ranta et al. (2011). The model is constructed by finding two functions describing the mid-curve and the width of a measured hysteresis loop. Thus, the Chua–Stromsmoe hysteresis model is only applicable for the input waveform that was used to generate the measurement and identifying the model.

In Schwartze et al. (2025a, 2025b), a control-oriented model is introduced which combines both advantages: it can be written in a state-space representation and is valid for arbitrary input signal forms. This model lays the basis for the present paper, which focuses on a more empirical investigation of the control-oriented model with practical applications in mind. The main contributions are as follows:

  • A data analysis based on a clustering approach for identifying fundamentally different transformers for a comprehensive model validation.

  • The experimental model validation based on the previously chosen transformers exhibiting distinct characteristics.

  • A runtime experiment comparing the computational effort between the control-oriented model and the Jiles–Atherton model.

  • An investigation concerning the numerical stability of the simulations.

In Section 2, an overview of the control-oriented model described in Schwartze et al. ( 2025b) is given. To validate the mathematical model on a wide range of possible dynamics while minimizing the testing effort, a data analysis is performed in Section 3 and Section 4. In Section 5, the control-oriented model and a model based on the Jiles–Atherton equations are compared to measurements. A performance analysis regarding the computational effort of both models is conducted in Section 6. A conclusion is given in Section 7.

This section gives an overview of the mathematical model derived in Schwartze et al. ( 2025b). The proposed control-oriented model describes one side of a single-phase transformer (i.e. either the primary or the secondary side). This configuration is especially important for model-based testing of transformers, which is the motivation and application considered in this article. It is assumed that a voltage u acts as the input to the system and the output current i⁠, on the same side, is regarded as the output of the system. Thus, the model reproduces the dynamics of the inductance seen from one side of the transformer. The corresponding equivalent circuit diagram is depicted in Figure 1.

Figure 1.
A schematic diagram illustrates an electrical circuit model with resistances, inductances, flux paths, and a transformer section.The diagram depicts an electrical equivalent circuit beginning at input terminals with voltage u of t and current i of t. A winding resistance R w is in series with a leakage inductance L k. The circuit then splits into two parallel branches. The first branch contains a core loss resistance R c. The second branch contains a magnetizing inductance L h in series with a hysteresis loss resistance R h, and a nonlinear element s of phi m. Arrows indicate magnetic flux directions, with phi m for the main flux and phi h for the hysteresis flux. On the far right, the circuit connects to a partial transformer representation.

Equivalent circuit diagram (Schwartze et al.,  2025b) representing the lumped parameters of the control-oriented model

Note(s): The input voltage u is applied to the system, which invokes a corresponding output current i⁠. The ohmic winding resistance is denoted by Rw and Lk describes the leakage flux. The parameters Rh and Lh capture the hysteresis effect, while Re is used to model the eddy current losses. The function sϕm⁠, where ϕm denotes the magnetic flux present in the core at any time instance, introduces the saturation. The primary side of the transformer is not included

Source: Authors’ own work and also compared to Schwartze et al. ( 2025b)

Figure 1.
A schematic diagram illustrates an electrical circuit model with resistances, inductances, flux paths, and a transformer section.The diagram depicts an electrical equivalent circuit beginning at input terminals with voltage u of t and current i of t. A winding resistance R w is in series with a leakage inductance L k. The circuit then splits into two parallel branches. The first branch contains a core loss resistance R c. The second branch contains a magnetizing inductance L h in series with a hysteresis loss resistance R h, and a nonlinear element s of phi m. Arrows indicate magnetic flux directions, with phi m for the main flux and phi h for the hysteresis flux. On the far right, the circuit connects to a partial transformer representation.

Equivalent circuit diagram (Schwartze et al.,  2025b) representing the lumped parameters of the control-oriented model

Note(s): The input voltage u is applied to the system, which invokes a corresponding output current i⁠. The ohmic winding resistance is denoted by Rw and Lk describes the leakage flux. The parameters Rh and Lh capture the hysteresis effect, while Re is used to model the eddy current losses. The function sϕm⁠, where ϕm denotes the magnetic flux present in the core at any time instance, introduces the saturation. The primary side of the transformer is not included

Source: Authors’ own work and also compared to Schwartze et al. ( 2025b)

Close Figure 1.

As the spatial resolution of the relevant field quantities is not important in this context of model-based controller design, a lumped element interpretation can be used.

The proposed control-oriented model consists of three coupled nonlinear ODE, which can be represented in a state-space form (Khalil, 2014) denoted by the following:

(1a)
(1b)

where u is the scalar input, y is the scalar output and x represents the state vector:

(2)

of the system.

According to the electrical equivalent diagram presented in Figure 1, the states can be interpreted in a physically meaningful way. The first state x1= ϕm denotes the magnetic flux that is present in the ferromagnetic core at any given time instance. The second state x2= ϕh represents a part of this flux, that is responsible for the width of the hysteresis. The output current that flows through the coil is x3=i⁠.

As derived in Schwartze et al. ( 2025b), the right-hand side fx,u of the ODE presented in (1a), is as follows:

(3)

where Rw, Lk, Re, Rh, Lh∈ R+ denote the model parameter. Each of these parameters belong to a corresponding physical effect. The ohmic resistance of the transformer windings is represented by Rw⁠. A common approach to model the leakage flux, i.e. magnetic flux that is induced by the coil but not applied to the ferromagnetic material, is to include a leakage inductance Lk in series to the main coil. A constant core-loss resistor Re is used to represent the eddy current losses in the ferromagnetic core. The combination of Rw and Lh is assumed to capture the material-specific hysteresis losses. For a more detailed theoretical analysis, see Schwartze et al. ( 2025b).

As the magnetic flux increases, the ferromagnetic material saturates, which is captured by the function sx1⁠. Theoretically, any appropriate function that best describes the saturation characteristic of the ferromagnetic material can be used, a comprehensive overview is given in Fischer and Moser (1956). The function:

(4)

presented in Ranta et al. (2011) is found to give an reasonable trade-off between complexity and adaptability to a wide range of different magnetic materials. The positive parameters β, S, Lm∈ R+ can be selected to adjust the saturation characteristic. Thus, the proposed model (1) consists of 8 constant parameters, which can be identified through a combination of measurements and optimization, which is explained in more detail in Schwartze et al. ( 2025b).

The characteristics of transformers can be very diverse, depending on their shape, size, ferromagnetic core material and application. To validate the model across a broad range of devices, a data set comprising measurements from 2,000 current transformers (CTs) is analysed.

A standard method for evaluating the operation and accuracy class of CTs is the model-based test (Jäger et al., 2014). A device capable of performing such a test is the CT Analyzer (OMICRON electronics GmbH, 2022). During this test, characteristics of the connected CT are identified. Some important parameters are as follows:

  • unsaturated inductance Lm⁠: inductance of the CT below the saturation of the material;

  • saturation flux ϕm,s⁠: material-specific limit for magnetic flux that can be present in the core;

  • remanence factor Kr⁠: ratio of remanence flux to saturation flux;

  • winding resistance Rw⁠: ohmic resistance of the windings in the coil; and

  • eddy current resistor Re⁠: constant core resistor representing the power losses due to the eddy currents induced in the ferromagnetic core.

These parameters are identified purely based on measurements and there is no need for any user input during the procedure. Further, they are connected to physical quantities of the CT such as shape, size and ferromagnetic material. These parameters are assumed to uniquely describe a CT and are therefore selected to form the feature-space for the clustering algorithm.

To identify groups, i.e. clusters, of CTs within the data set, an unsupervised clustering algorithm is used. As described in (Xu and Tian, 2015), selecting a suitable clustering algorithm is not trivial and strongly depends on the data set and the problem at hand. It can be observed that the data points are spread out continuously through the search space which leads to the assumption that the underlying clusters might be convex sets. The present data set is large with a smaller dimensional search space. For this clustering problem, the K-means algorithm (Thorndike, 1953) is an appropriate choice (Xu and Tian, 2015). The K-means algorithm is a well-established method that is computational efficient compared to other clustering algorithms. This property allows multiple reinitializations of the clustering process with reasonable time consumption. The K-means algorithm partitions the data set into k distinct clusters through an iterative process. The data points are assigned to the nearest cluster centroid. In a second step the centroids are updated based on the mean of the assigned points. This process continues until convergence is achieved. For this analysis, the scikit-learn implementation of the K-means algorithm (Pedregosa et al., 2011) is used in Python 3.8.13. The quality of the clustering heavily depends on the initialization of the starting values selected in the first iteration, here the “k-means++” initialization method is used. The algorithm is reevaluated with different starting points for 1,000 times, before the best clustering result is accepted.

The K-means algorithm needs to know the expected number of clusters ka priori. Different methods for selecting the k parameter are described and compared in Yuan and Yang (2019). A simple, yet powerful method is the “elbow-method”. There, the sum of the inertia of all clusters is calculated by the following:

(5)

where nj denotes the number of CT measurements per cluster j⁠. The variable dj,i represents the ith data point belonging to the jth cluster with the centroid μj⁠. This value J is calculated for clustering results with varying values k⁠, which leads to the elbow plot depicted in Figure 2. The elbow-method is a visual approach, with the goal of choosing an appropriate number of clusters k that gives a good clustering result (i.e. a small value for J⁠) while at the same time preventing over-fitting. The elbow-method is also criticized for being a subjective measure, as no formal definition of the “elbow-point” exists. Therefore, no statement about the optimality for the value of k can be drawn from this method. A detailed analysis of the elbow-method and other alternative approaches is discussed in Schubert (2023). While the elbow-method can be criticized, it is frequently utilized because of its intuitive simplicity and demonstrated effectiveness in practical applications. According to Figure 2, k=7 is selected for the considered data set. This marks the first point where the inertia starts to decrease approximately linearly with increasing number of clusters k⁠.

Figure 2.
A graph depicting the relationship between inertia and k values, showing a downward curve with points plotted. One point is highlighted in red.The image displays a graph indicating the relationship between inertia, represented on the vertical axis, and k values on the horizontal axis. The vertical axis ranges from zero to approximately one hundred and thirty, while the horizontal axis extends from zero to about eighteen. The plotted data points create a downward curve, demonstrating a decrease in inertia as k increases. One specific data point is highlighted in red, indicating its particular significance within the dataset, while the rest of the points are shown in blue. The graph includes an axis grid for better readability and analysis.

The elbow-method is used to select the number of clusters for the K-means algorithm

Note(s): As no formal definition of the “elbow-point” exists, the method is considered to be subjective. Nevertheless, due to the simplicity and intuition, this approach is often applied. The elbow-method selects k by finding a compromise between the quality of the clustering and over-fitting. The value is chosen to be k=7⁠, resulting in 7 distinct clusters

Source: Authors’ own work

Figure 2.
A graph depicting the relationship between inertia and k values, showing a downward curve with points plotted. One point is highlighted in red.The image displays a graph indicating the relationship between inertia, represented on the vertical axis, and k values on the horizontal axis. The vertical axis ranges from zero to approximately one hundred and thirty, while the horizontal axis extends from zero to about eighteen. The plotted data points create a downward curve, demonstrating a decrease in inertia as k increases. One specific data point is highlighted in red, indicating its particular significance within the dataset, while the rest of the points are shown in blue. The graph includes an axis grid for better readability and analysis.

The elbow-method is used to select the number of clusters for the K-means algorithm

Note(s): As no formal definition of the “elbow-point” exists, the method is considered to be subjective. Nevertheless, due to the simplicity and intuition, this approach is often applied. The elbow-method selects k by finding a compromise between the quality of the clustering and over-fitting. The value is chosen to be k=7⁠, resulting in 7 distinct clusters

Source: Authors’ own work

Close Figure 2.

A subspace of the full feature-space depicting the clusters is shown in Figure 3. It is demonstrated that the data set can be visually separated in the feature-space consisting of saturation flux ϕm,s⁠, unsaturated inductance Lm and remanence factor Kr⁠. This analysis allows to assign some characteristics to each identified group. It is assumed that by verifying the model for the centroid of the cluster, the model will work reasonably well to describe all CTs within this cluster. Three representative CTs of the clusters C0⁠, C2 and C3 are selected for the comparison presented in Section 5 and are described in more detail in Table 1.

Table 1.

Current transformer clustering selection

ClusterColorSymbolDescription
C0Purple◊Material with higher relative permeability μr and lower saturation flux ϕm,s
C2Light blueXMaterial with lower remanence factor Kr⁠, typically including transformers with thinner hysteresis
C3Turquoise□Larger cores, material with higher saturation characteristic ϕm,s⁠, thin hysteresis
Source(s): Authors’ own work
Figure 3.
A 3 D scatter plot displays data points categorised into seven classes, plotted across saturation flux, unsaturated inductance, and remanence factor.The plot presents a three-dimensional scatter chart with the horizontal axis for saturation flux in webers, the depth axis for unsaturated inductance in henries, and the vertical axis for remanence factor in percent. Data points are grouped into seven distinct classes, labeled from class 0 to class 6, each with a unique marker style. The distribution shows that certain classes cluster at lower inductance and higher remanence factor values, while others extend across higher saturation flux and inductance ranges. The spatial separation of these clusters indicates variations in the underlying properties represented by each class.

Identified seven clusters which are separable in the feature-space: saturation flux ϕm,s⁠, unsaturated inductance Lm⁠, remanence factor Kr

Note(s): Transformer-specific characteristics can be generalized for each group

Source: Authors’ own work

Figure 3.
A 3 D scatter plot displays data points categorised into seven classes, plotted across saturation flux, unsaturated inductance, and remanence factor.The plot presents a three-dimensional scatter chart with the horizontal axis for saturation flux in webers, the depth axis for unsaturated inductance in henries, and the vertical axis for remanence factor in percent. Data points are grouped into seven distinct classes, labeled from class 0 to class 6, each with a unique marker style. The distribution shows that certain classes cluster at lower inductance and higher remanence factor values, while others extend across higher saturation flux and inductance ranges. The spatial separation of these clusters indicates variations in the underlying properties represented by each class.

Identified seven clusters which are separable in the feature-space: saturation flux ϕm,s⁠, unsaturated inductance Lm⁠, remanence factor Kr

Note(s): Transformer-specific characteristics can be generalized for each group

Source: Authors’ own work

Close Figure 3.

The model’s performance is assessed by comparing measurement and simulation over one period, where the voltage u is supplied as an input and the output current i is measured.

Before the control-oriented model can be simulated, its parameters need to be identified during a separate step. This is done based on a combination of direct measurements for the parameters:

  • winding resistance Rw⁠;

  • eddy current losses Re⁠; and

  • leakage inductance Lk⁠.

and an optimization procedure for the parameters:

  • hysteresis losses Rh and Lh⁠; and

  • saturation characteristic sϕm⁠.

The identified parameters of the control-oriented model for the CTs from the clusters C0⁠, C2 and C3 are provided in Table 2.

Table 2.

Parameters of the identified control-oriented model

Control-oriented model parametersTransformer C0Transformer C2Transformer C3
Rw3.0370.41872.8255
Re2.3632e + 47.7384e-350e + 4
Lk0.84e-70.18e-68.3725e-4
Rh2.4402e + 47.1468e + 33.4538e + 03
Lh10.28685.35471.67e-2
Lu4.4452e + 31002.6315
S208.862.0002
β0.2230.03555.0719
Source(s): Authors’ own work

The control-oriented model is also compared to the well-established Jiles–Atherton model, which acts as a baseline and is described in more detail in Albert et al. (2024) and consists of 11 parameters. The ohmic resistance of the windings on the secondary side is reflected in the parameter Rsec⁠. The parameters L12σ and R12σ describe the leakage flux between the primary and the secondary side. The geometry of the transformer is modelled by the parameters Ac⁠, i.e. the cross-sectional area, and LFe⁠, the length of the ferromagnetic core. RFe describes the eddy current losses and corresponds to Re of the control-oriented model. All other parameters are the standard parameters of the Jiles–Atherton model where Ms is the magnetization saturation, a the domain wall density, α the interdomain coupling, k the pinning force and c the magnetization reversibility. The identified parameters of the Jiles–Atherton model are shown in Table 3. For a fair comparison, both, the control-oriented model and the Jiles–Atherton model, are implemented in MATLAB/Simulink 2023a (The MathWorks Inc, 2023a).

Table 3.

Parameters of the identified Jiles–Atherton model

Jiles–Atherton parametersTransformer C0Transformer C2Transformer C3
Rsec3.0370.41872.8255
R12σ13.42e-621.79e-62.5e-5
L12σ0.84e-70.18e-68.3725e-4
Ac1.5e-42.5e-425e-4
LFe0.125660.471240.8324
RFe2.3632e + 47.7384e-350e + 4
Ms26.5e + 55.6674e + 517.3e + 5
a2104.3433499
α2.3581e-42e-52.0e-8
k342.97792.1949e-8
c0.04090.081.2222e-6
Source(s): Authors’ own work

To show the general applicability of the control-oriented model, the experiment is performed over a wide range of different input voltage settings uo⁠. Here, o∈1,2,3,4,5,6 denotes the operation point:

  • u1=58 V square wave at 50 Hz;

  • u2=56 V sine wave at 53 Hz;

  • u3=8.75 V square wave at 50 Hz;

  • u4=14 V square wave at 50 Hz;

  • u5=118 V square wave at 5.6 Hz; and

  • u6=118 V square wave at 7.3 Hz.

Although the input voltage settings cover a large variety, no real statement can be made outside of these experiments. From Schwartze et al. (2025b), it is known that the model is frequency dependent. The exact frequency band within which the model reproduces the measurement up to a certain allowed error is derived and experimentally validated in Schwartze et al. ( 2025a).

This assessment is conducted on the CTs representative for the clusters C0⁠, C2 and C3⁠. While these CTs were deliberately selected to encompass a broad range real and technically feasible CTs, completeness cannot be guaranteed. CTs with uncommon or atypical characteristics that occur infrequently in the data set do not form distinct clusters are therefore not covered by this evaluation.

The qualitative visual comparison of the resulting current response and the corresponding hysteresis for the three representative CTs can be seen in Figure 4, Figure 5 and Figure 6.

Figure 4.
Four graphs compare current and magnetic flux for square and sine waves, each showing measurement, proposed model, and J A model over time and current.The image contains four graphs in a two-by-two arrangement. The top left graph plots current in amperes against time in seconds for a square wave at fifty-eight volts and fifty hertz, while the top right graph shows magnetic flux in webers versus current for the same square wave. The bottom left graph plots current against time for a sine wave at fifty-six volts and fifty-three hertz, and the bottom right graph shows magnetic flux versus current under the same sine wave conditions. Each graph features three lines: a solid line for measurement, a dashed line for the proposed model, and a dotted line for the J A model. The axes are clearly labeled with units, and each graph includes annotations specifying wave type, voltage, and frequency, enabling direct comparison of the models under the two waveform conditions.

Comparison of the current response of the control-oriented model with measured data and the Jiles–Atherton model on the representative CT of cluster C0

Note(s): The input voltages u1 and u2 are used

Source: Authors’ own work

Figure 4.
Four graphs compare current and magnetic flux for square and sine waves, each showing measurement, proposed model, and J A model over time and current.The image contains four graphs in a two-by-two arrangement. The top left graph plots current in amperes against time in seconds for a square wave at fifty-eight volts and fifty hertz, while the top right graph shows magnetic flux in webers versus current for the same square wave. The bottom left graph plots current against time for a sine wave at fifty-six volts and fifty-three hertz, and the bottom right graph shows magnetic flux versus current under the same sine wave conditions. Each graph features three lines: a solid line for measurement, a dashed line for the proposed model, and a dotted line for the J A model. The axes are clearly labeled with units, and each graph includes annotations specifying wave type, voltage, and frequency, enabling direct comparison of the models under the two waveform conditions.

Comparison of the current response of the control-oriented model with measured data and the Jiles–Atherton model on the representative CT of cluster C0

Note(s): The input voltages u1 and u2 are used

Source: Authors’ own work

Close Figure 4.
Figure 5.
Four graphs compare current and magnetic flux at two voltage conditions, showing measurement data alongside proposed and J A model predictions.The image contains four graphs arranged in a two-by-two grid. The top-left graph plots current in amperes against time in seconds for a voltage of 8.75 volts at 50 hertz, with a solid line for measurement, a dashed line for the proposed model, and a dotted line for the J A model. The top-right graph shows magnetic flux in webers versus current in amperes for the same voltage condition. The bottom-left graph mirrors the top-left but for 14 volts at 50 hertz, while the bottom-right graph presents magnetic flux versus current for this 14 volt case. Identical axis ranges and units are used across graphs to allow direct visual comparison of measurement and model outputs for both voltage conditions.

Comparison of the current response of the control-oriented model with measured data and the Jiles–Atherton model on the representative CT of cluster C2

Note(s): The input voltages u3 and u4 are used

Source: Authors’ own work

Figure 5.
Four graphs compare current and magnetic flux at two voltage conditions, showing measurement data alongside proposed and J A model predictions.The image contains four graphs arranged in a two-by-two grid. The top-left graph plots current in amperes against time in seconds for a voltage of 8.75 volts at 50 hertz, with a solid line for measurement, a dashed line for the proposed model, and a dotted line for the J A model. The top-right graph shows magnetic flux in webers versus current in amperes for the same voltage condition. The bottom-left graph mirrors the top-left but for 14 volts at 50 hertz, while the bottom-right graph presents magnetic flux versus current for this 14 volt case. Identical axis ranges and units are used across graphs to allow direct visual comparison of measurement and model outputs for both voltage conditions.

Comparison of the current response of the control-oriented model with measured data and the Jiles–Atherton model on the representative CT of cluster C2

Note(s): The input voltages u3 and u4 are used

Source: Authors’ own work

Close Figure 5.
Figure 6.
Four graphs compare current and magnetic flux over time and against current for different excitation conditions, with measurement, proposed model, and J A model distinguished by line styles.The image contains four graphs arranged in a two-by-two format. The top left graph shows current in amperes against time in seconds for square wave excitation at 118 volts and 5.6 hertz, with a solid line for measurement, a dashed line for the proposed model, and a dotted line for the J A model. The top right graph presents magnetic flux in webers versus current in amperes for the same excitation conditions. The bottom left graph depicts current in amperes against time in seconds for square wave excitation at 118 volts and 7.3 hertz, using the same line style coding. The bottom right graph plots magnetic flux in webers versus current in amperes for sine wave excitation at 118 volts and 7.3 hertz. All graphs maintain consistent axes and legend formatting, enabling direct comparison between measured and modelled data.

Comparison of the current response of the control-oriented model with measured data and the Jiles–Atherton model on the representative CT of cluster C3

Note(s): The input voltages u5 and u6 are used

Source: Authors’ own work

Figure 6.
Four graphs compare current and magnetic flux over time and against current for different excitation conditions, with measurement, proposed model, and J A model distinguished by line styles.The image contains four graphs arranged in a two-by-two format. The top left graph shows current in amperes against time in seconds for square wave excitation at 118 volts and 5.6 hertz, with a solid line for measurement, a dashed line for the proposed model, and a dotted line for the J A model. The top right graph presents magnetic flux in webers versus current in amperes for the same excitation conditions. The bottom left graph depicts current in amperes against time in seconds for square wave excitation at 118 volts and 7.3 hertz, using the same line style coding. The bottom right graph plots magnetic flux in webers versus current in amperes for sine wave excitation at 118 volts and 7.3 hertz. All graphs maintain consistent axes and legend formatting, enabling direct comparison between measured and modelled data.

Comparison of the current response of the control-oriented model with measured data and the Jiles–Atherton model on the representative CT of cluster C3

Note(s): The input voltages u5 and u6 are used

Source: Authors’ own work

Close Figure 6.

Both models are able to reproduce the measurements, while the control-oriented model is better at generalizing over a wider range of input settings. This is especially visible for CTC0 in Figure 4, where the sinusoidal input u2 evokes characteristic saturation currents in the Jiles-Atherton model, whereas the control-oriented model and the measurement show no saturation. This is caused by the selected saturation function. The Jiles–Atherton model uses a hyperbolic cotangent function with limited adaptability. The proposed saturation function (4) of the control-oriented model is more versatile and better at depicting the correct material characteristic.

The overall shape of the hysteresis curve for CTC2 in Figure 5 is accurately captured by both the Jiles–Atherton model and the control-oriented model. However, the Jiles–Atherton model overestimates the current amplitude at voltage setting u3⁠, while the control-oriented model slightly overestimates it at voltage setting u4⁠.

The hysteresis of transformers from cluster C3 are generally thinner compared to other clusters, as seen in Figure 6. This is also reflected in the identified parameters of both models, where smaller values for the parameters describing the width of the hysteresis are found. The control-oriented model is capable of accurately describing the current response in the linear range of the material as well as in saturation.

Additionally, to the qualitative visual evaluation, the model is compared quantitatively by calculating the normalised error between the measured current imes,o and the simulated current response io over one period of the input voltage. This quality measure is defined in Schwartze et al. ( 2025b) as follows:

(6)

It is calculated at each operation point o corresponding to a voltage setting uo∈1…6⁠. The squared error between the simulation and the measurement is integrated over one period. This error is normalized by the integral over the squared measured current. This ensures that the error metric depends less on the amplitude of the current and thus the values at different operation points are comparable. Table 4 shows the calculated error metric for both models at all operation points uo and CTs. It can be seen that the control-oriented model consistently achieves better values in the error metric – with the exception of u4 on CTC2⁠. The reason for the increased error value is that the model does not estimate the amplitude peak of the current response correctly. A potential solution could be a better fitting saturation function sx1 that can more accurately describe the saturation characteristic of this ferromagnetic core material.

Table 4.

Error metric for all CTs and input voltage settings

CT and operation pointError of control-oriented modelError of Jiles–Atherton model
C0⁠, u10.101350.21627
C0⁠, u20.196950.33138
C2⁠, u30.0346280.60834
C2⁠, u40.819340.012977
C3⁠, u50.00262420.29105
C3⁠, u60.00179070.033927
Source(s): Authors’ own work

In this simulation study, the computational effort of the control-oriented model and the transformer model based on the Jiles–Atherton equations from Albert et al. (2024) is compared. The simulation is performed on an Intel Core i5-8365U CPU @ 1.60 GHz with 16 GB RAM. Both models are implemented in MATLAB/Simulink R2023a (The MathWorks Inc, 2023a). The significant difference between the value of the leakage inductance Lk and the unsaturated main inductance Lm renders solving the parametrized model a stiff problem. This motivates the usage of a stiff solver such as the ODE15s variable step-size solver (The MathWorks Inc,  2023b). For a fair comparison, the same solver settings are selected. The maximum step-size is chosen to be 10-6 s, which means that the solver can automatically adjust the step-size, up to the specified limit.

To compare the computational effort of the simulation of different transformers, the relative simulation time:

(7)

is defined. It is calculated by dividing the computation time Tc that is needed to execute one run of the simulation, with the duration of the simulation Ts⁠, i.e. the time instance that is simulated. The value of Tr explains, how much computation time is needed to predict one second of simulation time. A value of Tr>1 means that the simulation takes longer than reality. Conversely, a value of Tr<1 indicates a simulation that is usable for real-time applications.

The bar chart in Figure 7 compares the mean relative simulation time Tr for the three selected representative CTs. It can be seen that the relative simulation time for the control-oriented model is significantly shorter compared to the Jiles–Atherton model. Each run-time experiment is performed for 10 consecutive times to minimize random effects introduced by the operating system. By solving the pre-compiled model, the time introduced by the compiler is excluded from the computation time measurement.

Figure 7.
A bar chart compares relative simulation times for current transformers C0, C2, and C3 using the Jiles-Atherton model and the proposed model, with error bars representing standard deviation.The image presents a bar chart with relative simulation time on the vertical axis and current transformer types C0, C2, and C3 on the horizontal axis. Two bars are shown for each transformer: one represents the Jiles-Atherton model and the other represents the proposed model. The Jiles-Atherton model shows higher relative simulation times for all three transformer types, reaching nearly 30 for C0, about 24 for C2, and around 22 for C3. The proposed model shows substantially lower times, approximately 10 for C0, 9 for C2, and 6 for C3. Error bars indicate the standard deviation for each value, with the proposed model showing slightly larger variability for C0 and C2 compared to C3.

Comparison of the relative simulation time Tr between the Jiles–Atherton model and the control-oriented model

Note(s): The ODE15s variable step-size solver with a maximum step-size of 10-6 s is selected. For all representative CTs (⁠C0 C2 and C3⁠), the computation time of the control-oriented model is significantly smaller

Source: Authors’ own work

Figure 7.
A bar chart compares relative simulation times for current transformers C0, C2, and C3 using the Jiles-Atherton model and the proposed model, with error bars representing standard deviation.The image presents a bar chart with relative simulation time on the vertical axis and current transformer types C0, C2, and C3 on the horizontal axis. Two bars are shown for each transformer: one represents the Jiles-Atherton model and the other represents the proposed model. The Jiles-Atherton model shows higher relative simulation times for all three transformer types, reaching nearly 30 for C0, about 24 for C2, and around 22 for C3. The proposed model shows substantially lower times, approximately 10 for C0, 9 for C2, and 6 for C3. Error bars indicate the standard deviation for each value, with the proposed model showing slightly larger variability for C0 and C2 compared to C3.

Comparison of the relative simulation time Tr between the Jiles–Atherton model and the control-oriented model

Note(s): The ODE15s variable step-size solver with a maximum step-size of 10-6 s is selected. For all representative CTs (⁠C0 C2 and C3⁠), the computation time of the control-oriented model is significantly smaller

Source: Authors’ own work

Close Figure 7.

The main advantage of the control-oriented model compared to the model based on the Jiles-Atherton equations is the representation of hysteresis. The Jiles–Atherton model is a special implementation of the Duhem hysteresis model (Visintin, 1994), which is given by the differential equation:

(8)

where the input to the model is denoted by u and x describes the output. Based on the time-derivative of the input u˙⁠, the hysteresis is described by the increasing function g1u,x or the decreasing function g2u,x⁠. Therefore, all Duhem-type hysteresis models need the derivative of the input, which is not necessarily available. Thus, for implementation purposes, the derivative is often approximated by a finite difference scheme such that:

(9)

Especially, for noisy inputs or inputs with a large time derivative, such as square waves, this approach can become numerically problematic. In Figure 8, the maximum solver step-size of ODE15s is varied. It is observed that increasing the maximum solver step-size introduces numerical issues which leads to a loss in quality of the predicted output current. To mitigate this numerical issue, the time steps need to be chosen very small, which in turn increases the computational effort and thus the computation time.

Figure 8.
A line graph shows current in amperes versus time in seconds, comparing measurement data with three simulation traces for different step sizes.The image presents a line graph with time in seconds on the horizontal axis, ranging from 0 to 0.18 seconds, and current in amperes on the vertical axis, spanning from negative 15 to positive 15. Four traces are plotted: one for the measurement, and three for simulation runs with step sizes of 1 × 10 to the power of negative 6, 1 × 10 to the power of negative 5, and 1 × 10 to the power of negative 4. The measurement trace follows a waveform pattern with pronounced peaks and troughs, particularly between 0.04 and 0.06 seconds. The simulation traces closely follow the measurement but show minor deviations depending on the step size, allowing for direct comparison of accuracy and stability across the different step settings.

Numerical stability investigation of the Jiles–Atherton model based on the C3 transformer

Note(s): The voltage u6 is applied. With increasing solver step-size the numerical artefacts become more visible. The simulation result with a step-size of 10-4 seconds is not feasible

Source: Authors’ own work

Figure 8.
A line graph shows current in amperes versus time in seconds, comparing measurement data with three simulation traces for different step sizes.The image presents a line graph with time in seconds on the horizontal axis, ranging from 0 to 0.18 seconds, and current in amperes on the vertical axis, spanning from negative 15 to positive 15. Four traces are plotted: one for the measurement, and three for simulation runs with step sizes of 1 × 10 to the power of negative 6, 1 × 10 to the power of negative 5, and 1 × 10 to the power of negative 4. The measurement trace follows a waveform pattern with pronounced peaks and troughs, particularly between 0.04 and 0.06 seconds. The simulation traces closely follow the measurement but show minor deviations depending on the step size, allowing for direct comparison of accuracy and stability across the different step settings.

Numerical stability investigation of the Jiles–Atherton model based on the C3 transformer

Note(s): The voltage u6 is applied. With increasing solver step-size the numerical artefacts become more visible. The simulation result with a step-size of 10-4 seconds is not feasible

Source: Authors’ own work

Close Figure 8.

In Figure 9, the same run-time experiment is performed for the control-oriented model. The maximum step-size of the ODE15s is increased to 10-3⁠, which can be done without any loss of quality in the simulated current response, as demonstrated in Figure 10. The resulting relative simulation time Tr is drastically improved compared to Figure 7. It can be observed that Tr falls below 1⁠, which means that the computation time is shorter that the simulation time. This allows to apply the control-oriented model in real-time applications where the computed current output i is time critical.

Figure 9.
A bar chart compares relative simulation times for current transformers C0, C2, and C3, with error bars representing standard deviation.The image shows a bar chart illustrating relative simulation times for three current transformers labeled C0, C2, and C3. The vertical axis, labeled relative simulation time, ranges from 0 to 30, and the horizontal axis lists the transformer identifiers. Each transformer has two bars, representing results from the Jiles-Atherton model and the proposed model. C0 shows the highest simulation times for both models, followed by C2, while C3 records the lowest. Error bars extend from the tops of all bars, indicating the standard deviation for each measurement. A legend identifies the bar patterns for the two models and denotes the standard deviation indicator.

Due to the structure of the control-oriented model, it is possible to significantly increase the maximum step-size of the ODE15s solver to 10-3⁠, without suffering from numerical issues

Note(s): This further decreases the relative simulation time Tr to values below 1, which showcases the usability of the model for real-time applications

Source: Authors’ own work

Figure 9.
A bar chart compares relative simulation times for current transformers C0, C2, and C3, with error bars representing standard deviation.The image shows a bar chart illustrating relative simulation times for three current transformers labeled C0, C2, and C3. The vertical axis, labeled relative simulation time, ranges from 0 to 30, and the horizontal axis lists the transformer identifiers. Each transformer has two bars, representing results from the Jiles-Atherton model and the proposed model. C0 shows the highest simulation times for both models, followed by C2, while C3 records the lowest. Error bars extend from the tops of all bars, indicating the standard deviation for each measurement. A legend identifies the bar patterns for the two models and denotes the standard deviation indicator.

Due to the structure of the control-oriented model, it is possible to significantly increase the maximum step-size of the ODE15s solver to 10-3⁠, without suffering from numerical issues

Note(s): This further decreases the relative simulation time Tr to values below 1, which showcases the usability of the model for real-time applications

Source: Authors’ own work

Close Figure 9.
Figure 10.
A graph shows the simulated current in amperes over time in seconds, with multiple line styles representing various step sizes.The graph plots current on the vertical axis in amperes and time on the horizontal axis in seconds. A solid black line represents the primary measurement, while additional lines indicate different step sizes: a dotted line for one times ten to the power of negative six, a dashed line for one times ten to the power of negative five, a dashed line for one times ten to the power of negative four, and a dashed line for one times ten to the power of negative three. The curves display the current's variation with time, showing a noticeable peak at approximately zero point zero four seconds, followed by fluctuations. Grid lines are present to support accurate value reading.

Numerical stability investigation of the Jiles–Atherton model based on the C3 transformer

Note(s): The voltage u6 is applied. Increasing the maximum step-size up to 10-3 s introduces virtually no differences between the simulated current responses

Source: Authors’ own work

Figure 10.
A graph shows the simulated current in amperes over time in seconds, with multiple line styles representing various step sizes.The graph plots current on the vertical axis in amperes and time on the horizontal axis in seconds. A solid black line represents the primary measurement, while additional lines indicate different step sizes: a dotted line for one times ten to the power of negative six, a dashed line for one times ten to the power of negative five, a dashed line for one times ten to the power of negative four, and a dashed line for one times ten to the power of negative three. The curves display the current's variation with time, showing a noticeable peak at approximately zero point zero four seconds, followed by fluctuations. Grid lines are present to support accurate value reading.

Numerical stability investigation of the Jiles–Atherton model based on the C3 transformer

Note(s): The voltage u6 is applied. Increasing the maximum step-size up to 10-3 s introduces virtually no differences between the simulated current responses

Source: Authors’ own work

Close Figure 10.

The proposed control-oriented model for inductive transformers can be written in a state-space representation which is especially beneficial for model-based controller design tasks. Additionally, the control-oriented model is valid for arbitrary input signal wave forms. It is validated on a wide range of different transformers, which were chosen based on a comprehensive data analysis. The simulated output current is in good agreement with measurements and other, well-established transformer models. Due to the structure of the proposed control-oriented model the simulation time can be significantly reduced compared to the Jiles–Atherton model, which opens the possibility for real-time applications.

The authors would like to acknowledge the support of the industrial cooperation partner, OMICRON electronics GmbH, for facilitating this study and the experimental analysis.

Albert
,
D.
,
Domenig
,
L.D.
,
Schwartze
,
N.
,
Hadzic
,
Z.
and
Roppert
,
K.
(
2024
), “
Current transformer hysteresis modelling for condition assessment under standard and non-standard operation
”,
CIGRE Paris Session 2024, presented at the CIGRE, Paris
.
Chua
,
L.
and
Stromsmoe
,
K.
(
1970
), “
Lumped-circuit models for nonlinear inductors exhibiting hysteresis loops
”,
IEEE Transactions on Circuit Theory, Presented at the IEEE Transactions on Circuit Theory
, Vol.
17
No.
4
, pp.
564
-
574
, doi: .
Fischer
,
J.
and
Moser
,
H.
(
1956
), “
Die nachbildung von magnetisierungskurven durch einfache algebraische oder transzendente funktionen
”,
Archiv f�r Elektrotechnik
, Vol.
42
No.
5
, pp.
286
-
299
, doi: .
Jäger
,
M.
,
Krüger
,
M.
,
Atlas
,
D.
,
Predl
,
F.
and
Freiburg
,
M.
(
2014
), “
Verfahren und vorrichtung zum testen eines transformators
”,
8 October
.
Jiles
,
D.C.
and
Atherton
,
D.L.
(
1986
), “
Theory of ferromagnetic hysteresis
”,
Journal of Magnetism and Magnetic Materials
, Vol.
61
Nos
1-2
, pp.
48
-
60
, doi: .
Khalil
,
H.
(
2014
),
Nonlinear Control
, (1st ed) .,
Pearson
,
Boston
.
Mörée
,
G.
and
Leijon
,
M.
(
2023
), “
Review of hysteresis models for magnetic materials
”,
Energies
, Vol.
16
No.
9
, p.
3908
, doi: .
OMICRON electronics GmbH
(
2022
), “
CT analyzer user manual
”,
User Manual No. 5.30, OMICRON electronics GmbH, Klaus
, p.
156
.
Pedregosa
,
F.
,
Varoquaux
,
G.
,
Gramfort
,
A.
,
Michel
,
V.
,
Thirion
,
B.
,
Grisel
,
O.
,
Blondel
,
M.
, et al. (
2011
), “
Scikit-learn: Machine learning in Python
”,
Journal of Machine Learning Research
, Vol.
12
No.
85
, pp.
2825
-
2830
.
Preisach
,
F.
(
1935
), “
Über die magnetische nachwirkung
”,
Zeitschrift f�r Physik
, Vol.
94
Nos
5-6
, pp.
277
-
302
, doi: .
Ranta
,
M.
,
Hinkkanen
,
M.
,
Belahcen
,
A.
and
Luomi
,
J.
(
2011
), “
Inclusion of hysteresis and eddy current losses in nonlinear time-domain inductance models
”,
IECON 2011 - 37th Annual Conference of the IEEE Industrial Electronics Society, presented at the IECON 2011 - 37th Annual Conference of the IEEE Industrial Electronics Society
, pp.
1897
-
1902
, doi: .
Rasilo
,
P.
,
Martinez
,
W.
,
Fujisaki
,
K.
,
Kyyrä
,
J.
and
Ruderman
,
A.
(
2019
), “
Simulink model for PWM-supplied laminated magnetic cores including hysteresis, eddy-current, and excess losses
”,
IEEE Transactions on Power Electronics
, Vol.
34
No.
2
, pp.
1683
-
1695
, doi: .
Schubert
,
E.
(
2023
), “
Stop using the elbow criterion for k-means and how to choose the number of clusters instead
”,
ACM SIGKDD Explorations Newsletter
, Vol.
25
No.
1
, pp.
36
-
42
, doi: .
Schwartze
,
N.
,
Albert
,
D.
,
Moschik
,
S.
and
Reichhartinger
,
M.
(
2025
a), “
Validation of the frequency dependency of an electromagnetic control-oriented current transformer model
”,
CIGRE Canada Symposium 2025
, presented at the CIGRE, Montreal.
Schwartze
,
N.
,
Moschik
,
S.
and
Reichhartinger
,
M.
(
2025
b), “
Control-Oriented modeling of single-phase transformers
”,
IEEE Transactions on Magnetics
, Vol.
61
No.
5
, pp.
1
-
12
, doi: .
The MathWorks Inc
(
2023
a), “
Matlab version: 9.14.0. (R2023a)
”,
The MathWorks Inc., Natick, Massachusetts, United States
.
The MathWorks Inc
(
2023
b), “
ODE15s documentation
”,
Ode15s Solve Stiff Differential Equations and DAEs — Variable Order Method, documentation
,
available at:
Link to ODE15s documentationLink to the cited article (
accessed
15 October 2024).
Thorndike
,
R.L.
(
1953
), “
Who belongs in the family?
”,
Psychometrika
, Vol.
18
No.
4
, pp.
267
-
276
, doi: .
Visintin
,
A.
(
1994
),
Differential Models of Hysteresis
,
Springer
,
Berlin, Heidelberg
, Vol.
111
, doi: .
Xu
,
D.
and
Tian
,
Y.
(
2015
), “
A comprehensive survey of clustering algorithms
”,
Annals of Data Science
, Vol.
2
No.
2
, pp.
165
-
193
, doi: .
Yuan
,
C.
and
Yang
,
H.
(
2019
), “
Research on K-Value selection method of K-means clustering algorithm
”,
J, Multidisciplinary Digital Publishing Institute
, Vol.
2
No.
2
, pp.
226
-
235
, doi: .
Zirka
,
S.E.
,
Moroz
,
Y.I.
,
Chiesa
,
N.
,
Harrison
,
R.G.
and
Høidalen
,
H.
(
2015
), “
Implementation of inverse hysteresis model into EMTP—part II: dynamic model
”,
IEEE Transactions on Power Delivery
, Vol.
30
No.
5
, pp.
2233
-
2241
, doi: .
Zirka
,
S.E.
,
Moroz
,
Y.I.
,
Harrison
,
R.G.
and
Chiesa
,
N.
(
2014
), “
Inverse hysteresis models for transient simulation
”,
IEEE Transactions on Power Delivery
, Vol.
29
No.
2
, pp.
552
-
559
, doi: .
Published by Emerald Publishing Limited. This article is published under the Creative Commons Attribution (CC BY 4.0) licence. Anyone may reproduce, distribute, translate and create derivative works of this article (for both commercial and non-commercial purposes), subject to full attribution to the original publication and authors. The full terms of this licence may be seen at Link to the terms of the CC BY 4.0 licenceLink to the terms of the CC BY 4.0 licence.

or Create an Account

Close subscription notice
Close access options