PRES22_0231.docx DOI: 10.3303/CET2294182 Paper Received: 15 April 2022; Revised: 13 July 2022; Accepted: 14 July 2022 Please cite this article as: Antonano H.C., Panganiban J.M.F., Yu J.V.T., Castro M.T., Ocon J.D., 2022, Multiphysics Modeling of High Temperature Planar Sodium-Sulfur Batteries, Chemical Engineering Transactions, 94, 1093-1098 DOI:10.3303/CET2294182 CHEMICAL ENGINEERING TRANSACTIONS VOL. 94, 2022 A publication of The Italian Association of Chemical Engineering Online at www.cetjournal.it Guest Editors: Petar S. Varbanov, Yee Van Fan, Jiří J. Klemeš, Sandro Nižetić Copyright © 2022, AIDIC Servizi S.r.l. ISBN 978-88-95608-93-8; ISSN 2283-9216 Multiphysics Modeling of High Temperature Planar Sodium Sulfur Batteries Hezekiah C. Antonano, Jay Mark F. Panganiban, John Vincent T. Yu, Michael T. Castro, Joey D. Ocon* Laboratory of Electrochemical Engineering, Department of Chemical Engineering, University of the Philippines Diliman, Quezon City 1101 jdocon@up.edu.ph Sodium-sulfur (NaS) batteries are a promising energy storage technology that features high energy density, high cycle life, and no self-discharge. One of the cell design considerations that can affect the performance and construction of NaS batteries is cell geometry. While the planar geometry has advantages in power output, cell packing, ease of assembly, and thermal management, it has thermo-mechanical issues in sealing, which makes it less commercialized than the tubular geometry. In this work, the first multiphysics model of a NaS cell with a planar geometry was developed by extending an existing multiphysics model for a NaS cell with a tubular geometry. A previously developed multiphysics model for the tubular NaS cell was first adapted into COMSOL Multiphysics®. The model’s results were then validated against experimental discharge profiles and surface temperatures. After this, a planar NaS cell was simulated using the dimensions of an experimental planar cell, with the mass, heat, and charge transfer equations and parameters previously used in the tubular multiphysics model. The results from the planar model were then validated through comparison to experimental discharge profiles for a planar cell. The developed model can be used in future research and development of NaS batteries by comparing performance parameters of the planar geometry such as discharge profiles, energy densities, and temperatures with a multiphysics model of the tubular NaS battery. 1. Introduction A major drawback of renewable energy such as wind and solar power is the inherent lack of stability of these resources (Ibrahim et al., 2008). As such, energy storage systems may be employed to store energy during peak times and release it when renewables are unable to meet the demand (Beaudin et al., 2010). The sodium- sulfur battery (NaS) is a promising energy storage technology with several key features such as a high energy density, high efficiency, low cost, versatility, and a high degree of technical maturity (Oshima et al., 2004). For energy storage applications, high-temperature NaS batteries are of greater interest than their room-temperature counterpart due to its higher capacity retention and better cycling properties (Nikiforidis et al., 2019). In a high-temperature NaS battery, molten sodium and sulfur serve as the battery’s anode and cathode (Oshima et al., 2004). Both are separated by beta-alumina solid electrolyte (BASE), which permits the selective transport of Na+ ions from the anode to the cathode. The cell may be built with either the tubular or planar geometry. The tubular geometry provides better mechanical stability for the BASE and the design can be easily scaled into high-capacity cells (Kim et al., 2016). On the other hand, the planar cell geometry has superior power output, cell packing, ease of assembly, and thermal management (Lu et al., 2016), but suffers from thermo-mechanical issues, which pose safety risks (Nikiforidis et al., 2019). The performance of different battery designs can be evaluated via multiphysics modeling. This technique incorporates various physical phenomena such as mass transfer, charge transfer, and heat transfer to describe systems where the properties of the entities in the system are coupled to other variables (Dickinson et al., 2014). Multiphysics modeling comes with several advantages, such as saving resources from constructing physical experimental setups while still allowing a comprehensive study and analysis of the system. It can also provide additional data and analysis that would otherwise be too difficult to obtain from experimental setups. 1093 Multiphysics models for the tubular geometry of the high temperature NaS cell have been previously developed. For instance, Schaefer et al. (2020) developed a 1D axisymmetric transient non-isothermal model of a tubular NaS cell. This study was based on earlier models developed by Kawamoto and Kusakabe (1989) and Kawamoto (1991). However, no models have been developed for planar geometry thus far. These models could be used to study the performance of high-temperature planar NaS batteries and compare it to the performance of a similarly sized tubular NaS cell. In this work, a multiphysics model of a high-temperature planar NaS battery was formulated based on existing models for tubular NaS cells. First, a multiphysics model for a tubular NaS battery was replicated in COMSOL Multiphysics®. The model was then validated by comparing operating voltages and surface temperatures during charging and discharging operations with experimental data in the work of Schaefer et al. (2020). Afterwards, the model was extended for planar geometries, and then validated with experimental data from Kim et al. (2016). 2. Methodology 2.1 Model description The discharge reaction for a sodium-sulfur battery is described by Eq(1) and Eq(2). The sodium metal in the anode liberates an electron to form Na+. The ion is then transported across the BASE and into the cathode, where it reacts with sulfur to form a polysulfide compound Na2Sx. 2Na → 2Na+ + 2e− (1) 2Na+ + 𝑥𝑥S + 2e− → Na2S𝑥𝑥 (2) The geometry of the tubular and planar batteries are presented in Figure 1. The tubular model was described using a 1D axisymmetric geometry, wherein the positive radial direction goes from the center of the tubular cell towards a radially outward direction until the outer surface of the molten sulfur cathode. The planar model assumes that 1D transport phenomena occur along its thickness. The positive 𝑥𝑥 direction in the planar cell starts from the outer surface of the sodium anode to the outer surface of the sulfur cathode. Figure 1: Schematic diagram and dimensions of tubular (a) and planar (b) geometries of the NaS cell The tubular model was developed based on the non-isothermal model of a NaS battery created by Schaefer et al. (2020), while the planar model is based on the experimental planar battery created and reported by Kim et al. (2016). These works were selected since they contained various data such as its operating voltage at different depth of discharge which can help in validating the accuracy of the model. The model was implemented in COMSOL Multiphysics® using the Battery with Binary Electrode, Heat Transfer in Solids interfaces, and Electrochemical Heating modules. The mesh was the default one generated by COMSOL Multiphysics®, and contained 102 vertices in the tubular model and 36 vertices in the planar model. Physical memory and virtual memory used during the computation were around 750 MB and 1450 MB, with a computation time of around one minute. 2.2 Governing equations Schaefer et al. (2020) developed a transient nonisothermal model of a tubular NaS cell based on the cell detailed in the works of Kawamoto and Wada (1987), Kawamoto and Kusakabe (1989), and Kawamoto (1991). The cathode melt of the NaS cell was considered by Newman and Thompson (1989) to be a binary melt consisting of sodium cations and sulfide anions in a neutral sulfur solvent. As such, the electrochemical equations for 1094 concentrated solution theory are used for this model. The electrochemical model includes the volumes of the sulfur electrode and the BASE only, with the former being models as a porous electrode and the latter modeled as a separator. While most of the parameters used are obtained directly from the work of Schaefer et al. (2020), the following parameters were modified for this work. The electrolyte conductivity, electrolyte salt diffusivity, and electrode conductivity of the cathode were modified with the porosity in the polysulfide melt 𝑝𝑝 as an empirical fitting parameter for the model. 𝜎𝜎l,eff = 𝑝𝑝𝑁𝑁𝜎𝜎l (3) 𝐷𝐷l,eff = 𝑝𝑝𝑁𝑁𝐷𝐷l (4) 𝜎𝜎s,eff = (1 − 𝑝𝑝)𝑁𝑁𝜎𝜎s (5) The transport number was modified from the work of Risch and Newman (1988). The empirical constants 𝛽𝛽1 and 𝛽𝛽3 are given in the same work. 𝑡𝑡+ = � 1 𝑀𝑀 ≥ 0.71 −[(𝛽𝛽3 𝛽𝛽1⁄ ) + 𝑋𝑋E] (1 − 𝑋𝑋E)⁄ 𝑀𝑀 < 0.71 (6) A flux condition is imposed on the sodium/BASE interface to simulate the anode compartment as a reservoir where sodium ions may enter, or exit based on the electrolyte current. Electric insulation is imposed on the sulfur/BASE interface. An electrolyte potential boundary condition was used to set the potential at the sodium/BASE interface to 0 as a reference potential. Finally, an electrode current boundary condition was used to set an applied current as the outer face of the cathode compartment. This applied current is varied in different simulations. The thermal model developed by Schaefer et al. (2020) considered the heat generation in each component through the electrochemical phenomena occurring in the component as well as heat transfer by conduction through the model and convection and radiation at the surfaces of the cell. The model only considered 1D heat transfer. The heat generation in the cell is the result of Joule heating, activation losses, and reversible heating due to reactions in the electrochemistry interface. COMSOL’s electrochemical heating multiphysics node was used to implement this for the model. The physical and chemical constants such as mass, density and other parameters are obtained from the work of Schaefer et al. (2020). Values of the polysulfide melt porosity’s exponent were studied from 0.9 to 1.2 in intervals of 0.1 to find the value that best fits the experimental data. The electrode volume fraction for the tubular cell was obtained from the work of Kawamoto and Wada (1987) while the porosity of the polysulfide melt was retained from the work of Schaefer et al. (2020). The heat transfer coefficient was estimated by fitting it to reduce the deviations with the experimentally determined surface temperature of the cell over time. Nevertheless, the estimated value of 45 W/m2·K was comparable to the 30 W/m2·K value used by Schaefer et al. (2020). Table 1: Input parameters to the model Parameters Value Reference Transport parameter correction factor, 𝑝𝑝 0.35 (Schaefer et al., 2020) Exponent to the correction factor, 𝑁𝑁 0.9 to 1.2 This work Electrode volume fraction, 𝜀𝜀s 0.0464 (Tubular) 0.2784 (Planar) (Kawamoto and Wada, 1987) This work Heat transfer coefficient, ℎ [W/m2·K] 45 This Work 2.3 Model validation To compare model results with experimental data, the model inputs listed in Table 2 were used to simulate the conditions in the experimental set-ups. The convention used for the current is (+) for discharge and (-) for charge. For comparing the voltage vs state of discharge (SOD) data and temperature vs time data with Schaefer et al. (2020) and Kim et al.’s (2016) experimental data, the simulation was run at a single current for a fixed amount of time. For comparing the voltage vs time data with Kawamoto and Kusakabe’s (1989) experimental data, the simulation was operated for 3 charge/discharge cycles to achieve a steady-state condition at high current density operation, with the third cycle’s data being used for the comparison. 1095 Table 2: Model inputs used in simulating the experimental set-ups Experimental Data Source Current [A] Runtime [h] Initial SOD Initial Temp. [°C] Schaefer et al. 3.0 / -3.0 16.0 0.18 / 1.00 335 (2020) 6.0 / -6.0 7.8 0.18 / 1.00 335 9.0 / -9.0 4.5 0.29 / 1.00 335 Kawamoto and Kusakabe 11.0 / -11.0 3.8 / 3.9 1.00 335 (1989) 16.5 / -22.0 2.6 / 1.9 1.00 335 22.0 / -33.0 1.45 / 1.0 1.00 335 Kim et al. (2016) 1.2 4.5 0.138 350 -1.0 4.5 0.781 350 3. Results and discussion 3.1 Model validation Figure 2 compares the model results to the experimental data from the work of Schaefer et al. (2020). The model appears to be able to approximate the experimental data. As the fitting parameter increases, the deviation of the charge and discharge curves from the equilibrium potential also increases. Figure 2: Comparison between model and experimental data for voltage with respect to SOD at (a) 3 A, (b) 6 A, and (c) 9 A at different values of the porosity exponent 𝑁𝑁 In Figure 3, the model can provide results which are able to follow the trend of the experimental data. However, larger deviations may be noted for high current density applications. According to Kawamoto and Kusakabe (1989), discrepancies at low SOD occur due to the formation of the sulfur insulating film at the surface of the BASE and discrepancies at high SOD may be attributed to the partial formation of Na2S2, which precipitates out of the polysulfide melt at the operating temperature (Kawamoto and Kusakabe, 1989). Figure 3: Comparison between model and experimental data for voltage with respect to time for high current density operation: (a) 11 A charge/11 A discharge, (b) 16.5 A charge/22 A discharge, and (c) 22 A charge/33 A discharge at different values of the porosity exponent 𝑁𝑁 In Figure 4, the model follows the trend of the experimental data, but the results are observed to have spikes and sudden drops. A spike may be observed at the transition zone between charge and discharge. A sudden drop is also observed at the start of the profiles as this is a transition point from discharging to charging. In general, the charging temperature is observed to be lower than the discharging temperature, which is to be 1096 expected as the reaction involved in the NaS cell is exothermic. Despite the ambient temperature being held constant at 335 °C (608.15 K), the temperature for high current density operation for the tubular model is very close to the safe temperature limit of NaS cells of 400 °C (673.15 K), with the highest current of 33 A discharge even surpassing this safety limit. Based on the comparison between the tubular model and the experimental data, a value of 𝑁𝑁 = 1 appears to be sufficient for the model. Figure 4: Comparison between model and experimental data for temperature with respect to time for high current density operation: (a) 11 A charge/11 A discharge, (b) 16.5 A charge/22 A discharge, and (c) 22 A charge/33 A discharge at different values of the porosity exponent 𝑁𝑁. Core temperature estimates were calculated in the work of Kawamoto and Kusakabe (1989) The model developed by Kawamoto and Kusakabe (1989), which formed part of the basis of the model developed by Schaefer et al. (2020), also exhibited similarly poor fit with experimental data. Kawamoto and Kusakabe (1989) calculated cell voltage as a function of the cell’s open circuit voltage and a pre-defined resistivity. Schaefer et al.’s (2020) model was able to improve the fit by obtaining the cell voltage as a consequence of 1D mass and charge transport, but otherwise still retained a noticeable error. As such, the inclusion of multidimensional transport phenomena, may be able to improve the fit. In Figure 5, the experimental data from the reference was compared to the model results (Kim et al., 2016). For the charge plot, the model is able to follow the experimental data well. However, for the discharge plot, the model’s discharge voltage has a relatively high curvature compared to the experimental data (Kim et al., 2016). This is able to verify that the planar multiphysics model developed is able to describe a real planar NaS cell. It would be ideal to compare the performance of the planar and tubular cells using the multiphysics models, but it should be noted that this is not a straightforward process. This is because the dimensions of the batteries for comparison (e.g., how to scale up/down the dimensions) must be carefully selected for a fair comparison between the planar and tubular cells. Figure 5: Comparison between the model and experimental data for voltage with respect to SOD for the planar geometry at different values of the porosity exponent 𝑁𝑁 4. Conclusions In this work, a multiphysics model for high temperature planar NaS batteries were formulated by extending models for tubular NaS batteries. A tubular NaS battery model was first implemented in COMSOL Multiphysics®, and then validated with experimental data by Schaefer et al. (2020). Afterwards, the tubular model was modified into a planar model by preserving the phenomenological equations while altering the geometry. The planar model was then validated with experimental data reported by Kim et al. (2016) 1097 While the models appear to be able to follow the trends in the experimental data, there are larger deviations at higher current densities. These may be resolved through the consideration of other phenomena, such as multidimensional mass, heat, and charge transport. In future work, the developed planar high temperature NaS battery model may be used to optimize the design of planar NaS batteries. Additionally, the developed planar multiphysics model may be compared to the existing tubular model in parameters such as charge-discharge voltage, maximum temperature, and energy density to identify the strengths and weaknesses of each geometry. Nomenclature 𝐷𝐷l – liquid phase diffusion coefficient, m2/s 𝑁𝑁 – porosity exponent, - 𝑝𝑝 – porosity of polysulfide melt, - 𝑡𝑡 – transport number, - 𝑋𝑋𝐸𝐸 – mole fraction of electrolyte in melt, - 𝛽𝛽1, 𝛽𝛽3 -empirical coefficients, V 𝜎𝜎l – liquid phase conductivity of melt, S/m 𝜎𝜎s – solid phase conductivity of melt, S/m Acknowledgments H.C.A., J.F.P., and J.V.Y. would like to thank the DOST-SEI Merit Scholarship Program. The authors would like to acknowledge the Department of Science and Technology (DOST) through the Niche Centers in the Regions for R&D (NICER) Program and The Commission on Higher Education – Philippine California Advanced Research Institutes (CHED-PCARI) through the CIPHER Project (IIID 2018-008). References Beaudin M., Zareipour H., Schellenberglabe A. Rosehart W., 2010. Energy storage for mitigating the variability of renewable electricity sources: An updated review. Energy for Sustainable Development, 14, 302–314. Dickinson E.J.F., Ekström H., Fontes E., 2014. COMSOL Multiphysics®: Finite element software for electrochemical analysis: A mini-review. Electrochemistry Communications, 40, 71–74. Ibrahim H., Ilinca A., Perronlu J., 2008. Energy storage systems: Characteristics and comparisons. Renewable and Sustainable Energy Reviews, 12, 1221–1250. Kawamoto H., 1991. Dynamic simulation of the charge-discharge characteristics of the sodium-sulphur cell. Journal of Applied Electrochemistry, 21, 409–414. Kawamoto H., Kusakabe Y., 1989. Performance and thermal behavior of sodium‐sulfur cell under high current density operations. Journal of The Electrochemical Society, 136, 1355–1361. Kawamoto H., Wada M., 1987. Two-dimensional distribution of electrochemical reaction rate in porous sulfur electrodes of NaS cells. Journal of The Electrochemical Society, 8, 280–285. Kim G., Park Y.C., Lee Y., Cho N., Kim C.S. Jung, K., 2016. The effect of cathode felt geometries on electrochemical characteristics of sodium sulfur (NaS) cells: Planar vs. tubular. Journal of Power Sources, 325, 238–245. Lu X., Li G., Meinhardt, K.D., Sprenkle V.L., 2016. Development of Na-beta alumina batteries at Pacific Northwest National Laboratory: From tubular to planar. Energy Storage Science and Technology, 5, 309– 316. Nikiforidis G., van de Sanden M.C.M., Tsampas M.N., 2019. High and intermediate temperature sodium-sulfur batteries for energy storage: Development, challenges and perspectives. RSC Advances, 9, 5649–5673. Oshima T., Kaijita M., Okuno, A., 2005. Development of sodium-sulfur batteries. International Journal of Applied Ceramic Technology, 1, 269–276. Risch T., Newman J., 1988. Transference number calculations for sodium polysulfides. Journal of The Electrochemical Society, 135, 1715. Schaefer S., Vudata S.P., Bhattacharyya D., Turton R., 2020. Transient modeling and simulation of a nonisothermal sodium–sulfur cell. Journal of Power Sources, 453, 227849. Thompson S.D., Newman, J., 1989. Differential diffusion coefficients of sodium polysulfide melts. Journal of The Electrochemical Society, 136, 3362. 1098 PRES22_0376.pdf Multiphysics Modeling of High Temperature Planar Sodium Sulfur Batteries