1. Introduction
In electrical power systems, transformers are one of the elements with the highest economic cost, accounting for around 60% of the investment in high-voltage substations. This requires a set of monitoring and diagnostic techniques that affect the life cycle of these important elements. These techniques include the following: dissolved gas analysis, oil quality test, infrared thermography test, power factor test, dielectric dissipation factor test and dielectric oil breakdown test, among others [
1,
2].
In this paper, an analysis of the quality of the oil is made using a combination of electrical, physical and chemical tests [
3]. The transformer oil corresponds to an oil sample without hours of use.
The most important and common tests are related to the dielectric breakdown voltage (BDV), water content, acidity and color. There are three main standards for the determination of the dielectric breakdown voltage of insulating liquids: ASTM D1816–12 (2019), ASTM D877/D877M–19 and IEC 60156:2018. In this paper, UNE-EN 60156, which is based on IEC 60156:2018, was followed.
The results of these tests are used to prevent incipient failures and evaluate preventive maintenance processes, such as transformer oil replacement or recovery [
4].
On the one hand, mineral oils in transformers play an important role as an element of electrical insulation between the parts under voltage and, on the other hand, they help to evacuate the heat generated due to hysteresis losses and eddy currents in iron, as well as the losses due to the Joule effect in the transformer coils. This last condition requires the oil to have a high thermal conductivity and a low coefficient of dynamic viscosity.
The breakdown strength of dielectric oils for transformers will depend on the nature of the impurities present in their solid or gaseous state. Oil analysis is important for extending the life of the transformer.
The state of knowledge of dielectric breakdown voltages in insulating liquids is less developed than in the case of gas and solid dielectrics. The studies carried out may in some cases be contradictory [
5].
Among these studies are those that explain the dielectric breakdown voltages of liquids based on an extension of the dielectric breakdown voltages in gases, in turn based on the avalanche ionization of the atoms caused by electron collisions in the applied field [
6].
Dielectric breakdown voltages in different temperature ranges show little dependence on temperature. This suggests that the cathode emission process is one of field emission rather than thermionic emission [
6].
Electronic theory predicts the relative magnitudes of the dielectric breakdown voltages well, but not the times at which these breakdowns occur in the insulating liquid. This last aspect—the temporal aspect—is partly explained by the presence of polluting particles inside the insulation. These give rise to local breakdowns, which in turn lead to the formation of small bubbles with a much lower dielectric strength and, hence, finally lead to breakdown.
Other phenomena that explain electrical breakdown include the electroconvection of dielectric breakdown, dielectric liquids subjected to high voltages and electrical conduction resulting mainly from charge carriers injected into the liquid from the electrode surface. The resulting space charge gives rise to Coulomb’s force which, under certain conditions, causes hydrodynamic instability, creating an eddy motion of the liquid which yields a convection current.
Thus, the charge transport will be largely via liquid motion and not ionic drift. The key condition for the instability onset is that the local low velocity exceeds the ionic drift velocity [
5].
In most of the studies referenced in the bibliography dealing with analyses of the dielectric breakdown voltage in dielectric oils, approximate analytical equations are used, or numerical methods are employed based on a differential formulation such as the finite elements method (FEM).
In the present study, we propose the finite formulation (FF) method [
7], together with the cell method (CM) [
8,
9] as an associated numerical method to analyze this type of device. Using this methodology, we consider the global magnitudes associated with space-oriented elements such as volumes, surfaces, lines and points of the discretized space, as well as temporal elements, instead of field magnitudes associated with independent variables, i.e., spatial and temporal coordinates [
7].
In addition, equations of the constitutive type—equations of the medium—are clearly differentiated from those of the topological type—equations of balance. In FF, the physical laws that govern the electromagnetic equations are expressed in their integral form. In this way, the final system of equations is obtained directly, without the need to discretize the equivalent differential equations [
8].
In the CM, the topological equations obtained directly from Maxwell’s laws are exact (balance equations), while the constitutive equations obtained from the discretization process are approximate. In the latter case, source-type quantities defined in the elements of the dual mesh must be related to the configuration quantities corresponding to the elements of the primal mesh [
10]. Field magnitudes and the physical properties of the medium are assumed to be constant, at least in the primal mesh. This ensures that the discrete equations are consistent with the continuous constitutive equations, in the sense that the discrete constitutive equations approximate the continuous constitutive equations with an error that decreases with the mesh size [
11].
Most of the research papers in the CM literature focus on the construction of discrete constitutive equations. Among those related to a quasi-electrostatic problem in 2D with plane symmetry, is [
12]. In this paper, an isotropic and anisotropic electrostatic field is studied by means of the CM. In [
13], the electrostatic problem is studied in 2D with plane symmetry. The constitutive equation is used with two approximations. The first approach assumes a uniform field with a triangular base inside each primal cell, and the second approach, which is more general, assumes the uniformity of the fields in subregions of each primal cell, with quadrilateral bases. In [
14], a 2D analysis with axial symmetry (axisymmetric) is performed for a quasi-electrostatic problem regarding a gas-insulated line for an ITER neutral beam injector.
In [
15], an electrostatic induction micro-motor is studied, using the CM in 2D.
In the literature, to the best of the authors’ knowledge, the 3D cell method has not hitherto been used to simulate dielectric breakdown voltage tests on transformer oils as a quasi-electrostatic problem. In this article, we propose the use of a geometric structure for the electrical permittivity constitutive matrix, analogous to the matrix that appears in [
16] for an electromagnetic problem in 3D to calculate eddy currents.
The advantage of using the same geometric structure in the constitutive matrices of conductivity and permittivity in the quasi-electrostatic problem in 3D is that it reduces the complexity of the programmed source code and the execution times. This is because the constitutive matrix is calculated in the assembly of the system of equations, only once for each tetrahedron. The electrical conductivity and electrical permittivity properties of each tetrahedron are multiplied by the common matrix. This is done, element by element, until the final system of equations is complete.
In the present study, a new constitutive matrix is formulated that, using the CM, relates the differences in electrical potentials (magnitudes of configuration) due to primal mesh edges with the electrical flux (magnitudes of source) due to dual mesh planes. The magnitudes of configuration are associated with the edges of a primal mesh made up of tetrahedra, and the source-type magnitudes are associated with the surfaces of a dual mesh (the control volume) obtained in a barycentric division of the primal mesh.
This paper presents an experimental study of the dielectric strength of transformer oil based on the IEC 60,156 standard [
17]. Our contribution consists of characterizing the behavior of the oil an instant before and after the electric arc rupture, combining a low-cost complementary metal oxide semiconductor (CMOS) imaging sensor and a new electrical permittivity matrix
, using the 3D CM [
7,
15]. In the standard test, only the effective value of the dielectric breakdown voltage is obtained. However, the information on the distribution of Kelvin forces [
18] an instant before the dynamic behavior of the arc begins is lost, as is the information on the gases that are produced an instant after the moment of breakdown via the electric arc in the oil.
This last aspect was analyzed by recording images of the movement of the gas bubbles that are produced within the oil. This also allowed the diameter of these bubbles to be measured. The measurement of their magnitudes was used to indirectly obtain the viscosity of the oil. The physical property of viscosity could be obtained by analyzing the post-arc images using an equation to predict the terminal velocity of the rise of isolated bubbles in Newtonian liquids [
19].
The data obtained with the sensors and the results of the simulations complement each other and offer information that would otherwise be lost when strictly following the standard test.
The use of low-cost camera systems in remote-sensing applications is not new. The use and study of low-cost cameras for engineering and scientific applications was addressed in detail in [
20]. In [
21] a study of the corona effect in aeronautical applications was performed, using low-cost Raspberry-Pi-type cameras. In [
22], an array of single-board computers produced by Raspberry Pi, and their associated 8 MP cameras, was used at the University of Cambridge to capture the images required for particle image velocimetry analysis or analysis of the correlation of digital images.
There are two main objectives in this paper. One is to present a new matrix of electrical permittivity, , that predicts the electric field before and after the electric rupture occurs. The other is to measure the kinematic viscosity of the dielectric oil using a low-cost CMOS imaging sensor to measure the distribution of bubbles, their diameters and their rates of ascent after the electric arc occurs. In addition, experiments were performed to estimate the dielectric breakdown voltage in order to obtain the boundary conditions for CM and FEM simulations. In this way, both objectives are related.
This paper has been divided into the following sections.
Section 2 determines the distribution of
and
in the dielectric strength test by applying the 3D CM using the new constitutive matrix
.
Section 3 describes the low-cost 8 MP CMOS imaging sensor used in the experimental studies.
Section 4 presents in detail the numerical results of the CM simulations with
. This matrix is verified by comparing the results obtained with those from the FEM analysis. Finally,
Section 5 presents the experimental setup of the oil testing device and the results that were obtained. In addition, it lays out the test procedure for the dielectric strength and kinematic viscosity of the transformer oil and establishes the theoretical basis of the new procedure for the determination of the dynamic viscosity of the oil. Finally, the data obtained on the kinematic viscosity using the proposed methods are verified by comparing them with the manufacturer’s data.
2. Distribution of and in the Dielectric Strength Test
The formation of the electric arc and the subsequent bubble formation are highly dependent on the estimation of the electric field distribution. Furthermore, its gradient determines the forces per unit volume that act on polluting particles and microbubbles within a dielectric [
23,
24].
It is important to determine the distribution of the field and the gradient of its square as these are factors that determine the DBV.
The following section explains the proposed method for obtaining this field distribution using the 3D CM as an alternative method to the FEM. In this section, a new electric permittivity matrix is proposed for use in the CM.
2.1. New Constitutive Matrix . Discrete Constitutive Equations of Transformer Oil in the Finite Formulation
The electrical constitutive equation in transformer oils is a complex equation based on the Fowler–Nordheim theory [
25]. In most dielectric materials, the conduction current of free carriers is relatively low, since their conductivity is usually several orders of magnitude lower than that of a metal or semiconductor. In new transformer oils at 50 °C it is usually of the order of 1 × 10
−13 S/m, and in used oils of the order of 1 × 10
−11 S/m [
26].
In this paper, a conductive-type model is considered and the conductivity is assumed to be the same throughout the volume of the oil, where the volumetric current density
is directly proportional to the electric field
[
25].
Taking into account the non-zero conductive properties of oil and the fact that it is subjected to an electric field, in the CM the current through the material can be described by its constitutive equation of current flow as a function of potential differences associated with the edges of the primal mesh
ei i = 1:6, as shown in
Figure 1. The constitutive equation shown in Equation (1) was developed in [
16].
The electrical constitutive matrix
is given by the expression
. The matrix
is a function of the electrical conductivity of each tetrahedron; its volume and the dot product of the surface vectors correspond to the dual planes (green planes) of primal edges, edges
ei and
ej i j = 1:6, as shown in
Figure 1.
In dielectrics, however, when the applied electric fields are variable with time, a new contribution to the free current appears, the so-called displacement current. This appears when there is variation in the electric flux with respect to time. The displacement current has two terms:
. The first term within the parentheses only depends on the potential difference in a vacuum. It is independent of the characteristics of the material. The second term depends on the insulating material used. This depends exclusively on the polarization of the dielectric—in this case the transformer oil—and contains the response of the material, which will be different according to the polarization mechanisms that occur for each stimulus of the net applied potential difference (due both to free charges and to those of polarization). We can group both terms with the total empty and material medium permittivities in
, with the constitutive equation given in Equation (2).
Equation (2) relates the differences in electric potential
associated with the edges of the primal mesh (magnitudes of configuration) with the electric fluxes
of the dual planes
, as shown in
Figure 1. In Equation (2), the matrix
is the new electrical constitutive matrix proposed in this work. Given the analogy with Equation (1), the following expression is proposed:
where
is a function of the electric permittivity of each tetrahedron,
, its volume
and the dot product of the surface vectors that correspond to the dual planes of primal edges, edges
ei and
ej. The value of
depends on the type of oil. At 20 °C, this lies between the values of 2.1 and 3.5 [
27]. Therefore, the displacement current is obtained from Equation (4).
The total current will be given by the sum of both contributions, according to Equation (5). The total current
is represented for the dual plane
in
Figure 1.
2.2. Maxwell’s Laws in Finite Formulation Applied to Transformer Oil
Maxwell’s laws, applied in their finite formulation in the dielectric strength test are, first of all, the laws corresponding to the configuration magnitudes. According to [
8], these are as follows.
- (a)
Gauss’s law for the magnetic field, Equation (6)
where
D is the volume–face incidence matrix of the primal mesh, which is equivalent to the standard divergence operator. The magnitude
represents a vector with all the magnetic fluxes associated with the four faces of the primal mesh tetrahedron if
i = 1:4, as shown in
Figure 1.
- (b)
Faraday’s law of induction, Equation (7)
where
C is the face–edge incidence matrix of the primal mesh, which is equivalent to the standard rotational operator.
is a vector of potential differences extended to all edges of the primal mesh and
t is time.
The next two laws correspond to the laws that operate with magnitudes of energy.
- (c)
Generalized Ampere’s law, Equation (8)
where
is the face–edge incidence matrix in the dual mesh, the vector
is a vector of magnetomotive force associated with all the edges of the dual mesh,
is a vector of electric currents extended to all planes of the dual mesh and, finally,
is the electric flux due to the polarization of the dielectric associated with the faces of the dual mesh.
- (d)
Gauss’s law of the electric field, Equation (9)
where
is the incidence matrix of the volumes–faces of the dual mesh and
is the charge contained in each dual volume.
Finally, if the divergence is applied to Equation (8), taking into account Equation (9), Equation (10) is obtained.
2.3. Maxwell’ Laws and Constitutive Equations
In this section, the constitutive equations and Maxwell’s laws given in
Section 2.1 and
Section 2.2 are combined to obtain the final equation in the time and frequency domain. The electrical scalar potential is used for this purpose, and it significantly reduces the number of unknowns.
As we mentioned in
Section 2.1, we assume that the electrical conductivity of the oil is low, of the order of magnitude of
S/m. The permittivity of oil is
F/m. Three time constants can be considered to fit the type of problem to be solved. The first is defined as the charge-relaxation time,
. The second is the electromagnetic time constant,
where
= 10 cm, the characteristic length of the domain. The constant
is the speed of light. The third constant is the magnetic time constant,
. We can consider that if the frequency is 50 Hz, τ ≈ 20 ms and the conditions for the field to be quasi-electrostatic are
and
. For further information on quasi-static laws and time rate expansions, see [
28]. In this way, the field can be considered to be quasi-electrostatic, and Equation (7) can be written as
Since the field is, in this situation, almost electrostatic, it is possible to work with a single electric potential
, where φ is an electric scalar potential. This is imposed on the surface of the electrodes, as shown in
Figure 2. Taking into account the constitutive Equations (2) and (9), this is as follows.
Taking into account Equation (10), and substituting in this equation the free volumetric electric charge
Qf from Equation (12) and the free current
from Equation (1), using
, Equation (13) is obtained.
Equation (13) corresponds to the differential equation derived in [
28]. If the electrodes work at a frequency of
f = 50 Hz, with an angular frequency of
ω = 2π50 s
−1, the final equation in the frequency domain will be Equation (14). It can be observed that this equation involves electrical permittivity and electrical conductivity. This is the equation to be programmed, together with the global electrode current
, which is calculated using Equation (15) where
is an incidence vector of the relative cut between the edges of the oil volume mesh and the surface of one of the electrodes, as shown in
Figure 3. The sum of all currents in that cut is equal to the total current entering or leaving at each of the electrodes.
The matrix representation of both equations is shown in Equation (16).
The unknowns are all the potentials of the primal mesh nodes and the global magnitude of the current associated with the surface of one of the electrodes.
Note that the total current has been defined twice by Equations (5) and (15) because Equation (5) shows the two main components of the total current, i.e., the displacement current and the conductive electric current, while Equation (15) details its explicit composition based on its physical and geometric properties.
In the matrix system (16), the first row is uncoupled from the second. However, we have preferred a compact equation system that includes the total intensity of the electrode. This avoids an additional post-processing calculus. It is also true that we have increased the dimensions of the system by one degree of freedom, corresponding to the total intensity through the electrode. In this way, by solving a single matrix system all the unknowns (degrees of freedom) are obtained at the same time, without post-processing.
Boundary conditions are stabilized on the electrode surfaces. All the nodes on one electrode surface have a value of the electric potential equal to zero and all the nodes on the surface of the other electrode have the dielectric breakdown voltage.
The electrodes have an axisymmetric geometry, but the oil container is cubic, making the global model non-axisymmetric. Hence, we must perform a 3D analysis. Furthermore, the new proposed electrical permittivity matrix is a 3D matrix, and these calculations give generality to the matrix .
2.4. Kelvin Polarization Forces in Dielectric Materials
The force suffered by a supposed spherical particle or a micro-bubble of gas in suspension of radius
r and relative permittivity
, in a liquid with relative permittivity
, in the presence of an electric field
, is calculated using the Kelvin polarization force formula, according to [
23], as shown in Equation (17).
If the
component of
is developed, Equation (18) is obtained.
It is developed in a similar way for
y coordinates as for
z coordinates. As
, Equation (17), according to [
24], is simplified and Equation (19) is obtained.
If
, this force tends to move the particles to the area where the field is stronger, aligning them and forming a bridge that makes it easier for the current to cross the liquid dielectric via that path. The field in the area of the particles increases and the breakdown value is reached [
5].
If the number of particles is not sufficient to bridge the gap, the particles give rise to a local field enhancement, and if the field exceeds the dielectric strength of the liquid, local breakdown will occur near the particles, thus resulting in the formation of gas bubbles which have a much lower dielectric strength and hence, finally, leading to breakdown.
It is important to highlight the dependence of the electric field strength and its gradient and hence the importance of the estimation of the distribution of this electric field and its gradient. This is achieved by applying the method proposed in this paper: the CM together with the proposed new electrical permittivity matrix .
4. Numerical Results of the Simulations in CM with vs. FEM
The system of equations shown in Equation (16), in a sinusoidal steady state, was programmed in C++. As the matrices are sparse and of large dimensions, the Krylov subspace was used. These algorithms were implemented in the PETSc software [
33].
In particular, the linear solver employed was the generalized minimal residual algorithm (GMRES).
The main characteristics of the computer used to make the simulations were: computer model: X399 AORUS PRO; architecture: x86_64; total memory: 128 GB; processors: 24; cpu: 2185.498 MHz; thread(s) per core: 2; core(s) per socket: 12.
The numerical validation of the proposed constitutive equation as shown in Equation (2) and the system of equations shown in Equation (16) was carried out by comparing it with the standard FEM and its implementation within the GetDP software [
34]. A very dense reference mesh was used in 3D with a number of tetrahedra equal to 2,790,589 volumes and 487,435 nodes. The Gmsh software [
35] was used for the mesh and the data visualization. The nodes determine the number of unknowns in the system of equations, as shown in Equation (16).
A section of this mesh and the solution of the potential distribution corresponding to a potential difference between the electrodes of 19,488 V for the dielectric breakdown voltage, can be seen in
Figure 7.
The results corresponding to the spatial distribution of the densities of the conductive and displacement currents are shown in
Figure 8 and
Figure 9, respectively.
The results of the distribution of the force density per unit volume scaled with the coefficient
are shown in
Figure 10 and
Figure 11. These graphs were obtained by cutting the simulation result with the planes 1 X + 0 Y +0 Z − 0.0489 = 0 and 0 X + 1 Y + 0 Z − 0.035 = 0, respectively. Note that the maximum force occurs in a ring around the center.
Figure 12 shows the capture of the electric arc recorded by cameras A, B and C, from the left, center and right, respectively. The images show the breakdown point of the arc in an area that coincides with the one estimated in the model of Equation (19).
Furthermore, the usefulness of calculating Kelvin forces an instant before the dielectric breakdown voltage is reached, using CM, helps in understanding where the dielectric breakdown voltage occurs.
Figure 10 and
Figure 11 show that maximum forces are located in a ring around the electrode, but not at its center. This is coherent with the experimental images shown in
Figure 12, obtained with the imaging sensors.
4.1. Validation of the Numerical Simulations
This section presents the numerical experiments to validate the results obtained by the proposed method, i.e., the cell method CM and the permittivity matrix from Equation (3). This is compared with the FEM in 3D, with a high number of elements (tetrahedra). The latter adapt well to the surface of the electrodes. The geometry of the problem has planes of symmetry. In the end, this always becomes a 3D problem. The experiments solved the entire problem in 3D.
Table 1 summarizes the geometric and physical properties of the dielectric oil used in the numerical experiments. Each type of experiment E1, E2 and E3 was subdivided into two data analyses corresponding to cut A and cut B, which are two characteristic areas of the electrodes, as shown in
Figure 7. Between these two zones, the potential gradients differ greatly and serve to compare the two numerical methods of CM and FEM–getdp using the same mesh.
The number of points in the cut used in all the experiments was 180. The analysis of the metrics was carried out on this basis.
Three numerical experiments were designed: E1, E2 and E3. These consisted of discretizing the electrodes and the volume around them, as shown in
Figure 7. Each experiment had a different mesh density, as indicated in
Table 1. The objective of these experiments was to check the convergence of the CM numerical method, as well as to determine the error produced when comparing it with another numerical reference method, FEM–getdp.
4.1.1. Results of Experiment E1
Figure 13 shows the results of the simulation using the CM and the FEM–getdp. These are the distributions of electric potential obtained in cuts A and B. The lines with a higher slope correspond to cut A and those with a lower slope to cut B.
Figure 14 shows the distribution of the electric field strength modules. A greater electric field is observed for cut A.
In experiment E1, a low-density mesh was used. In
Figure 13, it can be seen that the voltage values in the FEM–getdp and CM are almost coincident for cut A. However, they differ more for cut B. In
Figure 14 there is a notable difference between the values of the electric field given by the FEM–getdp and those given by the CM, in both cut A and cut B, because the mesh is not very dense.
4.1.2. Results of Experiment E2
Figure 15 and
Figure 16 represent the results for the electric potential and the modules of the electric field strength, respectively, corresponding to experiment E2.
In experiment E2, a medium-density mesh was used. In
Figure 15, it can be seen that the voltage values for the FEM–getdp and CM are practically the same. In
Figure 16, there are still differences between the electric field values given by the FEM–getdp and CM, in both cut A and cut B.
4.1.3. Results of Experiment E3
The maximum convergence corresponds to
Figure 17 and
Figure 18, with a total number of tetrahedra of 2,790,589 for experiment E3, as shown in
Table 1.
A high-density mesh was used in experiment E3. In
Figure 17 and
Figure 18, there are no differences in the values given by the FEM–getdp and CM for either the voltage or the electric field.
Therefore, it is observed that as the mesh density is increased, the two methods give coincident solutions. This confirms the validity of the new Mε matrix proposed in this paper for the CM.
4.2. Metrics of Numerical Experiments
To validate the proposed method ( using the CM), a series of comparisons was established between the results obtained with the CM and with the FEM–getdp in the numerical experiments performed.
Metrics are applied to the validation of a model against a reference or pattern. In this case, the model to be validated was the results obtained with the CM and the reference or pattern was the results obtained with the FEM. The CM and FEM are approximate numerical methods, therefore not exact. In both methods a tolerable error is pre-set.
There are various sources of error. One is the truncation of the figures and the accumulation of errors due to the numerical operations performed. Another, which is the one that most affects our problem, is the layout of the cuts. The proximity of the cut to the nodes of the mesh makes the calculated value at the cut more accurate. In contrast, when the cut moves further away from the node, it will be necessary to obtain interpolated values, thus producing a greater error. This happens in both the FEM and CM.
In this study, an analysis was carried out in the FEM with a very dense mesh. These results were used as a reference. Different mesh densities were established in the CM. With the results obtained using the CM, different comparisons were made with the results of the high-mesh-density FEM model. Various statistical indicators (metrics) were used to check the validity of the CM versus the FEM.
The comparisons made are shown in
Table 2.
Of the various statistics that can be used to measure the goodness of fit of a model, the following were chosen for the present study: the coefficient of determination (R
2), the root mean square percentage error (RMSPE), the mean absolute percentage error (MAPE) and the percentage bias (PBIAS).
Table 3 and
Table 4 show the range of these statistics, their optimal values and the value obtained in each the comparisons carried out. When the value obtained is closer to the optimum of the statistic, the analyzed mathematical model has a better goodness of fit.
Figure 19 shows various error histograms for some of the comparisons carried out. Following error theory, an ideally random error distribution should have its mean around zero and a Gaussian or normal distribution.
If the comparisons are grouped following the order {C1, C5, C9}, {C2, C6, C10}, {C3, C7, C11} and {C4, C8, C12}, it is observed, in a generalized way, that increasing the mesh density improves the model for any of the numerical experiments analyzed.
Comparing
Figure 13,
Figure 15 and
Figure 17 with
Figure 14,
Figure 16 and
Figure 18, it can also be seen that the voltages in {C3, C7, C11} and {C4, C8, C12} were modeled, in relative terms, much better than the electric fields in {C1, C5, C9} and {C2, C6, C10}.
The values of the coefficient of determination (R
2) indicate a good fit of the data for all the comparisons except those cases where the mesh was not very dense. The RMSPE, MAPE and PBIAS values are relatively high in the comparison. Nevertheless, all the indicators are in the optimal range. Even the largest error value of 0.36% for the RMPSE, as shown in
Table 3, is a more than acceptable value.
The distribution of the errors, following error theory, conforms to a normal distribution centered on the zero value, as shown in
Figure 19.
The detailed formulations of the metrics used appear in the Annex.
7. Annex
R
2: coefficient of determination.
Value: , the nearer to 1, the better.
Advantage: indicates the proximity to the regression line. The perfect regression line has slope of 1.
Disadvantages: it does not always indicate a linear correlation between the data. If the sample is small, the data may, when enlarged, indicate a nonlinear correlation.
RMSPE: root mean square percentage error.
Value: , the nearer to 0, the better.
Advantage: it is dimensionless and can be used to compare models.
Disadvantages: may underestimate the true measurement as it tries to reproduce the actual data.
MAEP: mean absolute percentage error.
Value:, the nearer to 0, the better.
Advantage: it is a dimensionless and robust measure of error.
Disadvantages: may underestimate the true measurement as it tries to reproduce the actual data.
Value: , the nearer to 0, the better.
Advantage: mean of the deviations, or difference in bias errors, or simply systemic error. It detects aberrant errors or extreme values. If the value is equal to or greater than 1 it is indicative of the existence of these extreme values. Measures the percentage tendency of the simulated data to be larger or smaller than the reference data.
Disadvantages: may underestimate the true measurement as it tries to reproduce the actual data.