A 2D Multi-Layer Model to Study the External Magnetic Field Generated by a Polymer Exchange Membrane Fuel Cell

: An original innovative two-dimensional (2D) multi-layer model based on the Maxwell– Fourier method for the diagnosis of a polymer exchange membrane (PEM) fuel cell (FC) stack is presented. It is possible to determine the magnetic ﬁeld distribution generated around the PEMFC stack from the (non-)homogenous current density distribution inside the PEMFC stack. Analysis of the magnetic ﬁeld distribution can indicate whether the FC is healthy or faulty. In this way, an explicit, accurate and fast analytical model can allow the health state of an FC to be studied. To evaluate the capacity and the efﬁciency of the 2D analytical model, the distribution of local quantities (i.e., magnetic vector potential and magnetic ﬁeld) in a PEMFC stack has been validated with those obtained by the 2D ﬁnite-element analysis (FEA). The comparisons demonstrate excellent results both in terms of amplitude and waveform. The validation of this 2D analytical model is essential for the subsequent generation of an inverse model useful for the diagnosis of a PEMFC.


Context of the Paper
Currently, the fuel cell (FC) system presents the advantage of combining high-energy efficiency (e.g., FC efficiency about 50% and heat engine efficiency about 32%) and environmental friendliness compared to other more conventional technologies [1]. This presents an interesting combination for conventional power generation systems, especially in stationary and portable FC applications [2][3][4]. Furthermore, a polymer exchange membrane fuel cell (PEMFC), which is the most commonly used FC in transport applications, is an interesting alternative to the internal combustion engine [5].
In order to transform fuel cells in to a mass-market technology, many challenges must be overcome. For that, it is necessary to estimate and comprehend the possible performance of PEMFC in automotive applications as analyzed in [6]. Hydrogen is a promising energy resource. In order to obtain an overview, various studies are based on economics and storage, while recommendations [7] or a critical view can be given on technologies, applications, trends and challenges [8].
The main application of FCs is in the automotive industry. A comprehensive review and research on FC electrical vehicles are presented in [9]. Many topics such as topologies, power electronic converters, energy management methods, technical challenges, marketing and future aspects are discussed.
Some research studies are conducted on the degradation mechanism and impacts on FC degradation with a focus in hybrid transport applications [10]. Similarly, FC diagnosis methods for embedded automotive applications are detailed in [11].
Moreover, studies were performed on the external magnetic field determination by imposing the (non-)homogeneous current density distribution inside a PEMFC with 3D numerical [35][36][37] and experimental [35] validations. It is interesting to note that a ferromagnetic circuit placed around a PEMFC emulator permits concentration of the external magnetic field and, consequently, the detection amplification of an FC failure [35]. The actual methods to determine the current density distribution can be globally described as redundant, time-consuming and have variable accuracy.

Objective of the Paper
In the literature, to the best of the authors' knowledge, there is no 2D multi-layer model based on the Maxwell-Fourier method for the diagnosis of a PEMFC stack. This method is based on the formal resolution of Maxwell's magnetostatic equations in Cartesian coordinates (x, y) by using the separation of variables method, the Fourier's series and the superposition principle. This purely analytical multi-layer model could give a perfect theoretical study of the magnetic field distribution with a rapid computation time. Moreover, this model could present more precision than the usual methods. Here, the magneto-tomography study is realized on a real bipolar plate of a PEMFC stack (without ferromagnetic circuit placed around the PEMFC stack) in order to diagnose the health state of this one, viz., (i) healthy FC, and (ii) FC failure. As reported in [35][36][37], the developed analytical model will allow the magnetic field distribution generated around the PEMFC to be determined by imposing the (non-)homogeneous current density distribution inside the FC. The reverse study, which is the second step of this scientific work, could then be easily performed on the 2D multi-layer model. However, this is beyond the scope of this paper. To evaluate the efficiency of the proposed analytical model, the distribution of local quantities (i.e., magnetic vector potential and magnetic field) for a healthy FC and an FC failure has been validated with those obtained by the 2D FEA [38]. The comparisons Mathematics 2022, 10, 3883 3 of 15 demonstrate excellent results both in terms of amplitude and waveform. The validation of this 2D analytical model is essential for the subsequent generation of an inverse model. This kind of model, which has never been previously developed for FC diagnosis application, must first be validated (e.g., by numerical validation) in order to finally obtain a powerful diagnostic tool of the FC operating mode.

Mathematical
, with a related constant current of I f c , consists of elements in parallel, i.e., compression of two end-plates, two current collectors and several unit cells. Each unit cell is composed of two bipolar plates with an inactive zone (seal) and an electrochemical active cell (of dimension w a × h a ), two gas diffusion layers and MEA. In electrochemical active cells, at the MEA level, the current density J = {0; 0; J z }, which is perpendicular to the cells with a constant value along the z-axis, is directly related to the flow of hydrogen protons travelling through the electrolyte membrane. Outside of the electrochemical active cell, J is neglected. An FC failure is characterized by a non-homogeneous distribution of J [12]. The considered bipolar plate for the diagnosis of a PEMFC stack, with the geometrical and physical parameters, is illustrated in Figure 1a. This system is surrounded by a vacuum via an infinite box (of dimension w × h). of this paper. To evaluate the efficiency of the proposed analytical model, the distribution of local quantities (i.e., magnetic vector potential and magnetic field) for a healthy FC and an FC failure has been validated with those obtained by the 2D FEA [38]. The comparisons demonstrate excellent results both in terms of amplitude and waveform. The validation of this 2D analytical model is essential for the subsequent generation of an inverse model. This kind of model, which has never been previously developed for FC diagnosis application, must first be validated (e.g., by numerical validation) in order to finally obtain a powerful diagnostic tool of the FC operating mode. , which is perpendicular to the cells with a constant value along the z-axis, is directly related to the flow of hydrogen protons travelling through the electrolyte membrane. Outside of the electrochemical active cell, J is neglected. An FC failure is characterized by a non-homogeneous distribution of J [12]. The considered bipolar plate for the diagnosis of a PEMFC stack, with the geometrical and physical parameters, is illustrated in Figure 1a. This system is surrounded by a vacuum via an infinite box (of dimen-  Table 1 for the various parameters).  Table 1 for the various parameters). Table 1. Geometrical and physical parameters of a PEMFC stack.

Symbol
Quantity Values

. Normal and Faulty Conductive Segments
To simulate (non-)homogeneous distribution of J, this electrochemical active cell of bipolar plates can be decomposed into Ny × Nx regular conductive segments (of dimension S cs = w cs × h cs with w cs = w a /Nx and h cs = h a /Nx), as shown in Figure 1b, where Ny ∈ N * and Nx ∈ N * are the segmentation number in the yand x-axis, respectively. In Cartesian coordinates (x, y), the starting and ending points of conductive segments in both axes can be defined by with j ∈ {N * , Nx} and i ∈ {N * , Ny}. Any conductive segments (i, j) may be faulty (not supplied) or normal (supplied). For a non-homogeneous distribution of J reflecting an FC failure, some conductive segments are not supplied (i.e., the current are assumed to be null in the conductive segment). For a homogenous distribution of J reflecting a healthy FC, all conductive segments are supplied by direct current I cs symbolized by ⊗ equal to with N n cs the number of conductive segments in a normal condition defined by where M cs i,j represents the element (i, j) of the power supply matrix of the conductive segments M cs (of dimension Ny × Nx), viz., 2.1.3. Modeling of Spatial Current Density Distribution in any Conductive Segments (i, j) Figure 2 illustrates the waveform of the spatial distribution of J zij in any conductive segments (i, j), which can be expressed as a Fourier's series with x ∈ [0, w]: with where J cs max = I cs /S cs is the maximum current density in each conductive segment, β k = kπ/w is the spatial frequency (or periodicity) of J zij with k ∈ {N * , K max } the spatial harmonic orders in which K max is the finite number of spatial harmonics terms.
is the maximum current density in each conductive segment, k k w β π = is the spatial frequency (or periodicity) of zij J with { } max k ,K * ∈  the spatial harmonic orders in which max K is the finite number of spatial harmonics terms.

List of Assumptions
The 2D magnetic field distribution generated around the PEMFC stack has been studied from multi-layer model by solving Maxwell's magnetostatic equations in Cartesian coordinates ( ) • the skin depth effect is not considered; • many industrial PEMFC stacks use graphite bipolar plates and non-magnetic stainless steels for rods and other mechanical parts [35]. Thus, the absolute magnetic permeability μ of all components used in bipolar plates of PEMFC stack are very close to the vacuum magnetic permeability H m

Problem Discretization in Regions
As shown in Figure 3, the bipolar plate for the diagnosis of a PEMFC stack is divided into three regions x 0, w ∀ ∈    , viz., • Vacuum above the electrochemical active cell: Region 1 (R1) for • The electrochemical active cell: Region 2 (R2) for e s Ny 1 Figure 2. Waveform of the spatial distribution of J zij in any conductive segments (i, j).

List of Assumptions
The 2D magnetic field distribution generated around the PEMFC stack has been studied from multi-layer model by solving Maxwell's magnetostatic equations in Cartesian coordinates (x, y). In this analysis, the simplifying assumptions of the analytical model are: • the end-effects are neglected (i.e., the magnetic variables are independent of z); • the magnetic vector potential and current density have only one component along the z-axis, i.e., A = {0; 0; A z } and J = {0; 0; J z }; • the skin depth effect is not considered; • many industrial PEMFC stacks use graphite bipolar plates and non-magnetic stainless steels for rods and other mechanical parts [35]. Thus, the absolute magnetic permeability µ of all components used in bipolar plates of PEMFC stack are very close to the vacuum magnetic permeability µ v = µ 0 = 4π × 10 −7 H/m.

Problem Discretization in Regions
As shown in Figure 3, the bipolar plate for the diagnosis of a PEMFC stack is divided into three regions ∀x ∈ [0, w], viz., • Vacuum above the electrochemical active cell: Region 1 (R1) for y ∈ y s 1 , h with The electrochemical active cell: Region 2 (R2) for y ∈ y e Ny , y s 1 , which is divided in the y-axis into sub-regions i (R2 i ) for y ∈ y e i , y s i , with µ 2 = µ v = µ 0 ; • Vacuum below the electrochemical active cell: Region 3 (R3) for y ∈ 0, y e Ny with

Governing Partial Differential Equations (PDEs) in Cartesian Coordinates
In quasi-stationary, the magnetostatic Maxwell's equations are represented by Maxwell-Ampère [39]

Governing Partial Differential Equations (PDEs) in Cartesian Coordinates
In quasi-stationary, the magnetostatic Maxwell's equations are represented by Maxwell-Ampère [39] ∇ × H = J (with J = 0 for the no − load operation), Maxwell-Thomson and the magnetic material equation where B, H and Mr are the magnetic flux density, the magnetic field and the remanent (or residual) magnetization vector (with Mr = 0 for other materials or Mr = 0 for the magnets), respectively. Using (6)~(8) with the general assumptions of the study, in Cartesian coordinates (x, y), the general PDEs of magnetostatic in terms of A = {0; 0; A z } inside an isotropic and uniform material (i.e., µ = C ste ) can be expressed by the Laplace's equation in (R1) and in (R3) and the Poisson's equation in (R2 i ) where J z2 i is the spatial distribution of J along the x-axis for each i, i.e., in (R2 i ), which is defined by For example, with Ny × Nx = 2 × 10 conductive segments, Figure 4 illustrates the waveforms of J z2 i for each i in (R2 i ) obtained by (10) Mathematics 2022, 10, x FOR PEER REVIEW 7 of 17 (a) (b)  Using (7b), the components of B = B x ; B y ; 0 whatever the region can be deduced from A z by From (8), and according to the materials assumptions, the components of H = H x ; H y ; 0 are defined by H = B/µ 0 .

Boundary Conditions (BCs)
In electromagnetic field problems, the general solutions of regions depend on boundary conditions (BCs) at the interface between two surfaces, which are defined by the continuity of parallel magnetic field H and A [40]. On the infinite box [see Figure 1], A z satisfies the Dirichlet's BC, viz., A z = 0. Figure 5 represents the respective BCs at the interface between the various regions.

General Solution of Magnetic Field in Each Region
Using the separation of variables method, the 2D general solution of A in each region can be defined by Fourier's series. The unknown coefficients (or integration constants) of series are determined analytically from a linear matrix system satisfying the BCs of Figure 5. To simplify the formal resolution of the Cramer's system, the superposition principle on each sub-region i (R2i) has been applied.

General Solution of Magnetic Field in Each Region
Using the separation of variables method, the 2D general solution of A in each region can be defined by Fourier's series. The unknown coefficients (or integration constants) of series are determined analytically from a linear matrix system satisfying the BCs of Figure 5. To simplify the formal resolution of the Cramer's system, the superposition principle on each sub-region i (R2 i ) has been applied.
In (R1) for y ∈ y s 1 , h and ∀x ∈ [0, w], the magnetic vector potential A z1 , which is a solution of (9a) satisfying the various BCs, is defined by where with in which a1 k = −sh(β k ·h cs )/ch(β k ·h cs ), (20) In (R3) for y ∈ 0, y e Ny and ∀x ∈ [0, w], the magnetic vector potential A z3 , which is a solution of (9b) satisfying the various BCs, is defined by (24) and the components of H 3 (i.e., H x3 and H y3 ) by where with In (R2 i ) for y ∈ y e i , y s i and ∀x ∈ [0, w], the magnetic vector potential A z2 i , which is a solution of (9c) satisfying the various BCs, is defined by with

Problem Description
The real bipolar plate of a PEMFC stack is illustrated in Figure 6. The electrochemical active cell is presented with the oxygen/hydrogen channels and surrounded by an inactive zone (seal). The geometrical and physical parameters of the PEMFC stack are reported in Table 1. conductive segments is considered. This discreti the electrochemical active cell to study health of the FC is illustrated in Figure 7. To simulate the (non-)homogeneous distribution of J, this electrochemical active cell is decomposed into Ny × Nx regular conductive segments. For example, the discretization with Ny × Nx = 25 × 20 = 500 conductive segments is considered. This discretization of the electrochemical active cell to study health of the FC is illustrated in Figure 7. To simulate the (non-)homogeneous distribution of J , this electrochemical active cell is decomposed into Ny Nx × regular conductive segments. For example, the discretization with Ny Nx 25 20 500 × = × = conductive segments is considered. This discretization of the electrochemical active cell to study health of the FC is illustrated in Figure 7. Such an important discretization is necessary due to the channels' size (viz.,~0.8 mm) as shown in Figure 6. The total current flowing through the PEMFC stack (of surface S f c = 196 cm 2 ) is considered constant and equal to I f c = 70 A. With an electrochemical active cell surface of S a = 100 cm 2 , for a healthy FC [see Figure 7a], this corresponds to J cs max = 0.7 A/cm 2 . An example of an FC failure, which can be a consequence of different operating problems, is exposed in Figure 7b. Considering 500 initial conductive segments, 144 segments are faulty (i.e., J cs max = 0 A/cm 2 ), so 356 conductive segments present a current density of J cs max = 0.983 A/cm 2 . To evaluate the capacity and the efficacy of the 2D analytical model, the two operating conditions (i.e., normal and faulty condition) are simulated and compared to 2D FEA. The geometrical and physical parameters described previously are simulated on FEMM software [38]. The FE calculations are performed with the same assumptions as the analytical model [see Section 2.2]. The infinite box surface has been imposed to S = 3·S f c = 1764 cm 2 (viz., w = 3·w f c and h = 3·h f c ) so that the equipotential lines of A z would not be perturbed by the Dirichlet's CB imposed on the edges. Figure 8 presents the mesh of the system (viz., 107,385 nodes) and different paths with their dimensions added around the PEMFC stack (in red, rated from A to D) in order to evaluate the magnetic field distribution. In this comparison, the analytical solution of A and H in each region have been computed with K max = 100.   Table 1 for the various parameters).

Healthy FC
The 2D multi-layer model is implemented in order to obtain values of z A in the bipolar plate for the diagnosis of a PEMFC stack. Figure 9 presents the equipotential lines (   Table 1 for the various parameters).

Healthy FC
The 2D multi-layer model is implemented in order to obtain values of A z in the bipolar plate for the diagnosis of a PEMFC stack. Figure 9 presents the equipotential lines (≈30 lines) of A z (or vector plot of H) for a healthy FC with the 2D analytical model and 2D FEA. As can be seen, a perfect accuracy between analytical and numerical results allows a first validation of the multi-layer model based on the Maxwell-Fourier method. Due to the homogenous distribution of J, and for a square FC, the field lines are symmetrical and circular along the centre of the bipolar plate. in the different regions. It is interesting to note that the waveform and the amplitude of magnetic field components are exactly the same in the four paths (i.e., A to D in Figure 8) independently of the axis.

FC Failure
With a non-homogeneous distribution of J reflecting an FC failure, the generated magnetic field is impacted. In Figure 10, the deformation of the equipotential lines of A z is clearly observed. It is also correctly predicted analytically. On the left, a distortion of the field lines (i.e., asymmetry of field lines) is widely visible due to the presence of a major defect on the right [see Figure 7b].    Figure 11 represents the distribution of the components of H = H x ; H y ; 0 for a healthy FC computed by the 2D multi-layer model and verified by the 2D FEA on the various paths. Highly accurate results are achieved between the analytical and numerical model, regardless of the paths and the components of H = H x ; H y ; 0 in the different regions. It is interesting to note that the waveform and the amplitude of magnetic field components are exactly the same in the four paths (i.e., A to D in Figure 8) independently of the axis.

FC Failure
As for the case of a healthy FC, Figure 12 shows the distribution of the components of

FC Failure
As for the case of a healthy FC, Figure 12 shows the distribution of the components of H = H x ; H y ; 0 for an FC failure obtained analytically and numerically on the different paths. Irrespective of the paths and regions, similar and accurate magnetic field results also permit validation of the 2D analytical model. A magnetic field variation reflects a failure in the FC. Around the most damaged part, the greatest decrease in magnetic field is observed (viz., on the path B and C), as shown in Figure 12b,c. As the overall current value I f c is maintained, the current density J cs max increases in the rest of the FC [see Figure 7]. Thus, the magnetic field variation is less observable for small surface defaults.