Broad Multi-Parameter Dimensioning of Magnetocaloric Systems Using Statistical Learning Classifiers

Prototyping innovative energy devices is a complex multivariable dimensioning problem. For the case of magnetocaloric systems, one aims to obtain an optimized balance between energy conversion performance, useful power generated, and power consumed. In these devices, modeling is entering a mature phase, but dimensioning is still time consuming. We have developed a technique that dimensions any type of magnetocaloric system by training statistical learning classifiers that are used to simulate the computation of a very large number of systems with different combinations of parameters to be dimensioned. We used this method in the dimensioning of a magnetocaloric heat pump aiming at optimizing the temperature span, heating power, and coefficient of performance, obtaining an f-score of 95%. The respective classifier was used to mimic over 940 thousand computed systems. The gain in computation time was 300 times that of computing numerically the system for each combination of parameters.


INTRODUCTION
One of the goals of the 2030 United Nations agenda for sustainable development is to "ensure access to affordable, reliable, sustainable and modern energy for all." One way for energy consumption to become more sustainable is to increase the energy efficiency of common household devices, such as water heaters. In this context, magnetocaloric systems have been pointed out as a reliable alternative to vapor compression technology for heat pumps and refrigerators (Gschneidner and Pecharsky, 2008;Yu et al., 2010;Kitanovski et al., 2015a). In fact, the coefficient of performance (COP) can attain high values when compared to the vapor compression systems (Yu et al., 2003). Magnetocaloric systems rely on magnetizing and demagnetizing magnetocaloric materials (MCM) for creating temperature gradients instead of compressing and expanding gases (Tishin and Spichkin, 2003;Kitanovski et al., 2015b). Therefore, the magnetocaloric technology does not require the use of nocive gases. One important breakthrough that increased the research on the development of new magnetocaloric systems was the use of the active magnetic regenerative cycle, where the magnetocaloric material acts as both the refrigerant and regenerator (Barclay and Steyert, 1982). This cycle allowed to increase the potential temperature span several times.
Because the regenerator is the core element of magnetocaloric systems, their thermal properties (Franco et al., 2012;Moya et al., 2014;Lyubina, 2017) and geometries are of paramount importance. In that respect, several geometries have been considered Lei et al., 2017). In particular, parallel plates regenerators and packed spheres have been pointed out as the most practical designs (Engelbrecht et al., 2013;Tušek et al., 2013;Aprea et al., 2017;. Although such designs show some level of simplicity, several parameters must be dimensioned. In that respect, several works aiming the dimensioning of prototypes have been performed. Li et al. (2008) optimized an active magnetic regenerative refrigerator with a secondlaw analysis. Nielsen et al. (2010) investigated the influence of the plate and channel thickness, cycle frequency and fluid motion on the performance of a magnetocaloric refrigerator. Roudaut et al. (2011) analyzed the implication of several parameters on the number of units, utilization factor, and axial conductivity. Li et al. (2012) dimensioned an active magnetic refrigerator by using an analytical model. Tagliafico et al. (2012) performed a parametric investigation on the frequency and fluid mass flow rate. The same authors, along with Aprea et al., used dimensioned models for the optimization (Aprea et al., 2013;Tagliafico et al., 2013). Risser et al. (2013) improved a numerical model aiming the optimization of the design of a magnetocaloric refrigerators. Trevizoli et al. (2014) and Trevizoli et al. (2016) included thermomagnetic phenomena and losses in their numerical model. Moreover, the same authors identified the most appropriate parameters and operating conditions following an entropy generation minimization analysis of active magnetic regenerators . Bouchekara et al. (2014) created a multiobjective optimization model. Niknia et al. (2016) included configuration losses in their parametric study. Finally, Roy et al. (2017) used a genetic algorithm in their multiobjective optimization model. Although these parametric investigations achieved the goal of dimensioning magnetocaloric systems, the development of a complete and reliable systematic dimensioning technique has been hindered by the computational cost of performing brute force computation, i.e., of systematically performing the numerical computation of systems with a large number of combinations of parameters to be dimensioned. For instance, for only 7 parameters, and only considering 10 different values for each one, it would be necessary to compute ten million different systems to figure out which values fulfills the requirements.
One way to reduce the number of computed systems for dimensioning procedures is by predicting, through statistical learning classifiers, the outcome of a dimensioned system, if such a method is computationally more efficient assuming trustworthy predictions. Generally speaking, statistical learning classifiers have recently been used in solving several problems in Physics, e.g., in predicting non-linear theories for multiscale modeling of heterogeneous materials (Matouš et al., 2017), in proposing strategies for systems with invariance properties (Ling et al., 2016), or in guiding the creation of improved closure models for computational physics applications (Parish and Duraisamy, 2016). In this work we present a new method for dimensioning magnetocaloric systems. We exemplify the method for a reciprocating heat pump with parallel plates regenerators. We optimize the method by testing with the statistical learning classifiers K-nearest neighbors (KNN), multilayer Perceptron (MLP), support vector machines (SVM), and random forests (RF). The performance of the method is showed for the temperature span, heating power and COP. FIGURE 1 | Geometry for the unidimensional active magnetic refrigerator model. The fluid is static, while the MCM, hot and cold heat exchangers (HHEX and CHEX, respectively) moves cyclically in the horizontal direction. Heat is transferred between the fluid and the 3 remaining elements (h).

MAGNETOCALORIC MODEL
The full process of obtaining the dataset by using the heatrapy python framework is described in this section (hea, 2018). The modeled system is unidimensional and is described in Silva et al. (2018) for hydraulic active magnetic regenerative systems. As depicted in Figure 1, it consists of four elements: one fluid, one MCM, one hot heat exchanger (HHEX), and one cold heat exchanger (CHEX). One should note that this process can easily be substitute by other numerical algorithms, e.g. by using other commercial softwares with similar numerical approaches.

Governing Equations and Boundary Conditions
The governing equations are the following (Petersen et al., 2008a,b;Nielsen et al., 2009Nielsen et al., , 2011Aprea and Maiorino, 2010;Lienhard, 2017): where T is the temperature, x the position, t the time, ρ the density of the material, C H the specific heat at constant H, k the thermal conductivity, V the velocity of the fluid, and the subscript s and f defines the type of material, solid or fluid, respectively. Q sf is the heat transfer power per volume between the solids and fluid. In the present model one considers the fluid framework so that the fluid is static (V = 0) and the other three components are moving forward and backward as shown in Figure 1. In this context, the convective term of Equation (2) vanishes and the overall problem is reduced to a heat conduction problem. The equations were solved by using the Crank-Nicholsen implicit finite difference method with the "implicit_k(x)" solver of the heatrapy package (Silva et al., 2018): where i and n stand for the space and time index, x and t stand for the space and time steps respectively and γ i = The used x and t were 0.001 m and 0.01 s, respectively. Note that the heatrapy package was validated for magnetocaloric fluidic systems (Silva et al., 2019b).
The magnetocaloric effect was considered by using a temperature step change (Silva et al., 2018): where -and + stand for the removal and application of H, so that T + ad and T − ad account for the adiabatic temperature change when applying and removing H, respectively. The T + ad (T) and T − ad (T) curves were obtained with the specific heat C H (T) curves for H=0 and H=1T, which were calculated with the Weiss, Debye and Sommerfeld models for Gd (Petersen et al., 2008b;Silva et al., 2018).
To obtain temperature span values, and since this work exemplifies magnetocaloric systems with the case of heat pumps, insulation was imposed for the HHEX, and the CHEX was kept at a fixed operating temperature. On the other hand, to obtain heating power quantities the temperature of both reservoirs were kept at fixed values. The fixed temperature boundary conditions were modeled by keeping constant the temperature of the boundary point: where b is the boundary position. Thermal insulation is modeled by keeping the temperature of the previous point equal to the boundary point: For more information about the implementation of these conditions consult (Silva et al., 2012(Silva et al., , 2014(Silva et al., , 2016(Silva et al., , 2019a.

Used Materials and Parameters to Be Dimensioned
The used values for C H , T ad , ρ, and k for Gd were extracted from Petersen et al. (2008b). While C H and T ad is strongly dependent on temperature, ρ and k can be considered as being fixed at 7,900 Kg/m 3 and 10.5 W/(mK). Water was chosen for the fluid and copper for both the heat exchangers (HEXs). All the physical properties of the used materials in the modeled system are listed in Table 1.
Seven different parameters were considered in the dimensioning of the magnetocaloric heat pump, and are the inputs of each performed simulation: 1. Operating frequency of the heat pump (ν), with values ranging from 0.01 to 4.51 Hz. 2. Thickness of each parallel plate (ǫ), with values ranging from 0.5 to 2.5 mm. 3. Distance between parallel plates (d p−p ), with values ranging from 0.5 to 2.5 mm. 4. MCM axial length (l MCM ), with values ranging from 1 to 10 cm.  5. HEX axial length (l HEX ), with values ranging from 1 to 10 cm. 6. Distance between the MCM and heat exchangers (d), with values ranging from 1 to 10 cm. 7. Stroke of the water pump, or motion amplitude of the fluid (s), with values ranging from 1 to 10 cm.

Performance Variables
In the present investigation three performance variables were considered: the temperature span of the heat pump, heating power, and COP. The temperature span of a magnetocaloric heat pump T is the temperature difference between the HHEX (T hot ) and CHEX (T cold ) when the CHEX is at a fixed temperature and the HHEX is insulated: These values are calculated when attaining the stationary state criteria, which in the present case occurs when the temperature change between two consecutive cycles is less than 10 −6 . The heating power per sectional area (P heat ) is the energy per time and area that, within one cycle in the stationary state, the system is pumping to the HHEX. To obtain this value one must fix both temperatures splitted below the no load temperature span. In the present case the zero temperature span condition was chosen, i.e., both HEX were fixed at the same temperature (ambient temperature).
One of the most important performance quantities is the COP of a magnetocaloric system, which relates the energy needed, in the form of work, to pump a useful quantity of heat within one cycle. For the case of heat pumps it follows the formula: where P work is the average power per sectional area within a cycle in the stationary state required to pump the heat from the cold to the hot reservoirs. This power can be approximated by where s max and ν max are the maximum stroke and maximum operating frequency considered, i.e., 10 cm and 4.51 Hz, respectively. P max is the total maximum power that need to be given to the motor that moves the magnets and to the water pump. Most of the losses are indirectly incorporated in Equation (9), e.g., viscous losses. In the present case it was considered P max =5 × 10 6 Wm −2 , which was the value to obtain a COP=3 for zero temperature span when using a frequency of 0.5 Hz and stroke of 2 cm. Note that, e.g., a P work of 10 6 Wm −2 means that a regenerator with a cross section of 1m 2 will require 10 6 W of working power. A more realistic regenerator of 2cm 2 , with the same P work , will require 200W of working power.

Dataset
The obtained dataset was extracted from 3,000 simulations, each one considering random values with a uniform distribution of the 7 input variables, known in statistic learning as predictors.
The computation time for the 3,000 simulations was 2 months using 25 cores (Intel(R) Xeon(R) CPU E5-2690 v4 @ 2.60GHz). This corresponds to approximately 4 years of computation if using a single core. The histograms for the performance variables are shown in Figure 2. Note that the count number axis is logarithmic. The large majority of temperature spans, cooling powers and COPs are below 5K, 5 × 10 5 Wm −2 and 5, respectively. The goal of the present work is to dimension a magnetocaloric system without the need to perform intensive computation. For that end a method of predicting the output of computed systems is here exposed by using classifiers. The first step is to classify the 3,000 simulations into acceptable or not acceptable for the defined requirements. This new variable, goal G, will determine the class of each performed simulation. If it complies with the requirements G=1, if not G=0. After classifying the dataset one must balance the data, so that there are equal numbers of simulations with G=0 and G=1. To that end, the exceeding number of the most popular G value is removed randomly. At the end of the first treatment, each predictor must be normalized to 1, and is divided into a training dataset (60% of the data) and a testing dataset (40% of the data).

CLASSIFIERS
A classifier is an algorithm that maps input data to a specific class. Classifiers must be trained with training data for that purpose. In the present case, the trained classifiers must characterize a given system with G=0 or G=1. Seven different classifiers were considered in this work to find the most reliable method for dimensioning the system described above. Each one will be described and tested below for a requirement of a minimum temperature span of 4K. The whole analysis in this work was performed by using the Python sklearn module (skl, 2018).

K-Nearest Neighbors
The K-nearest neighbors (KNN) classifier is known as one of the most simple and straightforward classifier to implement (James et al., 2013). In the present case, the K nearest input parameters of the training dataset is used to determine the G variable. For that end, the average G for the K nearest neighbors is calculated by where y is the vector that represents all the seven input parameters (predictors) of the machine that we want to predict the performance, and z represents the input parameters vector of one computed system of the training dataset. Taking Z as the set of all z vectors, Z K is a subset of Z that includes the K nearest z vectors to y. If Pr(x) > 0.5, than G=1, otherwise G=0. Depending on the K value, different prediction performances may arise, as depicted in Figure 3 for the error of the test dataset, i.e., fraction of bad predictions. For K values below ∼ 13, the model results in under-fitting, i.e., larger values for K improves the precision of the model. For K values above ∼ 13, the model results in overfitting, i.e., larger values for K increases the complexity of the model beyond reality. Therefore K = 13 shows an error of less than 12%: 88% of the predictions are correct.

Multilayer Perceptron
The Perceptron is one of the most famous classifiers (Bishop, 2006). It is well-known for its high accuracy, although its interpretation can be often complex. It consists of applying edge factors to the predictor nodes and inference of the class by using an activation function. By using hidden layers, instead of one for the predictor nodes as seen in Figure 4, the resulting classifier, known as multilayer Perceptron, is one of the most accurate.
The value of each of the hidden nodes and output nodes is determined by where N m,j is the node value in layer m and index j, β m,j is a bias value, N m−1 is the node set of layer m − 1, and α m,j,q are the coefficients for the determination of the node value N m,j . In this work a multilayer Perceptron was created for predicting the goal G. In that respect, a backpropagation algorithm was used and the best parameter for the perceptron FIGURE 3 | Error rate as a function of K neighbors for the KNN classifier.
was obtained by choosing the ones which showed the best performance according to a three-fold cross validation method (James et al., 2013). The parameter that were varied were the regularization term (α), which accounts for the a penalty to avoid the node coefficients α m,j,q to be very large, from 10 −8 to 10 −1 , the learning rate (τ ), which accounts for the weight of each training datum in the classifier, from 0.001 to 1, and the number (h) of nodes of the hidden layer. For this assessment only one hidden layer was used. The fitting procedure showed that if using h = 6, τ = 0.1, and α = 10 −6 the resulting average f-score of 95% was largest. The f-score is the harmonic mean of the precision, which is the fraction of true positives over true positives plus false positives, and the recall, which is the fraction of true positives over the sum of true positives and false negatives. The resulting classifier shows a precision, recall and f-score for G=0 of 96, 94, and 95%, respectively, and 95, 96, and 95% for G=1, respectively.

Support Vector Machines
Support vector machines are classifiers that rely on hypersurfaces that divide different classes (Bishop, 2006;James et al., 2013). This division is dictated by the so-called kernel K. Different kernels represent different types of hypersurfaces, and depends on the predictor vector y and on the set of training data vectors Z. This way, each predictor point y is classified with the following function:  where η 0 is a bias coefficient, and σ i are prefactors. Each σ i must be > 0 and < C, where C is the so-called tuning parameter (James et al., 2013). While η 0 and σ i values are determined during the training procedure, C must be pre-defined. In the present work four different kernels were considered, which will be described in the following. The respective precisions, recalls, and f-scores are in Table 2.

Linear
The linear kernel for the support vector machine is simply the internal product of y and z: K(y, z) = y, z .
Generally speaking, this classifier divides the two classes by a hyperplane. By using a three-fold cross validation, a fitting procedure was performed by varying C from 0.001 to 100. The fitting procedure showed that if using C = 100 the resulting average f-score of 92% was largest.

Polynomial
The polynomial kernel for the support vector machine is where λ is a prefactor, r is a static coefficient and u is the polynomial degree. In this case the hypersurfaces are no longer linear, but follows a polynomial approach. By using a three-fold cross validation, a fitting procedure was performed by varying C from 0.001 to 1,000, by varying λ from 0.0001 to 1, by varying r from 0 to 1, for u =2 and u =3. The fitting procedure showed that if using C = 10, λ = 0.1, u=3, and r = 1 the resulting average f-score of 95% was largest.

RBF
The radial basis function (RBF) kernel for the support vector machine is where φ is a prefactor. This classifier can use closed volumes for the classification. By using a three-fold cross validation, a fitting procedure was performed by varying C from 10 −10 to 10 10 and by varying φ from 10 −15 to 10 5 . The fitting procedure showed that if using C = 10, 000 and φ = 0.01, the resulting average f-score of 94% was largest.

Sigmoid
Finally, the sigmoid kernel for the support vector machine is where tanh is the hyperbolic tangent, ψ is a prefactor and w is a static coefficient. By using a three-fold cross validation, a fitting procedure was performed by varying C from 0.01 to 10,000, by varying ψ from 0.0001 to 1, and by varying b from 0 to 1. The fitting procedure showed that if using C = 1, 000, ψ = 0.01, and b = 0 the resulting f-score of 91% was largest.

Random Forest
A decision tree classifier consists in stratifying, or segmenting, the predictor space into a number of simple regions (James et al., 2013). Although for some very simple problems decision trees are very easy to interpret, it does not have the same level of predictive accuracy when comparing to other classifiers. A Random forest classifier is a tree-based model that consists in building decision trees on bootstrapped training samples (James et al., 2013).
In the present work a fitting procedure with three-fold crossvalidation was performed in order to obtain the best number of estimators, i.e., best number of predictor space segments. In that respect, it was considered a number of estimators ranging from 1 to 100. The best performance was obtained for 71 estimators. The resulting classifier showed a precision, recall and f-score for G=0 of 95, 91, and 93%, respectively, and 92, 96, and 94% for G=1, respectively. These numbers show that for the present case, random forests are comparable to other classifiers.

STATISTICAL MODELING
So far, only the requirement of a minimum temperature span of 4K was considered. In this section, other requirements will be taken into account (temperature span, heating power, and COP), which will allow us to conclude which of the classifiers are more suited for each performance variable and for each range of values. In that respect, a continuous range of performance variables was investigated for all the classifiers described above, by using a fitting procedure that uses a three-fold cross validation. The details of the fitting parameters will be described below.
The variables for the fitting procedure of the KNN classifier were the number of K neighbors. For the linear kernel of the SVM it was varied the C value from 0.001 to 100. For the polynomial kernel, the C value ranged from 0.001 to 1,000, γ from 0.0001 to 1, and the considered polynomial degree was choosen between 2 and 3. For the RBF and sigmoid SVMs, the C value was varied from 10 −10 to 10 10 and from 0.01 to 10,000, respectively. While FIGURE 5 | Precision and recall of the temperature span from 4K to 20K, for the two possibilities of the G variable (0 and 1), and the respective f-score average. In the f-score plot it is also shown the balanced dataset size. φ ranged from 10 −15 to 10 5 , ψ of the sigmoid kernel ranged from 0.0001 to 1, respectively. Moreover, the b value was varied from 0 to 1. Only the number of estimator was considered for the random forests fitting. Finally, for the multilayer Perceptron the α coefficient, the learning rate τ , and number of hidden nodes ranged from 10 −8 to 10 −1 , from 0.001 to 1 and from 1 to 19, respectively.

Temperature Span
A range of temperature spans from 4K to 20K was considered in this work. Figure 5 shows the precision and recall for the two possible G values, as well as the average of the f-score.
One can observe that the larger the temperature span the more fluctuations we have. This is due to the fact that the balanced dataset reduces with the temperature span, i.e., the larger the temperature span the less simulation from the dataset we have. Still, the accuracy of the predictions, i.e., the average f-score, is 90% for almost all classifiers and within the whole range of minimum temperature spans. Only the KNN classifier has an average f-score 90%. In fact, only for a temperature span of 16K and 20K other classifiers show smaller f-scores than the KNN one. From Figure 5 one can also note that, generally speaking, the precision and recall have similar values for all classifiers, so that the large f-scores does not result from only one of the precision FIGURE 6 | Precision and recall of the heating power from 0.25 × 10 6 to 0.95 × 10 6 Wm −2 , for the two possibilities of the G variable (0 and 1), and the respective f-score average. In the f-score plot it is also shown the balanced dataset size.
or recall values. One interesting result is that the random forest classifier behaves in a similar manner as the support vector machine and multilayer Perceptron classifiers. Although they have very similar results, the multilayer Perceptron classifier is the one with the largest average (also over the temperature span) f-score.

Heating Power
The used values for the heating power ranged from 0.25 × 10 6 to 0.95 × 10 6 Wm −2 . The precision and recall metrics can be seen in Figure 6. In this case the fluctuations already occur for the small heating power values, and increases with the heating power. These fluctuations have to do with the fact that the requirements are very strict for the considered range so that the balanced dataset size is very small. In fact, the balanced dataset has less than 250 simulations already for 0.25 × 10 6 Wm −2 , and reduces with the required heating power. By comparing with the temperature span requirement analysis, the balanced dataset decreases from 1,000 to 250 computed systems. Such low number of computed systems decreases the f-score, i.e., the predictor accuracy, from 90% to 80%. Moreover, with such small numbers one cannot conclude which classifier is more accurate FIGURE 7 | Precision and recall of the COP span from 2 to 10, for the two possibilities of the G variable (0 and 1), and the respective f-score average. In the f-score plot it is also shown the balanced dataset size.
in predicting the G variable. Either one needs to reduce the heating power requirement for the same dataset, or increase the number of computed systems in the dataset. This shows clearly the importance of having a robust number of balanced datasets after defining the requirements.

COP
The used COP values ranged from 2 to 10. Similarly to the previous cases, the high intensity of the f-score fluctuations increases with the COP, as observed in Figure 7. Nevertheless, it is clear that the KNN is the worse classifier. This can be seen in Figure 7 in the f-score plot, where for COP=2 the KNN classifier has an f-score of ∼ 85% while the other classifiers have an f-score above 93%. The balanced dataset size is between that of the heating power analysis and temperature span analysis, which explains why only above COP ∼ 7, the f-score starts to show values below 90% for SVM, random forests and MLP classifiers.

SYSTEM DIMENSIONING
In the previous section, the performance of different classifiers were analyzed for individual requirements, i.e., either a minimum temperature span, a minimum heating power or a minimum COP. In this section, one will impose a requirement involving all the three performance variables, through a three-fold filtering procedure, for a multilayer Perceptron and use the resulted classifier for dimensioning the 7 input parameters by computing simulations.

Three-Fold Filtering
To describe the overall process of dimensioning one will use an example with a requirement of 15K for the minimum temperature span, of 0.5 × 10 6 Wm −2 for the minimum heating power, and a 3 for the minimum COP. By using the 3,000 computed systems used in the dataset of this work, and following the procedure described in the dataset subsection, i.e., procedure of balancing and normalizing the dataset, one arrives at 36 simulations with G=0 and 36 with G=1. With this resulting dataset, a multilayer Perceptron was trained with 60% of the data. The testing used 40% of the data, which resulted in an average f-score of 95%. This f-score shows the high performance of the classifier even though only 72 simulations were used.

Dimensioning by Computing Simulations
The method here proposed involves the computation of a large number of simulations with the trained classifier described in the previous subsection. We highlight the difference between the computation of simulations and the computation of systems used for obtaining the dataset. For the example shown above, the computation of simulations was performed for 8 values of ν ranging from 0.625 to 5 Hz, 7 of ǫ ranging from 0.357 to 2.5 mm, 7 l p ranging from 0.357 to 2.5 mm, 7 l MCM ranging from 1.423 to 10 cm, 7 l HEX ranging from 1.423 to 10 cm, 7 d ranging from 1.423 to 10 cm, and 7 s ranging from 1.423 to 10 cm, i.e., a total of 941,192 different combinations. The computation of simulations resulted in 894,287 combinations of input parameters with G=0, and 46,905 combinations with G=1. Since the final goal is to pick up 1 combination, a filtering procedure must be handled. For the present example, the filtering of the 46,905 combinations with ν < 1Hz, l MCM ≤ 8cm, 3cm≤ s ≤ 5cm, d ≤ 2cm, l p ≥1mm, and d ≥1mm, resulted in 4 combinations. In all the possible combinations, the frequency ν was 0.625 Hz, s was 4.286 cm, l MCM was 7.143 cm, ǫ was 1.071 mm, l p was 1.071 mm, and d was 2.857cm. The only difference between the 4 possible combinations was the length of the heat exchangers that needs to be between 1.423 and 5.714 cm, i.e., the length of the heat exchangers does not influence considerably the performance of the system. Figure 8 summarizes the method presented in this work for the dimensioning of magnetocaloric systems with statistical learning classifiers. First, a robust dataset must be computed. Then the requirements, defined by a parameter G (goal) must be applied to the dataset. Then, the dataset must be balanced and normalized. The resulting data is used to train and test a classifier (in the present case a multilayer Perceptron). If the resulting f-score is not acceptable one needs to increase the initial dataset and proceed with the workflow. If the f-score is acceptable, the classifier is used to mimic a very large number of simulations to produce a dataset of computed simulations, i.e., a huge number of mimicked simulations defined by a very large number of combinations of input parameters. Finally, since a significant number of combinations of parameters results in G=1, one needs to filter the input parameters of the dataset of computed simulations.

Computational Cost
The statistical learning method for dimensioning magnetocaloric systems proposed in this work can reduce by orders of magnitude the time consumed in computing numerically the performance of systems. In the example treated in this section it was mimicked the computation of 941,192 different systems by using a dataset of 3,000 simulations. It took approximately 4 years×core of computation for obtaining the 3,000 computed systems. If computing numerically the 941,192 different systems it would take approximately 1,250 years×core of computation, which is unfeasible. In the present study 25 different cores were used at the same time, so that the final computation time was 2 months. Since the computation of simulations lasted only a few hours, i.e., contributes marginally for the overall computation time, the time consumed with the present method is approximately 1 300 of the time consumed if performing the computation of the 941,192 systems with a brute force approach.

CONCLUSIONS
We presented a method for the dimensioning of magnetocaloric systems by using statistical learning classifiers. The method consists in training a classifier and using it to mimic the numerical computation of a very large number of systems with different combinations of variables, e.g., frequency. With the resulting computed simulations, and by filtering some input variables, one is able to obtain a few combinations for the defined requirements. We tested several statistical learning classifiers, and found out that, apart from the K-nearest neighbors classifier, support vector machines, multilayer Perceptrons and random forests present similar behaviors with high f-scores. We applied the method to one example where the heatrapy package was used to produce 3,000 computed systems, and showed that a multilayer Perceptron resulted in an f-score of 95% for a requirement of temperature span above 4K, heating power above 0.5 × 10 6 Wm −2 and COP above 3. By mimicking 941,192 computed systems with the multilayer Perceptron classifier, and by filtering the 7 input parameters, we arrived at 4 different combinations of input variables. We emphasize that this method can be applied to any model (or software) for producing the dataset.

DATA AVAILABILITY STATEMENT
The datasets generated for this study are available on request to the corresponding author.

AUTHOR CONTRIBUTIONS
All authors provided substantial contributions to the conception or the work and interpretation of data. While DS drafted the work, JA and VA revised critically the manuscript. All authors provided approval for publication of the content.