Analysis of electrical and thermal models and modeling techniques for polymer electrolyte membrane fuel cells

*Corresponding author: [email protected] Analysis of electrical and thermal models and modeling techniques for polymer electrolyte membrane fuel cells F.J. Asensioa*, J.I. San Martína, I. Zamorab, G. Saldañaa, O. Oñederrab aElectrical Engineering Department, Engineering School of Gipuzkoa (Section of Eibar), University of the Basque Country UPV-EHU, Eibar 20600, Spain bElectrical Engineering Department, Engineering School of Bilbao, University of the Basque Country UPV-EHU, Bilbao 48013, Spain A B S T R A C T Polymeric Electrolyte Membrane Fuel Cell (PEMFC) modeling considering thermal and electrical behavior in a coupled manner is a key aspect when evaluating new designs, materials, physical phenomena or control strategies. Depending on the behavior to be emulated, it is important to choose the modeling technique that best suits the needs required. In this sense, this paper describes the most commonly used PEMFC modeling techniques in the context of analytical-mechanistic approach, semi-empirical approach based on theoretical formulation and empirical correlations, as well as empirical approach based on experimentation with a real system. In addition, an in-depth analysis of PEMFC models at the cell and stack level that emulate the thermal and electrical behavior of these systems in a coupled manner is carried out. A chronological classification of the most relevant models has been made based on the modeling technique used, purpose of the model, state and dimension of the model, and the real system, other developed models or experimental results that have been used to validate the proposed new model. Additionally, guidelines to improve the energy efficiency of PEMFC systems through the development of new models are given. Key words Polymer Electrolyte Membrane Fuel Cell (PEMFC), Modeling techniques, Electrical and thermal behavior, Model Development, Review. 1. Introduction In the current literature, the feasibility of PEMFC (Proton Exchange Membrane Fuel Cell) systems to be used in the field of stationary cogeneration, automotive applications and portable generation has been shown [1], [2], [3]. However, there is still room for improvement in the development of new membranes, catalysts, bipolar plates, collectors, etc., as well as in the management of these systems to be optimally integrated and to increase their energy efficiency and reduce the costs of the system in which they are integrated. On the one hand, to verify the effectiveness of new materials and designs for the elements that form a PEM type cell, it is necessary to develop models that reproduce with high accuracy the behavior of the desired objective (improvement of mass transport, activation losses, ohmic losses, concentration losses, heat transfer, membrane humidification, etc.) when operating with the new design or material. On the other hand, to optimally operate a PEMFC, it is essential to develop accurate models according to the optimization or operation strategy that will be implemented. In this sense, the modeling of PEMFC systems that This is the accepted manuscript of the article that appeared in final form in Renewable and Sustainable Energy Reviews 113 : (2019) // Article ID 109283, which has been published in final form at https://doi.org/10.1016/j.rser.2019.109283. © 2019 Elsevier under CC BY-NCND license (http://creativecommons.org/licenses/by-nc-nd/4.0/) exactly emulate its operation is a key aspect when studying operation strategies, its integration with other systems or the evaluation of control algorithms aimed at maximizing the energy efficiency of the system. To do this, the model must be able to predict accurately the efficiency of the system based on the energy production of the PEMFC. Besides, the choice of the appropriate modeling technique is a fundamental aspect when modeling a PEM cell or a complete stack. In this sense, to verify the effectiveness of new components, materials or phenomena, it usually involves the use of more complex models that require a higher computational cost. However, if the purpose of the model is to evaluate control strategies at system level or to analyze the behavior of the PEMFC when it operates simultaneously with other generation or storage technologies, the modeling technique will be completely different, so that the user can abstract from the processes that are not related to the variables that take part in the process to be controlled or analyzed. At present, there are a large number of models that reproduce with sufficient precision the electrical behavior of PEMFC systems in dynamic regime, as well as in steady state [4]-[5]. On the other hand, there is also a wide variety of models that allow studying the thermal behavior of PEMFC systems [6]-[7]. However, since the electrochemical and thermodynamic reactions that occur in a PEMFC are strongly coupled, to emulate the global behavior of a PEMFC, the model necessarily has to contemplate both aspects in a coupled manner. In this context, there are several research works that present non-isothermal models, which integrate the electrical and thermal behavior together and which can be taken as a reference so that the scientific community could improve what has been done so far in terms of evaluating new designs, materials or control strategies. Taking into account all mentioned above, section 3 of this paper describes the modeling techniques used to develop PEMFC models. Next, section 4 describes the current literature of the PEMFC models that consider coupled electrical and thermal behavior, which are classified in chronological order considering the modeling technique used, the purpose of the model, and the real system or published data used for its validation or comparison, among others. Additionally, section 5 describes the opportunities for improving the energy efficiency of PEMFC systems by developing new models aimed at developing new optimization strategies. Prior to the description of the models mentioned above, the following section 2 introduces general aspects of PEMFC modeling necessary to facilitate the understanding of the descriptions made in the following sections. 2. General aspects of modelling PEMFC systems Regarding the modeling of PEMFC systems, there are many characteristics that allow differentiating some models from others. In general, models of PEMFC systems can be classified according to the formulation used (semi-empirical, empirical or analytical-mechanistic), size of the model (zero, one, two or three dimensions), state of the model (stationary, transient) or the limit of the model (cell, complete stack, etc.), among other characteristics. Table 1 shows the general modeling characteristics of PEMFC systems. Table 1. General modeling characteristics of PEMFC systems Key aspects Options / Characteristics Approach Systemic approach (semi-empirical, empirical) or analytical-mechanistic approach (macroscopic / microscopic) State Stationary; dynamic (transient); real time Limit A component of a cell; a complete cell; stack; stack including auxiliary systems; Integrated PEMFC. Computational domain Single domain; multidomain Spatial dimensions Zero (0D); one (1D); two (2D); three (3D). It is also possible to combine the above. Purpose Thermal analysis, efficiency analysis, structural analysis, parametric studies, element design (gas and cooling channels, GDL, catalyst, collectors), etc. Integration in control systems, integration in cogeneration and trigeneration systems, hybridization with other technologies, etc. Complexity Effects of temperature, two-phase flow, porosity modeling, CO kinetics and poisoning, catalyst degradation, contact resistances, gravitational effects, etc. In the following subsections, the general aspects shown in Table 1 are introduced in order to facilitate the understanding of the description of the models that can be found in the current literature. 2.1. Approach of the model Depending on the requirements and the expected functionality, a systemic model (semi-empirical or empirical) or an analyticalmechanistic model can be used. Models based on an analytical-mechanistic approach are usually very accurate and provide minute details of the operation of the fuel cell at the microscopic or macroscopic level [8]. Its formulation is based on a series of elementary electrochemical and thermodynamic relationships that describe the processes that occur within a PEM (Proton Exchange Membrane) fuel cell. This type of formulation, commonly, addresses three main processes: the electrochemical reactions in the catalytic layers, the migration of protons in the polymeric membrane and the transport of heat and mass in all regions of the cell. The complexity of these models depends on the desired purpose and are generally very difficult to implement due to the highly non-linear nature that describes the behavior of a PEMFC. In this sense, in most cases, due to the highly complex equations on which they are based and the high computational cost required to solve them, these models are not very suitable for use in real-time control systems. Consequently, in most cases they are used as support for the design of components. Models based on a systemic approach focus on the general behavior of the fuel cell, without going into detail in the processes that occur at the atomic level. Depending on whether part or all of the formulation that describes the behavior of the model is based on empirical relationships, a distinction is made between semi-empirical or empirical models, respectively. This type of formulation is widely used when there are parameters that are not known a priori, such as the diffusivity or conductivity of the membrane, which are described by semi-empirical expressions. The empirical formulation covers a very extensive field, since it allows to identify the behavior of the system through a great variety of identification techniques, such as spectroscopy, voltammetry, interruption of the current, test and error methods, as well as implementation methods, such as the electrical circuits or the black-boxes [4], [9]. These types of models are more suitable for use in operation, control or optimization strategies. 2.2. State of the model The state of a model differs between steady state and dynamic (transient) state. In the steady state it is assumed that the variables do not change with time, that is, the volumetric flows, entropies, enthalpies, etc., are not a function of time. Commonly, these types of models are used to characterize the polarization curve of a fuel cell, including losses in the activation, ohmic and concentration region. They are also used to study in detail the behavior of the different parts of a PEM cell, such as the GDL (Gas Diffusion Layer), the CL (Catalyst Layer) or the membrane under different operating conditions. Transient state models include derivatives as a function of time to characterize the behavior of the system dynamically. Generally, time constants are used to characterize the electrochemical behavior of the double layer, the global transfer of heat and mass and the dehydration of the membrane, which usually vary several orders of magnitude [10], [11]. Normally, these types of models are used to analyze the behavior of the system in the event of changes in operating conditions. 2.3. Limit of the model Depending on the application, the model can be a component of a PEM cell (GDL, CL, membrane, etc.), a complete cell including the MEA (Membrane Electrode Assembly), a stack composed of several cells, or a stack including the auxiliary systems (cooling circuit, coolant pump, air compressor, etc.). In general, the models developed using an analytical-mechanistic approach focus on components or on a single cell, since simulating a complete stack would require a very high computational cost. On the contrary, models developed with a systemic approach are more suitable for modeling complete stacks. 2.4. Computational domain A PEMFC model can be based on a single domain or on multiple domains (multi-domain). When a single domain is used, only the source terms (oxygen and hydrogen flow, pressures, velocities and inlet temperatures, etc.) or sink (current density generated, cell voltage, etc.) vary according to the position inside the cell (without internal boundary conditions). In this sense, all the equations are written in the form of a generic convection-diffusion equation, and all the terms that do not conform to that format are coupled with the source or sink term [12]. Multidomain models use different modeling equations in each domain (GDL, CL, membrane, etc.) and require careful management with border boundary conditions, initial and internal (for example, continuity) and external (pressure, temperature, etc.). 2.5. Spatial dimension A model of a cell or a stack can be zero (0D), one (1D), two (2D) and three (3D) dimensions. In this sense, a zero-dimensional model is one in which it is assumed that all the magnitudes (temperature, flows, pressures, etc.) are homogeneous in space, considering only the average values of the inputs and outputs of the system. 1D models contemplate variations of physical magnitudes in only one direction, for example, the variation of the temperature through the GDL, the CL and membrane. Following the same line, 2D and 3D models, contemplate variations of magnitudes in a surface and space, respectively. In the early stages of fuel cell modeling, researchers used 1D models, with varying degrees of complexity, in cross direction assuming the cell arranged as a "sandwich". With this type of models, flows, concentrations, temperatures and electric potentials were analyzed to determine the limit conditions of the cells. These types of models provide enormous information, especially when modeling CLs. 2D models are presented as an improvement of1D models. In this sense, they offer a more realistic view of certain phenomena since spatial variations are considered. These models are usually implemented with the cell in a sandwich arrangement in the xy direction or in a domain along the channel (design known as along-the-channel) in the y-z direction. The sandwich models are mainly used for flow analysis, heat and mass transfer and concentrations that include the effects of bipolar plate and gas channels. The models with domains along the channel are used, mainly, to analyze the concentrations of the reactants along the gas channels. 1D and 2D models can include the same conservation equations as the 3D models, so they provide a lot of information with sufficient accuracy if the boundary conditions and initial conditions are carefully selected. 3D models (x-y-z directions) are the most appropriate when it is required to study the general behavior of a PEMFC. In practice, it can be interpreted that a 3D model is obtained from the combination of 2D domains, which allows studying the blocking effect of bipolar plates, the detailed distribution of current density, temperature or efficiency of a field flow design. As an example, Figure 1 shows a schematic illustration of the different dimensions that are contemplated when modeling a PEM cell, based on the selected coordinate axes (1D, direction y; 2D, direction x-y and y-z; 3D, direction x-y-z). x-y Anode supply (H 2 ) Cathode supply (O 2 ) y-z y Membrane (PEM) Catalyst Mayer (CL) Gas Difussion Layer (GDL) Bipolar Plates (BP) y x z Fig. 1. Schematic illustration of the different dimensions that are contemplated when modeling a PEM cell. In addition to the mentioned dimensions, it is also possible to combine the previous ones (1 + 1, 2 + 1, 2 + 1/2) when a higher dimension model is coupled with a lower dimension model or when one of the magnitudes of the model is calculated in a different dimension from the global dimension of the model. 2.6. Purpose of the model The type of modeling chosen is related to the purpose that will be given to the resulting model. In general, 1D designs are used to analyze the limitation of mass and load transport, while 2D or more dimensions are used for thermal analysis, flow analysis, structural analysis, etc. In addition, if the model is oriented to the design of elements (gas and cooling channels, GDL, catalyst, collectors, etc.), the model is implemented with an analytical-mechanistic approach. When the purpose of the model is to be used to implement a control or optimization strategy, be integrated into a cogeneration system or hybridized with other technologies, zero-dimensional models developed with a systemic approach are usually used (semi-empirical or empirical), although it is also possible to use more complex models depending on the level of details desired. However, for cases in which real-time computation is required, the models are limited to semi-empirical or empirical models of 0D, since the resolution of analytical models of more dimensions supposes a lot of computational time and its implementation would not be viable. for that type of applications [13]. In this context, in [14], to emulate the thermal behavior of a 600 W PEMFC, starting from a 3D computational model, a technique that has allowed converting the 3D model into a 0D model to be implemented in a real time control system has been developed. The fact of developing the 3D model and then converting it to a 0D model provides the combined advantage of both models. In this sense, the 3D model provides high accuracy in the prediction of the temperature profiles of the stack and the refrigeration circuit, and a high speed in computation is obtained by subsequently reducing it to a 0D model. 2.7. Complexity of the model Thanks to the current computational capacity, simulations with PEMFC models are becoming more complex and demanding and include a high level of detail. In this sense, depending on the level of detail required, computational models can contemplate the effects of temperature, efficiency, two-phase flows, porous media, kinetics and carbon monoxide (CO) poisoning, catalyst degradation, contact resistances, gravitational effects, etc. However, depending on the operating conditions or the desired purpose, a series of assumptions and simplifications are usually established that minimize the complexity of the model [15]. 3. PEMFC modeling techniques One way to classify modeling techniques is based on the method by which a model is implemented. The most complex models that incorporate more details are the models based on the technique of analytical-mechanistic formulation. The models obtained through this technique, also called white-box models, are formulated by means of very complex equations that describe the physic-chemical phenomena that intervene in the functioning of the PEMFC. When the purpose of a model is to be used in a control strategy or to be integrated with other systems, simplified models can be considered, so that they emulate only the necessary variables that allow their control or integration. These simplified models, commonly known as gray-box models, are developed based on analytical formulation complemented by a prior knowledge of the system (experimental data), so that some very complex mathematical equations are replaced by empirical equations or mapping tables. Finally, there is another modeling technique in which the relationships between the inputs and outputs of the system are not based on equations or physical laws, but are deduced only through physical experimentation with the real system or throughout experimental databases. Models developed using purely empirical techniques can be classified into two groups: electric models or models based on artificial intelligence. The latter are known by the name of black box models. Fig. 2 shows an outline of the alternatives of modeling approaches, implementation methods, as well as resolution strategies typically used. Fig. 2. Outline of the alternatives of modeling techniques, implementation methods, as well as resolution strategies for PEMFC modeling. In summary, three major groups are distinguished: implementation through analytical-mechanistic formulation, implementation through semi-empirical formulation and implementation through purely empirical methods. Next, these modeling techniques used to model a PEMFC are described. 3.1. Theoretical or analytical-mechanistic modeling techniques Commonly the models developed using this technique are implemented in a multidimensional domain, considering that physical magnitudes vary in two or three axes of the spatial domain. The analytical-mechanistic formulation typically contemplates the laws of Fick, Nernst-Planck and Butler-Volmer to reproduce the phenomena of charge transport (electric and ionic) and mass transfer [16]. The fist law of Fick relates diffusive flow to concentration assuming a steady state. It is postulated that the flow goes from regions of high concentration to regions of low concentration, with a magnitude proportional to the concentration gradient (spatial derivative), or in simplistic terms the concept that a solute will move from a region of high concentration to a low concentration region through a concentration gradient. Equation (1) shows the formulation of Fick's law, assuming a spatial dimension. 𝐽=−𝐷𝑑𝜑 𝑑𝑥 (1) where 𝐽 is the diffusion flow [mol/(m2·s)] of the substance flowing through an area during a time interval, 𝐷 is the diffusion coefficient [m2/s], 𝜑 is the concentration of the substance [mol/m3] and 𝑥 is the position [m]. The law of Nernst-Planck extends the law of diffusion of Fick assuming that the diffusing particles also move with respect to the fluid due to electrostatic forces. Its formulation is based on the conservation of the mass used to describe the movement of a chemical species charged in a fluid medium. Equation (2) shows the formulation of the Nernst-Planck law. 𝜕𝜑 𝜕𝑡=∇·𝐷∇𝜑−𝑢𝜑+𝐷𝑧𝑒 𝑘𝑇𝜑∇ϕ+∂A ∂t (2) where 𝑡 is the time [s], ∇ represents the gradient, 𝑧 is the valence of the ionic species, 𝑒 is the elementary charge [C], 𝑘 is the Boltzmann constant [J/K], 𝑇 is the temperature [K], 𝑢 is the velocity vector of the fluid [m/s], ϕ is the electric potential [V] and A is the vector of the magnetic potential [V·s/m]. The law of Butler-Volmer is one of the most fundamental relations in electrochemical kinetics. It describes how the electrical current in an electrode depends on the potential of the electrode, considering that in the same electrode a cathodic, as well as anodic reaction is produced. Equation (3) shows the formulation of the law of Butler-Volmer. 𝑗=𝑗·󰇧𝑒𝑥𝑝𝛼𝓏𝐹 𝑅𝑇𝐸−𝐸−𝑒𝑥𝑝𝛼𝓏𝐹 𝑅𝑇𝐸−𝐸󰇨 (3) where 𝑗 is the current density at the electrode [A/m2], 𝑗 is the exchange current density [A/m2], 𝐸 is the electrode potential [V], 𝐸 is the equilibrium potential [V], 𝓏 is the number of electrons involved in the reaction, 𝐹 is the Faraday constant [C/mol], 𝑅 is the universal gas constant [J/(mol·K)], and 𝛼 and 𝛼 are the coefficients of anodic and cathodic charge transfer, respectively. The way to solve this type of models based on the formulation described above, in multidimensional domains, is based on the use of advanced numerical methods. The first computer programs based on advanced numerical methods are known as Computational Fluid Dynamics (CFD) programs, which is a computerized tool to simulate the behavior of systems that contemplate fluid flow, heat transfer and other related physical processes. This software works by solving the equations that describe the flow of the fluid over a region of interest, with boundary conditions previously specified in the limit of that region. Currently, there is software that is based on advanced numerical methods and combines the simulation of mechanical, thermal, electrical and fluid-mechanical properties, allowing to completely model virtually any system, including all the parts involved in a PEMFC. Two of the most complete commercial software packages that include multiphysics packages are ANSYS and COMSOL Multiphysics. The methodology of calculation by means of this type of methods consists of the following points:  Pre-processing. During pre-processing, the geometry of the problem to be solved is defined by CAD (Computer-Aided Design) software. Once the geometry is developed, the domains of each element are established. For example, in the case of an MEA, the gas channels, the GDL, the CL, the PEM, etc. X XX X X X Ø Optimal hyperplane Ø ( ) Ø ( ) Ø ( ) Ø ( ) Ø ( ) X Ø ( ) X Ø ( ) X Ø ( ) X Ø ( ) X Ø ( ) X Ø ( ) Fig. 7. Illustration of the mapping of a nonlinear function to a linearly separable function. In the expression (9) the general formulation of a classifier based on SVM is shown. 𝑓(𝑥)=𝑤Φ(𝑥) + 𝑏 (9) where parameters 𝑤 and 𝑏 represent the weight vector and the bias, respectively, that are determined during the training process by minimizing a cost function, and Φ(·) represents the non-linear mapping function to map the input vector 𝑥 in a space of higher dimension, in order to be able to easily separate the data contained in 𝑥 by a linear hyperplane (Fig. 8). Starting from the fact that a training sample (𝑥,𝑦) is a support vector if it satisfies 𝑦𝑓(𝑥)≤1, and denoting the support vectors extracted by 𝑠∈[1,𝐾], the function of the SVM can be represented according to the expression (10). 𝑓(𝑥)=𝐾(𝑥,𝑠)+𝑏   𝐾(𝑥,𝑠)=Φ(𝑥)Φ(𝑠) (10) where 𝐾(·,·) represents the kernel function to be implemented to represent the non-linear mapping of Φ(·). The most commonly used kernel functions are linear, polynomial and RBF (Radial Base Function) [38]. Some examples of PEMFC models based on SVM are those presented in [39], [40]. B) Fuzzy Logic systems Fuzzy logic systems focus on fixed and approximate reasoning as opposed to fixed and exact reasoning. A variable in fuzzy logic can take a range of true values between 0 and 1, instead of taking "true" or "false" values as in traditional binary sets. Since the true value is a range of values of the total set, fuzzy logic systems can only partially handle the truth. A model based on fuzzy logic maps inputs to outputs combining three components: if-then rules, membership functions and logical operators, that is, AND and OR [41]. Numerical data are converted into linguistic variables through membership functions that define how well a variable belongs to the output, that is, an evaluation between 0 and 1 [41]. For the correct definition of these functions or rules, unlike the ANNs, prior knowledge of the user is required, which does not allow to abstract from the physical behavior of the problem to be modeled [9]. Some examples of PEMFC modeling based on fuzzy logic are those presented in [42], [43]. C) Artificial Neural Networks (ANN) Like SVMs and systems based on fuzzy logic, ANNs are systems that can be used to solve classification problems. ANNs are inspired by biological neural networks, and have proved to be a powerful tool for the modeling of non-linear systems [44]. Fig. 6 shows the standard model of an artificial neuron, which was described by D.E. Rumelhart and J.L. McClelland [45], [46]. Xn Xj X2 X1 ∑f () wi1 wi2 wij win -1 Synapse Cellular body yi Output bias Inputs dendrites axon i Neuron Fig. 6. Standard model of an artificial neuron. A biological neuronal system is composed of millions of neurons organized in layers. In the emulation of said biological neural system, a hierarchical structure similar to that existing in the brain can be established by means of an artificial neural system. The essential element will be the artificial neuron, which will be organized in layers. Several layers will constitute a neural network. Finally, a neural network together with the input and output interfaces will constitute the neuronal system (Fig. 7). ∑ F () O U T P U T S Algorithm I N P U T S Neuron Layer Network Neuronal system Fig. 7. Global structure of a neuronal system. There are several models of neural network (linear associator, simple perceptron, Adaline, multilayer perceptron, etc.), architectures (unidirectional, feedback, monolayer, multilayer, etc.) and learning algorithms (Hebbian algorithm, Rosenblatt algorithm, Widrow-Hoff algorithm, backpropagation algorithm, etc.). However, due to the limitations of some of the neural network models and architectures when emulating very complex non-linear systems, within the PEMFC context, the most used model and architecture corresponds to the MLP (Multi-Layer Perceptron) with feedback. The architecture of the MLP has become so popular because, with a single hidden layer, this neural network model can approximate any continuous function in a range up to the desired level [47]. The MLP model is usually trained through the BP (Back Propagation) algorithm. That is why in the current literature, this architecture can commonly be found under the name of backpropagation network. Fig. 8 shows the structure of the MLP with a single hidden layer, as well as the activation function used in the multilayer perceptron, which is a sigmoid function. f (x) x Input Hidden Output j  ij w kj w´ k ´  Target  i x  j y  k z  k t º Fig. 8. Structure (left) and activation function (right) of a MLP. The operation performed by a multilayer perceptron with a single hidden layer and activation functions of the hidden layer and output layer of sigmoid and linear type, respectively, can be defined according to the expression (11). 𝓏=𝑤 󰆒𝑦−𝜃󰆒   =𝑤 󰆒𝑓𝑤𝑥−𝜃   −𝜃󰆒   (11) where 𝑥 are the 𝑛 inputs of the network, 𝑦 are the 𝑜 outputs of the hidden layer and 𝓏 are the 𝑠 outputs of the output layer (those that have to be compared with the targets 𝑡). 𝑤 and 𝜃 represent the weights and biases of the hidden layer, respectively, and 𝑤 󰆒 and 𝜃󰆒 represent the weights and biases of the output layer, respectively. 𝑓 represents a function of sigmoid type. Considering a three-layer MLP (Fig. 4), that is, including a single hidden layer, and with the inputs, outputs, weights and biases of the neurons defined above, given an input pattern 𝑥(𝑟=1,…,𝑁), the global operation of this architecture for each of the 𝑘 output neurons with (𝑘=1,…,𝑠) can be defined by means of (12). 𝓏=𝑤 󰆒𝑦−𝜃󰆒   =𝑤 󰆒𝑓𝑤𝑥−𝜃   −𝜃󰆒   (12) The commonly used cost function is the MSE (Mean Square Error), being for the case of the MLP the expressed in (13) and (14). 𝐸(𝑤,…,𝑤,…,𝑤 󰆒,…,𝑤 󰆒,𝜃,…,𝜃,𝜃󰆒,…,𝜃󰆒) 𝑤ℎ𝑒𝑟𝑒 𝐸∈ℜ(×)(×) (13) 𝐸(𝑤,𝑤󰆒,𝜃,𝜃󰆒)=12(𝑡−𝓏)     (14) There are several methods for the minimization of (14), although the most effective are the Levenberg-Marquardt, Bayesian Regularization and Conjugate Gradient methods [48]. Some examples of PEMFC modeling based on ANNs are those presented in [44], [49]. 4. PEMFC systems modeling: electrical and thermal behavior This section includes a compilation of the models that can be found in the current literature, which emulate both electrical and thermal behavior of PEMFC systems. Models that at least emulate the voltage and temperature variation as a function of current density are considered. In this sense, it should be noted that only non-isothermal models that consider the temperature as an output variable and not as an input variable are considered. Besides, only stack or cell domain models are considered, As seen above, there are several ways to classify a model (approach, regime, domain, etc.). To classify all the models contemplated in this section, a classification based on the dimension of the model has been considered. 4.1. Zero-dimensional models (0D) Among the first research works contemplated in this section that present zero-dimensional models is the model proposed by S. Busquet et al. [50]. The model presented is a novel empirical model for accurately calculating the V-I characteristic of a PEM fuel cell, an electrolyzer or a reversible fuel cell. The empirical adjustment is made through experimentation with a PEMFC of 4 kW. Additionally, the new methodology presented allows to adjust the model through experimentation to any PEMFC or electrolyzer. The model is mathematically coherent and convergent near zero, allowing to interpolate the experimental results of the PEMFC when the model is operated in regions where the current density does not cross the mass transfer limitations. The same year, A. Kazim [51], presents a comprehensive exergy analysis of a 10 kW PEM fuel cell at variable operating temperatures, pressures, cell voltages and air stoichiometry. The analysis is carried out by varying the temperature and pressure, cell voltage and stoichiometry of the air, in order to determine the effect of these variations on the efficiency of the fuel cell. The results obtained highlighted the importance of the operating temperature, pressure, cell voltage and air stoichiometry on the exergy efficiency of the fuel cell. However, they recommend to operate the fuel cell at stoichiometric proportions below 4 to maintain the RH (Relative Humidity) level in the air and to prevent the membrane from drying out at high operating temperatures. One of the first zero-dimensional models implemented by electrical circuits is found in the work of X. Kong et al. [52], in which a fuel cell model that is able to characterize the steady state behavior of the fuel cell, as well as the transient behavior is presented. To make the empirical adjustment of the model, a Ballard Nexa of 1.2 kW of electric power is used. The proposed model shows good behavior (93% accuracy) when validating it with the experimentation results in steady and dynamic state. Later, in [53], an update of the model implemented with electrical circuits in [52] is proposed, in such a way that the dynamics of the system using an ANN based on the MLP is calculated. The ANN is composed of two hidden layers with 30 neurons per layer. By means of this technique it is managed to reduce to half the MSE of the model based on the Nexa Ballard equivalent circuit. In [54], K.C.S. Wang et al. propose a model based on electrical circuits developed in Pspice environment to model the dynamic behavior of a Nexa Ballard of 1.2 kW of electrical power. In the study, it is put special emphasis on the dynamic behavior during the cold start of the stack, as well as on the temperature evolution in time. The simulation results show consistency with the experimental results. In [55], Z. Zhang et al. present a dynamic model of an equivalent electrical circuit for the Nexa Ballard of 1.2 kW, considering the characteristics of the temperature and the equivalent internal resistance. The results of the model show that the developed model can accurately represent the experimental results in a wide range of load conditions. M. Miansari et al. [56], develop a model of a PEM cell to study the effect of different operating conditions, such as temperature, pressure and air stoichiometry on the exergetic efficiencies and the irreversibilities of the cell. The effect of the depth of the anode and cathode channels on efficiency is also calculated. In this sense, the highest efficiencies for a channel depth of 1.5 mm for the anode and 1 mm for the cathode are obtained. ME. Youssef et al. [57], propose a model of zero dimensions, oriented to the study of the effect of temperature, pressure, stoichiometry, thickness of the membrane and thickness of the gas diffusion layer on the behavior of the cell. The results obtained with the model are compared with the results presented by A.R. Maher et al. [58]. K. Hyun-il et al. [59], present a model to characterize the slow transient response of PEMFC. In the work, the estimation of the necessary parameters to obtain the behavior in static, as well as in dynamic regime is presented. The model is validated with experimental results carried out with a Nexa Ballard of 1.2 kW. R. da Fonseca et al. [60], present a stack-level model, oriented to applying a control strategy using the theory of differential flatness. The model responds to the control signals that regulate the most important variables in the air supply subsystem: oxygen stoichiometry and cathode pressure. The model is based on a stack of 5 kW of electrical power composed of 80 cells. In [61], I. San Martín et al. develop a model of the Ballard Nexa 1.2 kW in the MATLAB/Simulink environment, obtaining the parameters that emulate the behavior of the electrochemical and thermodynamic phenomena empirically. The model is validated in static and dynamic regime. The behavior of the model is validated forming a micro-grid with 4 fuel cells of the same type. R. Salim et al. [62], present a technique with a heuristic approach to estimate up to 18 parameters to model a Nexa Ballard of 1.2 kW. The identification algorithm is based on PSO (Particle Swarm Optimization). The resulting model shows good accuracy and requires few mathematical relationships. M.M. Barzegari et al. present in [63] a stack level PEMFC model to investigate the temperature effect on performance of dead-end cascade PEMFC stack with an integrated humidifier and separator. The equations are posed using a semi-empirical approach and are solved using a fourth-order Runge-Kutta method. The model can predict the bulk humidifier and PEMFC temperatures and the stack voltage. Besides, cascade PEMFC operation in a dead-end mode is compared with an open-end mode. Authors propose the obtained model for system identification and control purposes. F.J. Asensio et al. develop in [64] a 0D stack level model of a 600 W PEMFC to evaluate the electrical and thermal efficiency of the system, including power electronics. The model is implemented using ANNs and is developed in MATLAB/Simulink environment. The model provides the temperature and the hydrogen consumption as a function of electrical and thermal demand with good accuracy. Later, same authors present in [65] an improvement of the model by adding the cooling fluid flow rate as an input variable to the model. The temperature of the cooling fluid is calculated using a 3D model in COMSOL Multiphysics environment and a dynamic look-up table is developed to couple the thermodynamic model to the previously developed ANNbased model, resulting in a 0D model more accurate. Authors propose the model to develop real-time control, optimization strategies and to optimally manage the cooling system of the PEMFC. In [66], X. Chen et al. show a thermodynamic model of a PEMFC that includes the main auxiliary components. In this sense, the model contemplates a heat exchanger, a water tank, a cooling pump, as well as the input gas processing components (humidifier and compressor). A parametric study is carried out to study the electrical and thermal efficiency of the PEMFC and the efficiency of the total system. The PEMFC is controlled with MOEA/D (Multi-Objective Evolutionary Algorithm based on Decomposition) to optimize the operating parameters of the system, aimed at maximizing the efficiency and power of the system. The model is obtained by formulating semi-empirical relationships. In [67], A. Kheirandish et al. propose an AI-based model using FCMs (Fuzzy Cognitive Map) to describe the behavior of a fuel cell of 250 W for a power electric bicycle system. Authors use fuzzy rules to explain the cause and effect between concepts. C. Ziogou et al. present in [68] a PEM cell model aimed to apply MPC (Multivariable Predictive Control) strategies. The model is developed using semi-empirical formulation and gPROMS software. Authors show how using the model and strategy implemented fuel cell system operates economically and at a stable environment regardless of the varying operating conditions. J. Chen et al. propose in [69] a dynamic scalable model of PEMFC systems considering two-phase water flow. The model is developed in MATLAB/Simulink environment using the toolbox Simscape. Simulation results show reasonable distributions of current density, temperature, pressure, and two-phase water flow at both the steady state and dynamic operations. Authors highlight that the pattern reconfigurability and the segmentation scalability of the proposed model meet the requirements for both controller design and system analysis for fuel cells. K. Sankar et al. present in [70] a PEMFC system level 0D model aimed at evaluating a nonlinear MSMC (Multivariable Sliding Mode Control) strategy. In addition to the electrochemical and thermodynamic behavior of the PEMFC, the model developed is capable of emulating the behavior of the air compressor, air cooler, primary manifold, supply manifold, humidifier and return manifold. All equations involved are implemented in MATLAB/Simulink environment. 4.2. One-dimensional models (1D) One of the most relevant one-dimensional models can be found in the work presented by J.C. Amplhett et al. [71], in which a PEMFC model designed to predict electrical and thermal behavior, both in steady state and in transitory regime is presented. To do this, the transient characteristics of heat and mass transfer are incorporated into an electrochemical model to form a general model that predicts the transient responses of a PEMFC. The developed model is based on the Ballard Mark V system, which is a PEMFC system of 5 kW of electrical power, formed by 35 fuel cells. Through experimentation with this equipment, a thermal model of the stack based on the conservation of mass and energy balance was developed. The thermal characterization of the stack includes the determination of sensible heat changes in anode, cathode and water circulation flows, the theoretical energy derived from the reaction, the electrical energy produced by the fuel cell, and the heat released through the surface of the stack. The thermal model is coupled to an electrochemical model, relating the power produced by the stack and the temperature of the stack with the amount of heat that must be extracted from the stack. The electrochemical model calculates the electrical power produced by the stack by predicting the cell voltage based on a complex expression involving the operating current, the temperature of the stack, and partial flow rates and pressures of hydrogen and oxygen. A. Rowe et al. [72], present a one-dimensional model of a PEM cell on which the effect of designing and operating conditions on cell efficiency, thermal response and water management is researched. It is shown how the water phase change in the electrodes affects the temperature profile, especially for unsaturated reactants and at low operating temperatures. The simulation results obtained are compared with the experimental results presented by E.A. Ticianelli et al. in [73] and [74]. N. Djilali et al. [75], present a work in which a theoretical model of transport phenomena is formulated for a PEM fuel cell. The model considers the diffusion of humidified fuels and oxidizing gases through the porous electrodes, the transport of water through the electrodes and the membrane, as well as the gradients of heat transfer and gas pressure in the fuel cell. The micro-hydrodynamic phenomena associated with the low permeability of the electrodes are also considered. The model is implemented in a one-dimensional code, and a parametric study is performed, comparing the results obtained with those presented in the work carried out by D.M. Bernardi et al. [76] and E.A. Ticianelli et al. [77]. In this sense, it is verified that, unlike the isothermal and isobaric models, the non-uniform temperature and pressure distributions have a great impact on the flows of liquid water and in the form of simulated vapor in the anode and cathode diffusion layers. In particular, the results indicate that water management requirements (i.e., humidification or removal of water) to prevent possible dehydration of the membrane or flooding of electrodes are much more conservative than when assuming isothermal conditions. It is shown that, in the permeability range of the porous electrodes used in the PEMFC (10-16 – 10-17 m2), the Knudsen diffusion must be considered when modeling the gas transport. In the work presented by X. Xue et al. [78], a dynamic model at PEMFC system level that emulates the temperature, the gas flow through the channels and the capacitance formed by the double charge layer in the MEA is presented. To quantify the dynamic interactions, the PEMFC system is divided into three control volumes: the anode channel, the cathode channel and the fuel cell body; developing the respective dynamic models with grouped parameters. The resulting model is simulated in Simulink environment and validated by comparing the results obtained with those derived from [71]. It is concluded that the model is useful to be used in the optimization and real-time control of PEM fuel cells installed in automotive or stationary applications. Y. Shan et al. [79], propose a model that is constructed based on the layers of a PEM cell, considering the following factors: dynamics of the temperature gradient across the cell, dynamics in the redistribution of water concentration in the membrane, dynamics of the concentration of protons in the catalyst layer, and dynamics in the redistribution of the concentration of reactants in the GDL of the cathode. For the construction of the model, they are based on the parameters presented in [80], [81]. In the work, the results obtained in transitory regime, during the start of the PEMFC and in steady state are shown. C. Wang et al. [82], present the development of a dynamic model for PEM fuel cells using electrical circuits implemented in MATLAB/Simulink and Pspice environments. The model contemplates the effect of double layer charge and the thermodynamic characteristic within the cell. The model responses obtained in steady state and transient conditions are validated by experimental data acquired from the PEMFC Avista Labs SR-12 of 500 W of electric power. The authors propose the model to be used in PEM fuel cell control studies. In [83], S. Kjelstrup et al. present a model for a PEM cell aimed at studying the rate of production of local entropy in various parts of the fuel cell. Authors present five sets of transport equations for a one-dimensional heterogeneous steady-state cell (compatible with the second law of thermodynamics) and solve them by an iterative process. For the implementation of the model, data found in the literature on cells that use the membrane based on Nafion 115 is used. In [84], A.Z. Weber et al. present a one-dimensional model of a PEM cell based on the membrane Nafion 112, which is used to compare the behavior of it in an isothermal and non-isothermal situation. P. Sang-Kyun et al. [85], present a PEMFC model for studying the effect of water variation (contemplating one and two phases) and the temperature distribution along the stack at variable loads on the behavior of the stack. Authors include the cooling circuit in the model. For the construction and validation of the model, an PEMFC of 80 W composed of two cells of 140 cm2 is used. From the developed model, several starting strategies for a PEMFC composed of 20 cells are shown. A.A. Shah et al. propose in [86] a two-phase model that includes a complex kinetic mechanism to describe the electrode reactions. The model is aimed at studying the sulfide poisoning in PEMFCs and is developed in COMSOL Multiphysics environment. Obtained results are compared with experimental data published in [87]. Authors conclude that the kinetic mechanism in the anode is intimately linked with the temperature and the water activity, which yield a wider influence on performance, through the form of the reaction rates. These reactions cause a reduction on the water levels in the anode, which reduce the current density and restricts back diffusion of water via proton migration. P. Hu et al. [88], present an ANN-based model aimed at characterizing the non-linear dynamic behavior of a PEMFC of 1.5 kW of electrical power. For the implementation of the ANN authors use a hybrid algorithm based on PSO and LM (LevenbergMarquardt). The architecture of the ANN is based on MLP with feedback, and consists of 3 neurons in the input layer, 11 in the hidden layer and 3 in the output layer. The model shows good accuracy compared to the real system. S. M. Sharifi et al. [89] present a model aimed at emulating the dynamic response of a PEMFC to variations in the load. The innovation of the model is that it calculates the water content in the membrane and considers the presence of water vapor in the cathode channel. The model is validated with experimental results from several real systems (SR-12, Ballard Mark V and BCS 500), obtaining good correlation with them. In [90], F. Tiss et al. propose a non-isothermal model that takes into account the double layer effect, the geometric capacity and the temperature gradient. The model is developed to operate in a dynamic regime and the results obtained are compared with the isothermal model presented in the work of A. Haddad et al. [91]. Authors conclude that the effect of the temperature distribution significantly influences the cell voltage and the gas flow rate. N. Noguer et al. [92], present a one-dimensional and two-phase model, aimed at developing a method to evaluate the reliability of a PEMFC. The method combines physical modeling with statistical analysis. The model is developed with the ModelicaDymola software and allows analysis in transitory regime with time constants greater than 0.1s. In [93], J.A. Salva et al. develop a one-dimensional, two-phase model of a 50 cm2 PEM cell to emulate the cell voltage and water content in the membrane. For the resolution of non-linear equations, they use the software EES 9.705-3D. The model is validated with a novel technique based on neutron images. Subsequently, in [94], authors use the model to emulate the behavior of a stack composed of 3, 5 and 7 PEM cells. 4.3. Two-dimensional models (2D) The first modeling based on two dimensions, is found in the work presented by T.V. Nguyen et al. [95], in which a twodimensional and single-phase model of heat and mass transfer for a PEM fuel cell is presented. To developed the model, the electro-osmotic coefficient is assumed to be constant. The model is developed as a designing tool for the development of PEMtype cell humidification systems. Subsequently, authors present in [96] an update of the model with updated data on the commercial Nafion membranes of that time, in order to investigate several techniques of membrane humidification. In the same line, in [97] T.F. Fuller et al. develop a model that allows quantifying the amount of heat to be extracted from the cell. J.H. Lee et al. develop in [98] a technique to numerically model a MEA in two dimensions, in order to be integrated as part of a PEMFC dynamic model. The MEA model includes processes, losses and electrical characteristics. The equations used for the construction of the model are based on the previous works developed by J. Kim et al. [99] and J.H. Lee et al. [100]. For the development of the numerical models, authors use a MEA of 350 cm2 and a stack composed of 125 cells. The simulation results show that the model developed using the proposed numerical technique is especially useful to study the effect of temperature, pressure, humidity, and variations in oxygen concentration on the efficiency of the MEA. Based on a semi-empirical formulation technique, V. Gurau et al. present in [101] a model for a two-dimensional PEM cell with two-phase flow, to study the distribution of oxygen and water vapor in the GDL for various current densities. Authors also study the water content in the membrane and several aspects that influence the efficiency of the cell. The equations are solved by means of SIMPLE algorithm for the resolution of CFD systems, developed by S.V. Patankar [21]. The simulation results of the implemented model are compared with the results of the work presented by E.A. Ticianelli et al. [73]. In [102], M. Noponen et al. present a two-dimensional model in which the current density in the active layer of the cathode is modeled assuming an agglomerated material. The model is developed in the environment of COMSOL (FEMLAB 2.3) and validated using a segmented PEM cell. E. Biergersson et al. propose in [103] a two-phase two-dimensional model developed in COMSOL (FEMLAB 2.5), with which the effect of contact resistances between cell components and the effect of different capillary pressures are studied. The catalyst layer is treated as a reactive limit. Authors conclude that heat transfer by convection is negligible under the given operating conditions. S. Litster et al. present in [104] a model for studying heat and mass transfer on the cathode side of a PEM cell. The model, which includes the cooling system, is solved using the CFX software and the SIMPLEC algorithm. The computational domain consists of more than 30,000 mesh elements for the air domain and 580 mesh elements for both electrodes. In [105], J.J. Hwang develops a two-dimensional model to characterize the electrochemical behavior and heat transfer of a PEM cell in a coupled manner. The model contemplates a single phase of water status. The resulting mesh consists of 8,789 elements and the coupled equations are solved by the Newton-Rapshon algorithm. The polarization curve of the model is validated with the work presented in [106]. M. Acosta et al. present in [107] a two-dimensional two-phase model to study the performance of the flow field. In the work emphasis is placed on physical parameters and capillary saturation pressure, and the results are compared with experimental investigations. A saturation level of liquid water of 60% is predicted for low cell voltages. The computational domain consists of the GDL and the CL treated as a thin layer (6,400 elements of equal size). To solve it, authors use the software MUFTE_UG. In [108], H. Wu et al. present a two-dimensional model designed to emulate the dynamics of water transport (in a single phase) in PEM cells under non-isothermal temperature conditions. In the work done, authors focus on studying the effect of membrane width (Nafion 117) on the efficiency of the cell. Authors conclude that the thermal effect has a great impact on the transient behavior of the cell. In [109], Y. Zhang et al. develop a model in two dimensions to study the atmospheric air intake in a PEM cell with an active area of 6 cm2. The effects of the orientation of the cell, the operating conditions and the geometrical parameters are analyzed. Fluent software is used to solve the model together with user-defined subroutines, consisting of 5,748 elements. The results of the model are compared with the experimental results presented in the technical report presented in [110]. In the same research line, B.P.M. Rajani et al. present in [111] a two-dimensional model aimed at studying the cell respiration. Despite being a twodimensional model, the catalyst layer is considered very thin and is treated as a single dimension. Authors also use Fluent and user-defined subroutines to study various effects. In the work authors show that most of the dynamic response processes are within few seconds. [54] 2005 0D Empirical Electrical circuit Stationary/ Transient Characterize the dynamic behavior during the cold start and the evolution of the temperature. Pspice Nexa Ballard 1,2 kW [83] 2005 1D Theoretical Non-linear equations Stationary Determine the local entropy rate in a PEM cell MATLAB 6.0.088 Nafion 115-based membrane cell [127] 2006 3D Semi-empirical FEM Stationary Study the phenomena of transport, saturation of the liquid phase of water, fields of velocities and temperature gradients. Ansys Fluent 6.0.1.2 SIMPLE Algorithm Validation considering balance of species [104] 2006 2D Theoretical FEM Stationary Study heat and mass transfer on the cathode side of a PEM cell CFX SIMPLEC Algorithm Validation considering local Nusselt numbers [128] 2006 3D Semi-empirical FEM Stationary Large-scale study of a 200 cm2 MEA Ansys Fluent 6.0.1.2 SIMPLE Algorithm 200 cm2 MEA [84] 2006 1D Semi-empirical Non-linear equations Stationary Compare the thermal and non-isothermal effects on the behavior of water transport in the cell Nafion 112-based PEM cell [105] 2006 2D Theoretical FEM Stationary Characterize electrochemical behavior and heat transfer in a coupled manner Newton Raphson Algorithm and Gaussian elimination Comparison with [106] [53] 2006 0D Empirical ANN Stationary/ Transient Update the proposed model in [52] implementing the dynamics of the PEMFC by using ANNs DSpace Nexa Ballard 1,2 kW (Comparison with [52]) [55] 2006 0D Empirical Electrical circuit Stationary/ Transient Emulate the dynamic behavior of the Nexa Ballard 1.2 kW by means of the internal equivalent resistance Nexa Ballard 1,2 kW [107] 2006 2D Semi-empirical FEM Stationary Oriented to simulate PEM cells with conventional and interdigitated gas distributors MUFTE_UG Célula PEM con electrodo ELATDS de E-TEK inc. [129] 2007 3D Theoretical FEM Transient Describe the transient process and the dynamic characteristics of a PEM cell with a fluid channel in the form of a coil. Ansys Fluent and SIMPLE Algorithm 21,32 cm2 PEM cell [130] 2007 ~ 3D Semi-empirical Non-linear equations Stationary/ Transient Oriented to apply control strategies in PEM-based stacks MATLAB/ Simulink Nafion 112-based 25 cm2 PEM cell [58] 2007 3D Semi-empirical FEM Stationary Use as a computer-assisted tool to optimize fuel cells with high power density and lower cost Software CFD Comparison with [124] [108] 2007 2D Semi-empirical FEM Stationary/ Transient Study the dynamic characteristics of the cell and the influence of the width of the membrane on the efficiency COMSOL Multiphysics Nafion 117-based PEM cell [109] 2007 2D Semi-empirical FEM Stationary Investigate the behavior of a PEM cell that runs on hydrogen fed at the anode and air supplied by natural convection at the cathode Ansys Fluent Comparison with [110] [111] 2007 2D Semi-empirical FEM Stationary Study the respiration of a PEM cell under atmospheric conditions of 23 °C and 20% of RH. Ansys Fluent 6 cm 2 PEM cell [112] 2007 2D Theoretical FEM Transient Study the effects of liquid water transport and heat transfer phenomena on the transient responses of a PEM cell during a change in cell voltage Ansys Fluent PEM cell [131] 2007 2D 3D Semi-empirical FEM Stationary/ Transient Study of the dynamics of GDL dehydration and its impact on the efficiency of a PEM cell Ansys Fluent 6.0.12 PISO Algorithm Comparison with [104] [85] 2008 1D Semi-empirical Non-linear equations Stationary/ Transient Study the effect of water variation (two phases) and heat to variable loads on the behavior of the stack PEMFC of 80 W composed of two cells of 140 cm2 [86] 2008 1D Semi-empirical FEM Stationary Investigate the effects of hydrogen sulfide contaminant on performance of PEMFCs COMSOL Multiphysics Comparison with [87] [56] 2009 0D Semi-empirical Non-linear equations Stationary Study the effect of stoichiometry, temperature, pressure and air on the exergetic efficiency and irreversibilities of a PEM cell 25 cm2 PEM cell [132] 2009 3D Semi-empirical FVM Transient Study the cold start processes in PEMFCs Ansys Fluent 6.3 UDFs PISO Algorithm Comparison with [133] [88] 2010 1D Empirical ANN Transient Characterize the nonlinear dynamic behavior of a PEMFC MATLAB R2008a PEMFC of 1,5 kW composed of 28 cells of 232 cm2 [89] 2010 1D Semi-empirical Non-linear equations Stationary/ Transient Study the dynamic response of the PEMFC to variations in the load MATLAB/ Simulink SR-12 Ballard Mark V BCS-500 [57] 2010 0D Semi-empirical Non-linear equations Stationary Study of the effect of temperature, pressure, stoichiometry, thickness of the membrane and thickness of the gas diffusion layer on the behavior of a PEM cell MATLAB Comparison with [58] [59] 2010 0D Semi-empirical Non-linear equations Stationary/ Transient Oriented to characterize the slow transient response of a PEMFC. Excel Runge-Kutta of order 4 Nexa Ballard 1,2 kW [135] 2010 3D Semi-empirical FVM Stationary/ Transient Study the effect of non-equilibrium phase transfer Ansys Fluent 6.3.26 Experimental results presented in [136]. [90] 2013 1D Semi-empirical Electrical circuit Transient Study of the behavior of a PEM cell in non-isothermal and transitory conditions Mathcad Comparison with [91] [60] 2014 0D Empirical Electrical circuit Transient Oriented to apply control techniques on the air subsystem (oxygen stoichiometry and cathode pressure) MATLAB/ Simulink PEMFC of 5 kW composed of 80 cells [61] 2014 0D Empirical Electrical circuit Stationary/ Transient Emulate the stationary and transient behavior of a microgrid formed by 4 PEMFC MATLAB/ Simulink Nexa Ballard 1,2 kW [62] 2015 0D Empirical Non-linear equations Stationary/ Transient Experimental methodology Identify up to 18 parameters of the model experimentally MATLAB fourth-order RungeKutta method Nexa Ballard 1,2 kW [92] 2015 1D Semi-empirical Non-linear equations Stationary/ Transient Oriented to develop a method to evaluate the reliability of a PEMFC Modelica-Dymola 220 cm 2 PEM cell manufactured by CEA LITEN [113] 2015 2D Semi-empirical FVM Transient Study the effect of permeability on the dynamic field in the PEM fuel cell FVM-based software Generic PEM cell [63] 2016 0D Semi-empirical Non-linear equations Transient Investigate the temperature effect on performance of deadend cascade PEMFC stack with integrated humidifier and separator. Model suitable for system identification and control purposes. Fourth-order RungeKutta method PEMFC stack composed of 4 cells of 225 cm2 [93] 2016 1D Semi-empirical Non-linear equations Stationary Emulate the cell voltage and the water content in the membrane (new validation technique based on neutron images) EES 9.705-3D 50 cm2 PEM cell [94] 2016 1D Semi-empirical Non-linear equations Stationary Emulate static behavior at stack level EES 9.705-3D 50 cm 2 PEM cell [137] 2017 3D Semi-empirical FVM Stationary Investigate the effect of a blockage in the flowfield channel of a PEM fuel cell on mass transfer of reactant gas from the channel into the catalyst layer Ansys Fluent SIMPLE Algorithm Previously published results in [138] [64] 2017 0D Empirical ANN Stationary/ Transient Emulate the hydrogen consumption and temperature as a function of the electrical and thermal load MATLAB/Simulink PEMFC of 600 W [139] 2017 3D Semi-empirical FVM Transient Investigate the water removal processes in a PEMFC during the gas purging prior to its shutdown Ansys Fluent and UDS equations Validated using experimental data [140] and numerical data [141] [66] 2017 0D Semi-empirical Non-linear equations Stationary Apply a multi-objective optimization algorithm (MOEA/D) aimed at maximizing the energy efficiency of the system PEMFC of 5 kW composed of 75 cells [142] 2017 3D Theoretical FVM Stationary Study the cooling flow fields effect on PEMFC performance Ansys Fluent PEM cell presented in [144] [67] 2017 0D Empirical Fuzzy logic Transient Determine de behavior of a fuel cell electric bicycle system Fuzzy Cognitive Map Rule-based FCM 22 cell air-cooled stack of 250 W [143] 2017 3D Theoretical FVM Stationary Study the cooling flow fields effect on PEMFC performance Ansys Fluent PEM cell presented in [144] [145] 2017 3D Semi-empirical FVM Stationary Study reactants distribution, current density and final power in PEMFCs for four square tubular configurations (simple, DPIE, DBIE and TPIE) Ansys Fluent 14 SIMPLE algorithm Validated using data reported in [124] [146] 2017 3D Semi-empirical FVM Stationary/ Transient Evaluate the performance of a PEM fuel cell stack with variable inlet flows under simulated driving cycle conditions Ansys Fluent 15 UDFs Scalars in C code PEMFC of 320 cells and 1,600 cm2 of active catalyst area [65] 2018 0D Theoretical FEM Empirical ANN Stationary/ Transient Provide a tool to develop new optimization strategies in real time applications COMSOL Multiphysics MATLAB/Simulink PEMFC of 600 W [68] 2018 0D Semi-empirical Non-linear equations Transient Apply a multivariable model predictive control (MPC) strategies for PEM fuel cells. gPROMS software 6 W PEM cell [147] 2018 3D Semi-empirical FVM Stationary Investigate the performance of a PEMFC. Ansys Fluent 16.2 Validated using data reported in [124] [148] 2018 3D Semi-empirical FVM Stationary Study the effects of agglomerate model parameters on transport characterization and performance of PEM fuel cells Ansys Fluent and UDS equations Validated using data reported in [150] [114] 2018 2D Semi-empirical FEM Stationary Investigate the competition between current collection and oxygen supply. Used for rib/channel design. COMSOL Multiphysics Comparison with results from AutoStack-CORE project [115] [69] 2018 0D Semi-empirical Non-linear equations Stationary/ Transient Develop a dynamic scalable segmented model of PEMFC systems with two-phase water flow MATLAB/Simulink Simscape [151] 2018 3D Semi-empirical FVM Stationary Investigate the steady polarization curves and long-term stability of poisoned PEMFC, as well as evaluate the performance of various commercial GDLs. Ansys Fluent and UDF equations Validated using data reported in [152] and [153] and a 25 cm2 PEM cell [149] 2018 3D Semi-empirical FVM Stationary Determine the effects of GDL deformation and obtain the transport characteristics in PEMFCs with interdigitated flow fields Ansys Fluent and UDS equations Validated using data reported in [150] [70] 2018 0D Semi-empirical Non-linear equations Stationary/ Transient Develop a Nonlinear multivariable sliding mode control for a reversible PEM fuel cell integrated system MATLAB/Simulink PEMFC composed of 35 cells [116] 2018 2D Semi-empirical FEM Stationary Investigate the reactant gas crossover effect in a PEM cell FreeMem++ PEM cell 5. Opportunities to improve energy efficiency through new model developments An important aspect when modeling a PEMFC system is that the developed model must contemplate the effect that temperature has on the efficiency of the system. In this context, it would be interesting if models would consider the electrical and thermal behavior of the PEMFC emulate as precisely as possible the effect that temperature has on the speed at which the electrochemical reactions occur. Since the operating temperature is directly related to the amount of heat extracted through the cooling circuit, a key aspect to keep into account is that the model has to consider the thermal management of the system. In this sense, among all research works that can be found in the current literature, only a small group of them focuses on models that include the thermal management of the system through the refrigeration circuit. However, most models that include in detail the effect of heat extraction through the cooling system are developed by a complex formulation, which does not allow to use the model to evaluate real-time control strategies or optimization strategies, for which a very reduced calculation time is required. Considering all mentioned above, the need to develop zero-dimensional models of PEMFC that contemplate the regulation of the flow rate of the cooling fluid as an input variable to the model is detected. This will allow to establish the optimal reference of electrical and thermal production of the fuel cell in order to maximize the electrical efficiency of the system and to minimize the production costs. This aspect is especially important if it is considered that a PEMFC could operate connected to the power grid and that other auxiliary electrical and / or thermal generation devices could be integrated with the PEMFC in a microgrid. In this way, the model developed could be used to predict the behavior of the system in different situations (for predicted electrical and thermal load profiles), which in turn will allow establishing the optimum operating reference that maximizes efficiency and minimizes operating costs. of the system based on the hourly prices of electricity, prices of the fuels used by other generators, as well as all the operation and maintenance costs of all the devices 6. Conclusions In this paper the general modeling aspects of the PEMFCs have been introduced. Likewise, a search was made of the state of the art of the PEMFC models that emulate in a coupled manner the thermal and electrical behavior, at cell and stack levels, contemplating a total of 78 research works. It has been possible to verify how approximately 80% of the models developed in the last 5 years (2014-2018) correspond to zero-dimensional and three-dimensional models, practically in the same proportion. Due to the existing computational capacity, it has been proven how simulations with PEMFC models are becoming more complex and exigent and include a high level of detail. In this sense, the global modeling of PEMFC has evolved from steady to dynamic, from one-dimensional to complex three-dimensional models, from isothermal to non-isothermal, from single-phase to multiphase and recently from straight channels to field structures with more complex flow, such as serpentine or interdigitated flow fields. It has been detected that, to use a PEMFC model in real-time control or optimization strategies, it is usually necessary to resort to zero-dimensional system models to avoid the high computational cost involved in the inclusion of details and phenomena that occur at a microscopic level inside a PEM cell. In this sense, it has been detected that MATLAB is the most used engineering tool by the research community in terms of modeling PEMFCs in zero dimensions. In terms of multidimensional modeling, Ansys Fluent (based on FVM) has been shown as the most widely used multiphysics software, followed by COMSOL Multiphysics (based on FEM). In the context of zero-dimensional systemic modeling, the modeling approach based on empirical techniques has been more widely used, due to the good relationship between simplicity and precision that it provides. Among these techniques, the identification of parameters through EIS and their subsequent implementation using equivalent electrical circuits has been shown as one of the most suitable alternatives to emulate the dynamic characteristics of the electrical behavior of the PEMFCs. However, to contemplate the coupled behavior, electrical and thermal, it has been found that artificial intelligence based on ANNs show a great potential in the modeling of these electrochemical devices of a highly non-linear character. Finally, it has been verified that, among all the models that have been identified in the scope of the application of operation and optimization strategies, the works that have considered the flow of the refrigerant as a control variable from the point of view of energy efficiency are very limited. Since the operating temperature of the PEMFC is directly related to the efficiency of the stack, and which in turn, the operating temperature is related to the amount of heat extracted from the PEMFC, it is detected the need to implement optimization strategies to control PEMFC systems that will take the regulation of the cooling flow as a control variable to maximize the electrical and thermal efficiency of the system. In addition to the efficiency of the system, the optimization strategy could be implemented as a multi-objective strategy that would contemplate the minimization of operating and maintenance costs and harmful emissions. Acknowledgement The authors thank the support from the Spanish Ministry of Economy, Industry and Competitiveness (project ENE2016-79145R AEI/FEDER, UE), the Basque Government (project ELKARTEK KK-2017/00083 and GISEL research group IT1083-16), as well as from the University of the Basque Country UPV/EHU (project EHUA15/25 and research group funding PPG17/23). Reference List [1] Isa NM, Tan CW, Yatim AHM. A comprehensive review of cogeneration system in a microgrid: A perspective from architecture and operating system. Renewable and Sustainable Energy Reviews 2018;81:2236-2263. [2] Lü X, Qu Y, Wang Y, Qin C, Liu G. A comprehensive review on hybrid power system for PEMFC-HEV: Issues and strategies. Energy Conversion and Management 2018;171:1273-1291. [3] Shimpalee S, Ohashi M, Van Zee JW, Ziegler C, Stoeckmann C, Sadeler C, et al. Experimental and numerical studies of portable PEMFC stack. Electrochimica Acta 2009;54:2899-2911. [4] Saadi A, Becherif M, Hissel D, Ramadan HS. Dynamic modeling and experimental analysis of PEMFCs: A comparative study. Int J Hydrogen Energy 2017;42:1544-1557. [5] Alberro M, Marzo FF, Manso AP, Domínguez V, Barranco J, Garikano X. Electronic modeling of a PEMFC with logarithmic amplifiers. Int J Hydrogen Energy 2015;40:3708-3718. [6] Oh K, Chippar P, Ju H. Numerical study of thermal stresses in high-temperature proton exchange membrane fuel cell (HTPEMFC). Int J Hydrogen Energy 2014;39:2785-2794. [7] Özden E, Tolj I, Barbir F. Designing heat exchanger with spatially variable surface area for passive cooling of PEM fuel cell. Appl Therm Eng 2013;51:1339-1344. [8] Lile JRD, Zhou S. Theoretical modeling of the PEMFC catalyst layer: A review of atomistic methods. Electrochim Acta 2015;177:4-20. [9] Petrone R, Zheng Z, Hissel D, Péra MC, Pianese C, Sorrentino M, et al. A review on model-based diagnosis methodologies for PEMFCs. Int J Hydrogen Energy 2013;38:7077-7091. [10] Weber AZ, Newman J. Modeling Transport in Polymer-Electrolyte Fuel Cells. Chemical Reviews 2004;104:4679-4726. [11] Wang C. Fundamental Models for Fuel Cell Engineering. Chem Rev 2004;104:4727-4766. [12] Cheddie D, Munroe N. Review and comparison of approaches to proton exchange membrane fuel cell modeling. J Power Sources 2005;147:72-84. [13] Cheddie DF, Munroe NDH. Semi-analytical proton exchange membrane fuel cell modeling. J Power Sources 2008;183:164-173. [14] Asensio FJ, San Martín JI, Zamora I, Oñederra O. Model for optimal management of the cooling system of a fuel cellbased combined heat and power system for developing optimization control strategies. Applied Energy 2018;211:413-430. [15] Siegel C. Review of computational heat and mass transfer modeling in polymer-electrolyte-membrane (PEM) fuel cells. Energy 2008;33:1331-1352. [16] Wu H. A review of recent development: Transport and performance modeling of PEM fuel cells. Appl Energy 2016;165:81-106. [17] Lim YI, Le Lann JM, Joulia X. Accuracy, temporal performance and stability comparisons of discretization methods for the numerical solution of Partial Differential Equations (PDEs) in the presence of steep moving fronts. Comput Chem Eng 2001;25:1483-1492. [18] Patankar SV. Numerical Heat Transfer and Fluid Flow CRC Press; 1980. [19] Botte GG, Ritter JA, White RE. Comparison of finite difference and control volume methods for solving differential equations. Comput Chem Eng 2000;24:2633-2654. [20] Vinokur M. An analysis of finite-difference and finite-volume formulations of conservation laws. J Comput Phys 1989;81:1-52. [21] S.V. P. Heat Trannsfer and Fluid Flow Hemisphere: Nueva York; 1980. [22] Logan DL. A First Course in the Finite Element Method. Sixth edition ed. Boston, MA : Cengage Learning; 2017. [23] Surana KS, Allu S, Tenpas PW, Reddy JN. k-version of finite element method in gas dynamics: higher-order global differentiability numerical solutions. Int J Numer Methods Eng 2007;69:1109-1157. [24] Springer TE, Zawodzinski TA, Gottesfeld S. Polymer Electrolyte Fuel Cell Model. Journal of The Electrochemical Society 1991;138:2334-2342. [25] Amphlett JC, Baumert RM, Mann RF, Peppley BA, Roberge PR, Harris TJ. Performance Modeling of the Ballard Mark IV Solid Polymer Electrolyte Fuel Cell: II . Empirical Model Development. Journal of The Electrochemical Society 1995;142:9-15. [26] Pisani L, Murgia G, Valentini M, D’Aguanno B. A new semi-empirical approach to performance curves of polymer electrolyte fuel cells. J Power Sources 2002;108:192-203. [27] Maggio G, Recupero V, Pino L. Modeling polymer electrolyte fuel cells: an innovative approach. J Power Sources 2001;101:275-286. [28] S. C. Page, A. H. Anbuky, S. P. Krumdieck, J. Brouwer. Test Method and Equivalent Circuit Modeling of a PEM Fuel Cell in a Passive State. IEEE Transactions on Energy Conversion 2007;22:764-773. [29] Kim J, Lee S, Srinivasan S, Chamberlin CE. Modeling of Proton Exchange Membrane Fuel Cell Performance with an Empirical Equation. Journal of The Electrochemical Society 1995;142:2670-2674. [30] Chu D, Jiang R, Walker C. Analysis of PEM fuel cell stacks using an empirical current–voltage equation. Journal of Applied Electrochemistry 2000;30:365–370. [31] Han I, Shin HK. Modeling of a PEM fuel cell stack using partial least squares and artificial neural networks. Korean Chem Eng Res 2015;53:236-242. [32] Meiler M, Schmid O, Schudy M, Hofer EP. Dynamic fuel cell stack model for real-time simulation based on system identification. J Power Sources 2008;176:523-528. [33] Springer TE, Zawodzinski TA, Wilson MS, Gottesfeld S. Characterization of Polymer Electrolyte Fuel Cells Using AC Impedance Spectroscopy. Journal of The Electrochemical Society 1996;143:587-599. [34] Zhiani M, Majidi S. Effect of MEA conditioning on PEMFC performance and EIS response under steady state condition. Int J Hydrogen Energy 2013;38:9819-9825. [35] Büchi FN, Marek A, Scherer GG. In Situ Membrane Resistance Measurements in Polymer Electrolyte Fuel Cells by Fast Auxiliary Current Pulses. Journal of The Electrochemical Society 1995;142:1895-1901. [36] J. R. J. Larminie. Current interrupt techniques for circuit modelling. IEE Colloquium on Electrochemical Measurement 1994:12/1-12/6. [37] Shen X, Niu L, Qi Z, Tian Y. Support vector machine classifier with truncated pinball loss. Pattern Recognit 2017;68:199210. [38] Ren J. ANN vs. SVM: Which one performs better in classification of MCCs in mammogram imaging. Knowledge-Based Syst 2012;26:144-153. [39] Zhong Z, Zhu X, Cao G. Modeling a PEMFC by a support vector machine. J Power Sources 2006;160:293-298. [40] Han I, Chung C. A hybrid model combining a support vector machine with an empirical equation for predicting polarization curves of PEM fuel cells. Int J Hydrogen Energy 2017;42:7023-7028. [41] Dash S, Rengaswamy R, Venkatasubramanian V. Fuzzy-logic based trend classification for fault diagnosis of chemical processes. Comput Chem Eng 2003;27:347-362. [42] Kishor N, Mohanty SR. Fuzzy modeling of fuel cell based on mutual information between variables. Int J Hydrogen Energy 2010;35:3620-3631. [43] Hissel D, Péra MC, Kauffmann JM. Diagnosis of automotive fuel cell power generators. J Power Sources 2004;128:239246. [44] Saengrung A, Abtahi A, Zilouchian A. Neural network model for a commercial PEM fuel cell system. J Power Sources 2007;172:749-759. [45] Rumelhart DE, MacCLelland JL. Parallel Distributed Processing. Vol 1. Foundations MIT Press; 1986. [46] MacClelland JL, Rumelhart DE. Parallel Distributed Processing. Vol 2. sychological and Biological Models MIT Press; 1986. [47] Funahashi K. On the approximate realization of continuous mappings by neural networks. Neural Networks 1989;2:183192. [48] Dennis JE, Schnabel RB. Numerical Methods for Unconstrained Optimization and Nonlinear Equations Prentice Hall; 1983. [49] Asensio FJ, San Martín JI, Zamora I, Garcia-Villalobos J. Fuel cell-based CHP system modelling using Artificial Neural Networks aimed at developing techno-economic efficiency maximization control systems. Energy 2017;123:585-593. [50] Busquet S, Hubert CE, Labbé J, Mayer D, Metkemeijer R. A new approach to empirical electrical modelling of a fuel cell, an electrolyser or a regenerative fuel cell. J Power Sources 2004;134:41-48. [51] Kazim A. Exergy analysis of a PEM fuel cell at variable operating conditions. Energy Conversion and Management 2004;45:1949-1961. [52] Xin Kong, A. M. Khambadkone, Soy Kee Thum. A hybrid model with combined steady-state and dynamic characteristics of PEMFC fuel cell stack. Fourtieth IAS Annual Meeting Conference Record of the 2005 Industry Applications Conference, 2005 2005;3:1618-1625 Vol. 3. [53] X. Kong, W. Yeau, A. M. Khambadkone. ANN Modelling of Nonlinear Subsystem of a PEMFC Stack for Dynamic and Steady State Operation. IECON 2006 - 32nd Annual Conference on IEEE Industrial Electronics 2006:4322-4327. [54] K. C. S. Stanton, Jih-Sheng Lai. A Thermally Dependent Fuel Cell Model for Power Electronics Design. 2005 IEEE 36th Power Electronics Specialists Conference 2005:1647-1651. [55] Zhang Z, Huang X, Jiang J, Wu B. An improved dynamic model considering effects of temperature and equivalent internal resistance for PEM fuel cell power modules. J Power Sources 2006;161:1062-1068. [56] Miansari M, Sedighi K, Amidpour M, Alizadeh E, Miansari M. Experimental and thermodynamic approach on proton exchange membrane fuel cell performance. J Power Sources 2009;190:356-361. [57] Youssef ME, Al-Nadi KE, Khalil MH. Lumped Model for Proton Exchange Membrane Fuel Cell (PEMFC). Int J Electrochem Sci 2010;5:267-277. [58] Sadiq Al-Baghdadi MAR, Shahad Al-Janabi HAK. Parametric and optimization study of a PEM fuel cell performance using three-dimensional computational fluid dynamics model. Renewable Energy 2007;32:1077-1101. [59] Kim H, Cho CY, Nam JH, Shin D, Chung T. A simple dynamic model for polymer electrolyte membrane fuel cell (PEMFC) power modules: Parameter estimation and model prediction. Int J Hydrogen Energy 2010;35:3656-3663. [60] da Fonseca R, Bideaux E, Gerard M, Jeanneret B, Desbois-Renaudin M, Sari A. Control of PEMFC system air group using differential flatness approach: Validation by a dynamic fuel cell system model. Appl Energy 2014;113:219-229. [61] San Martín I, Ursúa A, Sanchis P. Modelling of PEM Fuel Cell Performance: Steady-State and Dynamic Experimental Validation. Energies 2014;7:670-700. [62] Salim R, Nabag M, Noura H, Fardoun A. The parameter identification of the Nexa 1.2 kW PEMFC's model using particle swarm optimization. Renewable Energy 2015;82:26-34. [63] Barzegari MM, Dardel M, Ramiar A, Alizadeh E. An investigation of temperature effect on performance of dead-end cascade H2/O2 PEMFC stack with integrated humidifier and separator. International Journal of Hydrogen Energy 2016;41:3136-3146. [64] Asensio FJ, San Martín JI, Zamora I, Garcia-Villalobos J. Fuel cell-based CHP system modelling using Artificial Neural Networks aimed at developing techno-economic efficiency maximization control systems. Energy 2017;123:585-593. [65] Asensio FJ, San Martín JI, Zamora I, Oñederra O. Model for optimal management of the cooling system of a fuel cellbased combined heat and power system for developing optimization control strategies. Applied Energy 2018;211:413-430. [66] Chen X, Li W, Gong G, Wan Z, Tu Z. Parametric analysis and optimization of PEMFC system for maximum power and efficiency using MOEA/D. Appl Therm Eng 2017;121:400-409. [67] Kheirandish A, Motlagh F, Shafiabady N, Dahari M, Khairi Abdul Wahab A. Dynamic fuzzy cognitive network approach for modelling and control of PEM fuel cell for power electric bicycle system. Applied Energy 2017;202:20-31. [68] Ziogou C, Voutetakis S, Georgiadis MC, Papadopoulou S. Model predictive control (MPC) strategies for PEM fuel cell systems – A comparative experimental demonstration. Chemical Engineering Research and Design 2018;131:656-670. [69] Chen J, Huang L, Yan C, Liu Z. A dynamic scalable segmented model of PEM fuel cell systems with two-phase water flow. Mathematics and Computers in Simulation 2018. [70] Sankar K, Jana AK. Nonlinear multivariable sliding mode control of a reversible PEM fuel cell integrated system. Energy Conversion and Management 2018;171:541-565. [71] Amphlett JC, Mann RF, Peppley BA, Roberge PR, Rodrigues A. A model predicting transient responses of proton exchange membrane fuel cells. Journal of Power Sources 1996;61:183-188. [72] Rowe A, Li X. Mathematical modeling of proton exchange membrane fuel cells. J Power Sources 2001;102:82-96. [73] Ticianelli EA, Derouin CR, Srinivasan S. Localization of platinum in low catalyst loading electrodes to to attain high power densities in SPE fuel cells. Journal of Electroanalytical Chemistry and Interfacial Electrochemistry 1988;251:275-295. [74] Ticianelli EA, Beery JG, Srinivasan S. Dependence of performance of solid polymer electrolyte fuel cells with low platinum loading on morphologic characteristics of the electrodes. J Appl Electrochem 1991;21:597-605. [75] Djilali N, Lu D. Influence of heat transfer on gas and water transport in fuel cells. International Journal of Thermal Sciences 2002;41:29-40. [76] Bernardi DM, Verbrugge MW. Mathematical model of a gas diffusion electrode bonded to a polymer electrolyte. AIChE J 1991;37:1151-1163. [77] Ticianelli EA, Derouin CR, Redondo A, Srinivasan S. Methods to Advance Technology of Proton Exchange Membrane Fuel Cells. Journal of The Electrochemical Society 1988;135:2209-2214. [78] Xue X, Tang J, Smirnova A, England R, Sammes N. System level lumped-parameter dynamic modeling of PEM fuel cell. J Power Sources 2004;133:188-204. [79] Shan Y, Choe S. A high dynamic PEM fuel cell model with temperature effects. J Power Sources 2005;145:30-39. [80] Ceraolo M, Miulli C, Pozio A. Modelling static and dynamic behaviour of proton exchange membrane fuel cells on the basis of electro-chemical description. J Power Sources 2003;113:131-144. [81] Pukrushpan JT, Stefanopoulou AG, Peng H. Modeling and control for PEM fuel cell stack system. Proc Am Control Conf 2002;4:3117-3122. [82] Caisheng Wang, M. H. Nehrir, S. Shaw. Dynamic models and model validation for PEM fuel cells using electrical circuits. IEEE Power Engineering Society General Meeting, 2005 2005:2115 Vol. 3. [83] Kjelstrup S, Røsjorde A. Local and Total Entropy Production and Heat and Water Fluxes in a One-Dimensional Polymer Electrolyte Fuel Cell. J Phys Chem 2005;109:9020-9033. [84] Weber AZ, Newman J. Coupled Thermal and Water Management in Polymer Electrolyte Fuel Cells. Journal of The Electrochemical Society 2006;153:A2205-A2214. [85] Park S, Choe S. Dynamic modeling and analysis of a 20-cell PEM fuel cell stack considering temperature and two-phase effects. J Power Sources 2008;179:660-672. [86] Shah AA, Walsh FC. A model for hydrogen sulfide poisoning in proton exchange membrane fuel cells. Journal of Power Sources 2008;185:287-301. [87] Mohtadi R, Lee W-, Van Zee JW. The effect of temperature on the adsorption rate of hydrogen sulfide on Pt anodes in a PEMFC. Applied Catalysis B: Environmental 2005;56:37-42. [88] Hu P, Cao G, Zhu X, Li J. Modeling of a proton exchange membrane fuel cell based on the hybrid particle swarm optimization with Levenberg–Marquardt neural network. Simulation Modelling Practice and Theory 2010;18:574-588. [89] Sharifi Asl SM, Rowshanzamir S, Eikani MH. Modelling and simulation of the steady-state and dynamic behaviour of a PEM fuel cell. Energy 2010;35:1633-1646. [90] Tiss F, Chouikh R, Guizani A. Dynamic modeling of a PEM fuel cell with temperature effects. Int J Hydrogen Energy 2013;38:8532-8541.