Electrochemical-Thermal Modelling and Optimisation of Lithium-Ion Battery Design Parameters Using Analysis of Variance

A 1D electrochemical-thermal model of an electrode pair of a lithium ion battery is developed in Comsol Multiphysics. The mathematical model is validated against the literature data for a 10 Ah lithium phosphate (LFP) pouch cell operating under 1 C to 5 C electrical load at 25 ◦C ambient temperature. The validated model is used to conduct statistical analysis of the most influential parameters that dictate cell performance; i.e., particle radius (rp); electrode thickness (Lpos); volume fraction of the active material (εs,pos) and C-rate; and their interaction on the two main responses; namely; specific energy and specific power. To achieve an optimised window for energy and power within the defined range of design variables; the range of variation of the variables is determined based on literature data and includes: rp: 30–100 nm; Lpos: 20–100 μm; εs,pos: 0.3–0.7; C-rate: 1–5. By investigating the main effect and the interaction effect of the design variables on energy and power; it is observed that the optimum energy can be achieved when (rp < 40 nm); (75 μm < Lpos < 100 μm); (0.4 < εs,pos < 0.6) and while the C-rate is below 4C. Conversely; the optimum power is achieved for a thin electrode (Lpos < 30 μm); with high porosity and high C-rate (5 C).


Introduction
The creation of battery models is essential for better understanding of the impact of design variables as well as operating conditions on battery performance, in terms of efficiency, degradation and safety.They are also cost effective tools for determination of an optimised battery design that can reduce the need for experimental trial and error.Beyond that, most variables that underpin the operation of a battery are not directly measurable, whereas they can be predicted through a validated model [1].
Battery models can be broadly categorised into four groups, namely: empirical models (equivalent circuit and neural network), electrochemical engineering models, multiphysics models, and molecular/atomistic models [1,2].These battery models include different levels of details and differ in terms of complexity, computational cost and reliability.They can be chosen based on the particular needs for a specific application.Due to the complexity of batteries, electrochemical models are often viewed as the best approach for investigating the impact of design variables on battery performance during a charge and discharge process.
The most common electrochemical models in the literature [3] are single particle models (SPM) [4,5] porous-electrode models [6,7] and pseudo-two-dimensional models (P2D) [8,9].The main difference is the level of complexity and the computational time associated with their use.Once an Energies 2017, 10, 1278 2 of 22 efficient model is developed, it can be employed to address a number of real-world challenges, such as the identification of transport and kinetic parameters, the occurrence of capacity fade, improving life time, and improving energy/power density [1].The energy/power density can be improved by manipulating either the design parameters or operating protocols [1].Much of the previous work on improving the performance of batteries has focused on battery packs rather than a single battery cell design [10].The important design parameters which can be varied within the manufacturing process to achieve the optimal battery performance are known to be: electrode thickness [11][12][13][14], porosity [6,[10][11][12]15], particle size [14,16], electrode surface area, geometry and the dimensions of current collectors [6,17,18].
Given the importance of battery design, Newman et al. and co-workers [11,12] developed an analytical model of a lithium ion battery to optimise porosity and thickness of the positive electrode for maximum specific energy, while holding other parameters constant.In another study, Singh et al. [13] investigated experimentally the amount of energy that may be extracted from a cell manufactured using thick electrodes (320 µm) compared to cells that employ a thinner electrode (70 µm) for a graphite/lithium nickel manganese cobalt oxide (Gr/NMC) chemistry.They observed a significant capacity loss when using the thicker cells at C-rates of C/2 due to poor kinetics.The authors suggest that the proposed thick electrodes could be advantageous for certain applications where a continuous low C-rate is required.Research presented by Wu et al. [14], employed an electrochemical-thermal model of a lithium ion battery for a Li y Mn 2 O 4 chemistry.The model was used to investigate the impact of particle size and electrode thickness on heat generation and the performance of the battery.They found that a battery containing a thin electrode showed a better performance in terms of temperature rise and material utilisation, but the effect of particle size was not monotonic across the discharge rates.In a recent study published by Ramadesigan et al. [6] a multi-layered porosity distribution was investigated rather than optimising a uniform porosity in a lithium ion battery design.They developed a simple electrochemical porous electrode model, which did not include the solid-phase intercalation mechanism.The model was applied for a cathode made of lithium cobalt oxide.For a fixed value of active material, optimal multi-layered porosity distribution across the positive electrode was found.The authors managed to decrease the ohmic resistance by circa: 15-33%, employing the optimal porosity distribution.Golmon et al. [15] applied a gradient-based optimisation method on a multi-scale battery model to maximise the usable capacity.The electrochemical-mechanical multi-scale model was an extension of Doyle and Newman's electrochemical battery model [19,20].The variables were particle size and porosity, with constraints placed on the stress levels in the cathode particles.Other examples can be seen in [16,18].Darling et al. [16], developed a 1D theoretical model to investigate the influence of different particle size distributions on the operation of porous intercalation electrodes.Chen et al. [18] developed a technique to enhance the achieved capacity of Li-based batteries through improvements in both electronic and ionic conductivity of materials.
Parameter sensitivity study is another area of research performed in modelling approaches for finding the influential design variables which could potentially improve battery performance [14,[21][22][23][24][25].Zhang et al. [21,23] developed a coupled P2D electrochemical-thermal model for a 2.3 Ah cylindrical Li-ion battery by improving the open source FORTRAN code maintained by Newman's research group.They established a parameter sensitivity matrix and used clustering theory to group the parameters according to their average sensitivity.They performed a sensitivity study of up to 30 parameters under different operating conditions.Out of the 30 parameters, 10 were found to be highly sensitive, seven of those were sensitive, along with 10 low sensitivity and three insensitive parameters.The most sensitive parameters were found to be anode particle radius, diffusion coefficient in the negative electrode, stoichiometry of the anode, volume fraction of the electrolyte and active material in the negative electrode, contact resistance, reaction rate of the negative and positive electrode as well as the activation energy of the electrolyte ionic conductivity.Moreover, it was found that the sensitivity of the parameters is strongly dependent on the operating conditions.Edouard et al. [22] developed a single particle (SP) electrochemical-thermal model applying the pseudo-2D mathematical structure.The model was applied to evaluate the sensitivity of the key parameters that are involved in battery aging.They pointed out that the limiting mechanism in a battery can switch depending on the operating conditions and battery design [14].Moreover, the combined effect of design parameters can be even more significant than their individual effect.This presents a great challenge for experimentally optimising a battery design.Du et al. [24] introduced a surrogate modelling framework to map the effect of design parameters, such as cathode particle size, diffusion coefficient and electrical conductivity on battery performance.They quantified the relative impact of various parameters through global sensitivity analysis employing a cell-level model in conjunction with tools such as kriging, polynomial response, and radial-basis neural networks.Ghaznavi et al. [25] applied a mathematical approach to conduct a sensitivity study on a lithium-sulfur cell.They focused on the effects of discharge current and conductivity of the positive electrode over a wide range of values.
To date, much progress has been made in modelling and design optimisation of lithium ion batteries to map the trade-off between power and energy density.However, in most studies only a few design variables have been selected and the optimisation is limited to a narrow range of operating conditions.Moreover, the interaction effect of influential parameters has not been comprehensively studied.Therefore, there exists a critical need to stablish a framework to access both the individual and the interaction effect of various parameters on the energy and power of a cell.This study attempts to quantify the strength of design factors and the combined effect of variables on specific energy as well as specific power of a battery.
Within this paper, a 1D electrochemical-thermal model of an electrode pair of a lithium ion battery is developed in Comsol Multiphysics.Each pair is assumed to be a sandwiched model of different layers, a negative current collector, a negative electrode, a separator, a positive electrode and a positive current collector.The anode is made of graphite and the cathode material is lithium phosphate (LFP).The mathematical model is validated against the literature data for a 10 Ah LFP pouch cell operating under 1 C to 5 C electrical load at 25 • C ambient temperature.It is a commercial cell and can be employed for vehicle applications as it is a large format cell.The validated model is used to conduct statistical analysis of the most influential parameters that dictate cell performance, i.e., particle size (r p ), electrode thickness (L pos ), volume fraction of the active material (ε s,pos ) and C-rate, and their interaction on the two main responses, namely; specific energy and specific power.This is to achieve an optimised window for energy and power within the defined range of design variables.The design factors are chosen in a way that they can be varied during the manufacturing process of a cell, in order to make the developed statistical model more applicable for industry.
In Section 2, the mathematical modelling approach is explained and contains the derivation of the electrochemical-thermal model along with the statistical analysis.Section 3 presents the model validation accompanied by simulation results of analysis of variance (ANOVA), which elaborates the main and combined effect of the factors as well as the optimised design of the cell.Further work and conclusions are discussed in Sections 4 and 5, respectively.

Electrochemical-Thermal Model
A single cell is made up of a number of pairs, connected in parallel by current collectors of both electrodes, and the collective behaviour of these pairs represents the overall cell performance.Each pair is assumed to be a sandwiched model of different layers, a negative current collector, a negative electrode, a separator, a positive electrode and a positive current collector.The electrodes are composed of porous materials, which are filled with electrolyte.The electrolyte can be liquid, solid or polymer.The porous separator is a membrane between the negative and the positive electrodes, which allows the transport of lithium ions from one electrode to the other one.Therefore it needs to be Energies 2017, 10, 1278 4 of 22 ionically conductive and electrically insulating to prevent the cell from short circuiting [26,27].Current collectors are made from conductive materials.Typical materials for negative and positive current collectors are copper and aluminium, respectively.During a discharge event of a lithium ion battery the lithium is deintercalated from the anode, resulting in an equal number of electrons and ions at the particles/electrolyte interface.The chemical reaction only takes place at the surface of electrode particles that are in contact with the electrolyte.Lithium ions pass through the electrolyte and the separator and migrate to the cathode side where they react with the positive electrode.Diffusion of ions is driven by the lithium ion concentration.Electrons pass through the current collectors and the external circuit and generate electricity.
In this study, a 1D electrochemical-thermal model is developed for a single electrode-pair of a lithium ion battery with a LFP cathode.The model is based on a similar, experimentally validated model presented within Li et al. [28].The inputs to the model are load current, geometrical design parameters, material properties and ambient operating temperature, while the outputs from the model are the responses of the cell to the load, i.e., terminal voltage (V), generated heat (Q), temperature profile ( • C) and state of charge (SOC).Additional internal variables include: lithium concentration in the electrodes and the separator, potential distribution of different phases, reaction current, electronic and ionic current.The primary motivation for using Li et al.'s model [28] for verification, is because the authors have presented a complete set of electrochemical parameters for the cell along with a definition of parameter variations with temperature.Since the electrochemical model contains temperature dependent variables, coupling it with a thermal model yields more accurate results.Furthermore, the model presented in [28] has been validated using experimental data.The schematic geometry of a single electrode-pair is shown in Figure 1.
Energies 2017, 10, 1278 4 of 21 battery the lithium is deintercalated from the anode, resulting in an equal number of electrons and ions at the particles/electrolyte interface.The chemical reaction only takes place at the surface of electrode particles that are in contact with the electrolyte.Lithium ions pass through the electrolyte and the separator and migrate to the cathode side where they react with the positive electrode.Diffusion of ions is driven by the lithium ion concentration.Electrons pass through the current collectors and the external circuit and generate electricity.
In this study, a 1D electrochemical-thermal model is developed for a single electrode-pair of a lithium ion battery with a LFP cathode.The model is based on a similar, experimentally validated model presented within Li et al. [28].The inputs to the model are load current, geometrical design parameters, material properties and ambient operating temperature, while the outputs from the model are the responses of the cell to the load, i.e., terminal voltage (V), generated heat (Q), temperature profile (°C) and state of charge (SOC).Additional internal variables include: lithium concentration in the electrodes and the separator, potential distribution of different phases, reaction current, electronic and ionic current.The primary motivation for using Li et al.'s model [28] for verification, is because the authors have presented a complete set of electrochemical parameters for the cell along with a definition of parameter variations with temperature.Since the electrochemical model contains temperature dependent variables, coupling it with a thermal model yields more accurate results.Furthermore, the model presented in [28] has been validated using experimental data.The schematic geometry of a single electrode-pair is shown in Figure 1.The current flowing through a single pair is calculated as follows: where is the capacity of a battery cell (including all pairs) and represents the number of electrode-pairs in the cell.The dynamic performance of the cell is characterized by the solution of four partial differential equations describing the time evolution of the lithium concentration profile in the electrode and electrolyte phases.The field variables to be calculated are , , , [26,28].and are lithium concentration in the solid and the electrolyte phase, and represents the potential in the solid and electrolyte phase respectively.The variables are evaluated through solving electric charge conservation, mass conservation, electrochemical kinetics and energy conservation equations [26,[29][30][31][32][33]].

Mass Balance
Lithium in the Solid Phase Identical particle size of active materials is assumed in the model.The distribution of lithium in the electrode is described by Fick's second law: The current flowing through a single pair is calculated as follows: where δ cell is the capacity of a battery cell (including all pairs) and N pairs represents the number of electrode-pairs in the cell.The dynamic performance of the cell is characterized by the solution of four partial differential equations describing the time evolution of the lithium concentration profile in the electrode and electrolyte phases.The field variables to be calculated are C s , C e , φ s , φ e [26,28].C s and C e are lithium concentration in the solid and the electrolyte phase, φ s and φ e represents the potential in the solid and electrolyte phase respectively.The variables are evaluated through solving electric charge conservation, mass conservation, electrochemical kinetics and energy conservation equations [26,[29][30][31][32][33].

Mass Balance
Lithium in the Solid Phase Identical particle size of active materials is assumed in the model.The distribution of lithium in the electrode is described by Fick's second law: where C s , is the lithium ion concentration at the surface of the electrode and D s is the lithium diffusion coefficient in the solid phase.No species source exists at the centre of the electrode particles, hence the boundary condition is defined as follows: where r p is the radius of the active material, a s is the reaction surface area and j Li represents the reaction current density.

Lithium in the Electrolyte Phase
The lithium ion distribution in the electrolyte is dependent on the lithium diffusion in the electrolyte (C e ), the electrode porosity (ε e ), and the reaction current density (j Li ).This phenomenon is expressed by Fick's second law: where F is Faraday's constant and t 0 + represents the initial transference number.Lithium cannot diffuse through the current collectors, as set by the following boundary condition:

Electronic Charge Balance
During a charge and discharge event lithium ions and electrons flow in reverse directions, charge conservation states that the amount of lithium ions has to be equivalent to the amount of electronic charge transfer.

Potential in the Solid Phase
The potential in the solid phase is a function of electrode conductivity, current collector conductivity and reaction current density as presented by the following equation: At the electrode/current collector interface, the charge flux represents the external current, as expressed by: Also, there is no charge flux at the electrode/separator interface as stated below: Energies 2017, 10, 1278 Potential in the Electrolyte Phase The potential in the electrolyte phase is dependent on the reaction current density as well as the local concentration of lithium as expressed in the following equation: where k e f f D is the effective diffusional conductivity of the species.A zero gradient boundary condition is imposed at the electrode/current collector interfaces:

Electrochemical Kinetics-Reaction Current Density
The lithium concentration and charge distribution in the electrode and electrolyte phases, which are described through Equations ( 1)- (10), are coupled through the Butler-Volmer equation: The reaction surface area (a s ), is the interfacial area between the two phases, i.e., the solid active material and the liquid electrolyte.The interfacial area can be calculated by treating the solid phase as a collection of uniform spheres, as displayed by the following equation: where N p is the number of spheres per unit volume.The volume fraction of the active material (ε s ) is given by: ε s = N p 4/3πr p 3 (13) By combining the two equation the reaction surface area is evaluated as follows [34]: where ε f is the volume fraction of the fillers.The exchange current density (i 0 ) is a function of temperature and SOC and is expressed as: k i is the reaction rate which is temperature dependent.α is a dimensionless parameter, called the symmetry factor and defines the ratio between oxidation and reduction.The overpotential (η), is defined as the difference between the open circuit voltage (OCV) and the operating voltage of the positive/negative electrode.
The potential of the positive and negative electrodes are a function of local lithium concentration and are defined by using empirical expressions.The parameters and equations employed in this study for developing the electrochemcial-thermal model are elaborated in Appendix A.

Statistical Analysis of Design Variables
A statistical model is developed by ANOVA of the numerical data in a full factorial design framework [35][36][37][38][39]. ANOVA is primarily designed for the analysis of experimental data.However it is noteworthy that a number of comparable studies have also employed ANOVA on the numerical results from verified simulation models, for example: [24,[40][41][42][43][44][45].A full factorial design methodology is undertaken to analyse the results of the 1D electrochemical-thermal model and to determine the optimum energy and power by manipulating key design variables of the positive electrode.As highlighted in recent studies, cells that employ a LFP cathode material, the operation of the cell is more sensitive to the design of the cathode rather that the anode.This is due to the poor electronic conductivity and low solid diffusion coefficient of the LFP material [46].For this reason, this study focused on understanding the optimisation of the cathode rather than the anode.It is noteworthy that the methodology employed is transferable to the anode.

Full Factorial Design
The most influential parameters that dictate cell performance are known to be particle size (r p ), electrode thickness (L pos ), volume fraction of the active material (ε s,pos ) and C-rate.The effects of these four factors and their interaction on the two main responses, namely; specific energy and specific power are investigated.The design factors are chosen in a way that they can be varied during the manufacturing process of a cell, in order to make the developed statistical model more applicable for industry.In a 3-level full factorial design, each factor is varied over three levels within the determined ranges, as summarised in Table 1.The range of variation of the design variables for LFP lithium ion battery is determined based on the data from literature [11,28,29,33,[47][48][49], which is: r p : 30-100 nm, L pos : 20-100 µm, ε s,pos : 0.3-0.7,C-rate: 1-5.The number of required simulations in a 3 level full factorial design for four design variables is equal to 3 4 , i.e., 81 numerical case studies.It is clear that investigating 81 test cases, having different designs is not feasible through traditional physical design and experimentation.A simulation based approach highlights one of the advantages of ANOVA in terms of cost and time.As discussed within [48], for the experimentalist it is very challenging to optimise the cell when facing with numerous variables.The model can assist in identifying limiting cell properties and once identified, the experiments can be designed.

Analysis of Variance
The contribution of each factor to the response (energy/power) is evaluated by the ANOVA, based on which the response is expressed as follows: A semi-empirical regression model is employed for analysing the responses, as given by [35,36]: Energies 2017, 10, 1278 where y is the predicted response, the β's are the coefficients to be determined, the x values are the factors and defines a random error term.The quality of the fitness of the regression model is defined by the value of the correlation coefficient (R 2 ) value as well as the normal probability plot of the residuals.For an adequate model R 2 should be close to 1 [50].In general for a well fitted model R 2 value should be above 0.9 [35].Other important outputs from ANOVA are F-values (Fisher variation ratio) and p-values (probability value) [36].p-Values define whether the effect of a term is significant or not.The calculations are based on the confidence level of 95%, meaning that if the (p-value < 0.05), the factor is significant.In addition, the order of significance for the influential factors is found by comparing the F-values.A higher F-value indicates that the factor is more significant.

Validation of the Coupled Electrochemical-Thermal Model
A 1D electrochemical-thermal model is developed for a single electrode-pair of a lithium ion battery with LFP cathode.The model is based on a similar, experimentally validated model presented within Li et al. [28].The input to the model is load current, geometrical design parameters, material properties and ambient operating temperature.The parameterisation data and experimental results presented within [28] allow for the creation of a reference model to ensure the accuracy and robustness of the model, before proceeding with the statistical analysis.The design parameters of the reference case are those reported by [28] and are summarised in Table 2.
The output from the model shows a high degree correlation with the results published within [28] in terms of the terminal voltage response when the cell is discharged, the surface temperature gradient and other internal parameters (e.g., local SOC, lithium concentration, potential and internal current).The discharge curve for the 10 Ah LFP pouch cell under 1 C to 5 C at 25 • C is presented in Figure 2. It is seen that there is a good agreement between the simulation results and the experimental data.The simulation results are generated by the 1D COMSOL model used in this study and the experimental data are those reported by Li et al. [28].presented within [28] allow for the creation of a reference model to ensure the accuracy and robustness of the model, before proceeding with the statistical analysis.The design parameters of the reference case are those reported by [28] and are summarised in Table 2.
The output from the model shows a high degree correlation with the results published within [28] in terms of the terminal voltage response when the cell is discharged, the surface temperature gradient and other internal parameters (e.g., local SOC, lithium concentration, potential and internal current).The discharge curve for the 10 Ah LFP pouch cell under 1 C to 5 C at 25 °C is presented in Figure 2. It is seen that there is a good agreement between the simulation results and the experimental data.The simulation results are generated by the 1D COMSOL model used in this study and the experimental data are those reported by Li et al. [28].The validated model is used as a tool for simulating the 81 case studies (defined within Table 1) which are required for conducting the statistical analysis.For each case study the specific energy ( ) and specific power ( ) of the cell are defined as [51]: The validated model is used as a tool for simulating the 81 case studies (defined within Table 1) which are required for conducting the statistical analysis.For each case study the specific energy (E cell ) and specific power (P cell ) of the cell are defined as [51]: where m cell is the total weight of the cell and t dis is the final time when the cell reaches the cut off voltage of 2.5 V.

Results of ANOVA
The results of the fitted model derived from ANOVA for the selected responses, i.e., specific energy and specific power are shown in Table 3.In this study the R 2 values for the specific energy and power are equal to 96.04% and 100% respectively, verifying the high accuracy of the model.
As mentioned earlier, another important issue to check for the adequacy of the model is the normal probability plot of the residuals, as illustrated in Figure 3 for the specific energy and the specific power.The residuals follow a straight line, revealing that the fitted regression models are valid [50,52].
The p-values provide a cut-off beyond which it is asserted that the terms are statistically significant.Based on the displayed results for R 1 in Table 3, except for (r p × L pos ), for which the p-value is higher than 0.05, the rest of the factors are significant and the achieved specific energy of the cell is highly influenced by those factors.The ranking of the effective factors in the regression model, within the determined range, based on their impact on specific energy, R 1 is as follows: In respect to the specific power, based on the presented results for R 2 in Table 3, all the factors except for (r p × C) are significant and the order of significance is as follows: In respect to the specific power, based on the presented results for R2 in Table 3, all the factors except for ( × C) are significant and the order of significance is as follows:

Specific Energy: The Main and Interaction Effects of Factors
The main effect plot reveals the sole effect of the factors individually on a response, in which the mean value of the response is displayed for each factor level. Figure 4 shows the main effect plot of the factors on the specific energy.It is seen that as increases from 30 nm to 100 nm the specific energy decreases significantly.This can be attributed to the diffusion time ( ), which is defined as: where , is the lithium ion diffusion coefficient in the positive electrode [24].Using smaller particles can reduce the diffusion time, which in turn increases the intercalation rate.Consequently,

Specific Energy: The Main and Interaction Effects of Factors
The main effect plot reveals the sole effect of the factors individually on a response, in which the mean value of the response is displayed for each factor level. Figure 4 shows the main effect plot of the factors on the specific energy.It is seen that as r p increases from 30 nm to 100 nm the specific energy decreases significantly.This can be attributed to the diffusion time (t di f ), which is defined as: where D s,pos is the lithium ion diffusion coefficient in the positive electrode [24].Using smaller particles can reduce the diffusion time, which in turn increases the intercalation rate.Consequently, it can enhance the rate of the electrochemical reaction.On the other hand smaller particles have a larger surface area per unit volume which leads to a lower kinetic resistance [53].
Energies 2017, 10, 1278 10 of 21 it can enhance the rate of the electrochemical reaction.On the other hand smaller particles have a larger surface area per unit volume which leads to a lower kinetic resistance [53].The proper cell design requires a right balancing of the electrodes and it is known that the capacity of the electrodes is proportional to the weight of the active materials which itself is function of electrode thickness and volume fraction of the active materials.Hence, in this study the thickness of the negative electrode was varied proportional to the film thickness of the positive electrode to keep the right balancing.By increasing the thickness, the energy keeps increasing as more active material is available.As shown in Figure 4, by increasing from 20 μm to 60 μm the energy increases quite rapidly, but from 60 μm to 100 μm the rate of energy improvement reduces significantly which can be attributed to the higher transport resistance of the thick electrodes.
By increasing , from 0.3 to 0.7, the energy increases first until it reaches to the point at which the capacity of the positive and the negative electrodes are identical.After that, further increasing of The proper cell design requires a right balancing of the electrodes and it is known that the capacity of the electrodes is proportional to the weight of the active materials which itself is function of electrode thickness and volume fraction of the active materials.Hence, in this study the thickness of the negative electrode was varied proportional to the film thickness of the positive electrode to keep the right balancing.By increasing the thickness, the energy keeps increasing as more active material is available.As shown in Figure 4, by increasing L pos from 20 µm to 60 µm the energy increases quite rapidly, but from 60 µm to 100 µm the rate of energy improvement reduces significantly which can be attributed to the higher transport resistance of the thick electrodes.By increasing ε s,pos from 0.3 to 0.7, the energy increases first until it reaches to the point at which the capacity of the positive and the negative electrodes are identical.After that, further increasing of the active volume fraction of the positive electrode cannot improve the energy any more as the capacity of the negative electrode will be the limiting factor.Higher C-rates replicate a higher current, and it is known that the rate of the electrochemical reaction is higher under this condition.Hence, due to limitations in the ionic and electrical conductivities, higher kinetic resistances is observed, which in turn causes a significant energy loss [54].This is in agreement with the results in presented Figure 4, which illustrates by increasing the C-rate from 1 C to 5 C, the energy drops significantly.
One of the key findings from ANOVA is, clarifying the interplay between the different factors, and hence establishing a process window for a desired range of a certain predefined response, specific energy and power in this study.Figure 5 illustrates the second order interactions of any two factors on the specific energy.Based on the ANOVA results in Table 3, among the second order interactions, ( × , ) was found to be the most significant factor, and it is clearly displayed in Figure 5. From the main effect plot, Figure 4, has a positive effect on the specific energy.However, ( × , ) plot reveals that at = 0.3-0.5, by increasing , higher energy can be achieved, whereas the energy does not follow the same pattern at higher values.For ( > 0.5), increasing both and at the same time leads to an energy drop and the reduction rate is higher for higher active volume fraction.
The next important factor is the interaction of ( × ).As seen has a negative impact while does not show a monotonic effect on the energy.The interaction plot shows that by increasing the particle size (from 30 to 100 nm) at low , the energy drops rapidly.However, the sensitivity of the energy to the particle size is less pronounced at higher , , where a higher portion of active material is available and it will compensate the energy loss.and rate both have a negative effect on the energy.The same trend is displayed in their interaction plot, ( × ).Moreover, the derived results show at higher C-rates the energy is more sensitive to the particle size.Meaning that for high C-rate application small particle are more desirable while at lower C-rates there is a bigger window for the optimum .Lower increases the intercalation rate, which in turn it takes a shorter time for the lithium ions to diffuse through the solid particles.This is particularly advantageous for high C-rate application where the discharge time is quite short.Based on the ANOVA results presented in Table 3, among the second order interactions, (L pos × ε s,pos ) was found to be the most significant factor, and it is clearly displayed in Figure 5. From the main effect plot, Figure 4, L pos has a positive effect on the specific energy.However, (L pos × ε s,pos ) plot reveals that at ε pos = 0.3-0.5, by increasing L pos , higher energy can be achieved, whereas the energy does not follow the same pattern at higher ε pos values.For (ε pos > 0.5), increasing both ε pos and L pos at the same time leads to an energy drop and the reduction rate is higher for higher active volume fraction.

Specific Power: The Main and Interaction Effects of Factors
The next important factor is the interaction of (r p × ε pos ).As seen r p has a negative impact while L pos does not show a monotonic effect on the energy.The interaction plot shows that by increasing the particle size (from 30 to 100 nm) at low ε s,pos the energy drops rapidly.However, the sensitivity of the energy to the particle size is less pronounced at higher ε s,pos , where a higher portion of active material is available and it will compensate the energy loss.r p and C rate both have a negative effect on the energy.The same trend is displayed in their interaction plot, (r p × C).Moreover, the derived results show that at higher C-rates the energy is more sensitive to the particle size.Meaning that for high C-rate application small particle are more desirable while at lower C-rates there is a bigger window for the optimum r p .Lower r p increases the intercalation rate, which in turn it takes a shorter time for the lithium ions to diffuse through the solid particles.This is particularly advantageous for high C-rate application where the discharge time is quite short.

Specific Power: The Main and Interaction Effects of Factors
Figure 6 presents the main effect plot of the design parameters on the specific power.It is shown that the particle size of the positive electrode does not have a significant effect on the power.Moreover, the power reduces for thicker electrodes, which is due to longer transport length of the ions and electrons.In addition to that, the electrode volume fraction does not have a significant effect on the power, however by increasing ε s,pos from 0.3 to 0.7 the power is slightly reduced.That is because of the lower ionic conductivity of the positive electrode at high ε s,pos .Finally, it is clearly observed that, the power linearly increases by the C-rate, which is in agreement with the Ohm's law.The same kind of combined influence as for energy, is shown for the specific power (see Figure 7).has a negative effect on the power, whereas C-rate has a positive effect.Their interaction, ( × ), shows that the thickness of the electrode becomes a more influential factor as the C-rate increases, whereas at 1C the thickness does not have a significant impact on the power.The behaviour of the power to the interaction of and C-rate, ( × ), is not monotonic.
For example by varying from 0.3 to 0.5, the power does not change irrespective of the C-rate.However, their combined effect is more pronounced at high active volume fraction when the C-rate is high.By increasing from 0.5 to 0.7 under 5 C discharge load, the power drops with a higher rate in comparison to lower C-rates.The same kind of combined influence as for energy, is shown for the specific power (see Figure 7).L pos has a negative effect on the power, whereas C-rate has a positive effect.Their interaction, (L pos × C), shows that the thickness of the electrode becomes a more influential factor as the C-rate increases, whereas at 1C the thickness does not have a significant impact on the power.
The behaviour of the power to the interaction of ε pos and C-rate, (ε pos × C), is not monotonic.For example by varying ε pos from 0.3 to 0.5, the power does not change irrespective of the C-rate.However, their combined effect is more pronounced at high active volume fraction when the C-rate is high.By increasing ε pos from 0.5 to 0.7 under 5 C discharge load, the power drops with a higher rate in comparison to lower C-rates.
increases, whereas at 1C the thickness does not have a significant impact on the power.
The behaviour of the power to the interaction of and C-rate, ( × ), is not monotonic.
For example by varying from 0.3 to 0.5, the power does not change irrespective of the C-rate.However, their combined effect is more pronounced at high active volume fraction when the C-rate is high.By increasing from 0.5 to 0.7 under 5 C discharge load, the power drops with a higher rate in comparison to lower C-rates.), volume fraction of the active material ( , ) and C-rate) at different levels on the specific power.

Optimised Battery Response
The optimised response is derived from the developed regression model through statistical analysis with much less simulation effort [45].Such results help to understand the existing interactions between the variables and to achieve the optimum design for a particular application.Figure 8 presents the contour plot of the specific energy and the specific power as a function of electrode thickness and particle size (Figure 8a,b), volume fraction of active material and particle size (Figure 8c,d), as well as C-rate and particle size (Figure 8e,f).

Optimised Battery Response
The optimised response is derived from the developed regression model through statistical analysis with much less simulation effort [45].Such results help to understand the existing interactions between the variables and to achieve the optimum design for a particular application.Figure 8 presents the contour plot of the specific energy and the specific power as a function of electrode thickness and particle size (Figure 8a,b), volume fraction of active material and particle size (Figure 8c,d), as well as C-rate and particle size (Figure 8e,f).
Discrete regions represent different energy or power levels as stated in the plots.The plots highlight the impact of two factors on the response simultaneously.By analysing the combined effect of different factors on the response the optimum range can be achieved.Having r p and L pos as variables (Figure 8a,b), the maximum energy is attained for 75 < L pos < 100 and r p < 40 nm.For r p > 40 nm the energy constantly reduces.From the power contour, it is observed that as the thickness of the electrode increases the power decreases.Hence, it is important to select a range of design variables in which the energy to power ratio is satisfactory for a specific application.
Besides r p and L pos , ε pos has also a significant effect on the energy.As displayed for porous electrodes, the particle size should be very small, below 40 nm, to achieve the optimum energy and yet the optimum region is quite restricted.On the other hand, from the power graph it is observed that, the power is not very sensitive to the interaction of these factors (r p and ε pos ).On the other hand, from the power graph it is observed that, the power is not very sensitive to the interaction of these factors (r p and ε pos ), and irrespective of the particle size, a higher power is obtained as the porosity increases.C-rate is another important factor to be considered when designing a battery.In general, the design of batteries operating at high C-rates is more crucial for achieving high power and energy.For example, at 1 C-rate, by changing the particle size from 30 nm to 100 nm, the energy drops by 21.6%, whereas the energy reduction reaches 55% when operating under a 5C discharge.In contrast to the energy, the power is not highly affected by the particle size.However, C-rate has a greater impact on it and as the C-rate increases the power is monotonically increasing.
The combined effect of (L pos , ε s,pos ) as well as (L pos , C) on the batteries responses is presented in Figure 9.The energy plot based on (L pos and ε s,pos ) shows a circle for the optimum energy where 70 < L pos < 100 and 0.4 < ε s,pos < 0.62, and outside this area the energy decreases.Moreover, from the power contour it is observed that for low porous electrode, 0.6 < ε s,pos < 0.7, a lower power is obtained.In general, for having a reasonable power at high ε s,pos , the electrode should be quite thin.Otherwise the power reduction is substantially high, especially at higher C-rates.
Energies 2017, 10, 1278 13 of 21 Discrete regions represent different energy or power levels as stated in the plots.The plots highlight the impact of two factors on the response simultaneously.By analysing the combined effect of different factors on the response the optimum range can be achieved.Having and as variables (Figure 8a,b), the maximum energy is attained for 75 < < 100 and < 40 nm.For > 40 nm the energy constantly reduces.From the power contour, it is observed that as the thickness of the electrode increases the power decreases.Hence, it is important to select a range of design variables in which the energy to power ratio is satisfactory for a specific application.Besides and , has also a significant effect on the energy.As displayed for porous electrodes, the particle size should be very small, below 40 nm, to achieve the optimum energy and yet the optimum region is quite restricted.On the other hand, from the power graph it is observed that, the power is not very sensitive to the interaction of these factors ( and ).On the other hand, from the power graph it is observed that, the power is not very sensitive to the interaction of these factors ( and ), and irrespective of the particle size, a higher power is obtained as the porosity increases.C-rate is another important factor to be considered when designing a battery.In general, the design of batteries operating at high C-rates is more crucial for achieving high power The maximum energy is achieved for C-rates less than 4. It is observed that as the C-rate increases the energy decreases and the optimum range for the energy becomes restricted.At 5C, to achieve the maximum energy, L pos has to be higher than 60 µm.In contrast to the energy, the power is higher at higher C-rates.However, it is very sensitive to the electrode thickness.For example the highest power can be achieved for L pos < 20 µm and as the thickness increases the cell will end up to a lower power region.
Energies 2017, 10, 1278 15 of 22 , in Figure 9.The energy plot based on ( and , ) shows a circle for the optimum energy where 70 < < 100 and 0.4 < , < 0.62, and outside this area the energy decreases.Moreover, from the power contour it is observed that for low porous electrode, 0.6 < , < 0.7, a lower power is obtained.In general, for having a reasonable power at high , , the electrode should be quite thin.Otherwise the power reduction is substantially high, especially at higher C-rates.The maximum energy is achieved for C-rates less than 4. It is observed that as the C-rate increases the energy decreases and the optimum range for the energy becomes restricted.At 5C, to achieve the maximum energy, has to be higher than 60 μm.In contrast to the energy, the power is higher at higher C-rates.However, it is very sensitive to the electrode thickness.For example the highest power can be achieved for < 20 μm and as the thickness increases the cell will end up to a lower power region.

Further Work
In this paper, we employed a 3-level full factorial design on the numerical results obtained from an electrochemical-thermal model.The influence of four design variables were investigated to achieve the optimal responses (here specific energy and power).However, the developed model is not limited to the defined design variables and responses.By changing the factors, or adding new design variables to the existing model, a new set of simulation studies can be run and used for further analysis.As an example, it would be of interest to investigate whether the optimal design for achieving the highest energy and power is influenced by the operating temperature of the cell.

Further Work
In this paper, we employed a 3-level full factorial design on the numerical results obtained from an electrochemical-thermal model.The influence of four design variables were investigated to achieve the optimal responses (here specific energy and power).However, the developed model is not limited to the defined design variables and responses.By changing the factors, or adding new design variables to the existing model, a new set of simulation studies can be run and used for further analysis.As an example, it would be of interest to investigate whether the optimal design for achieving the highest energy and power is influenced by the operating temperature of the cell.Another area of interest is to investigate the impact of particle size as well as porosity distribution on the obtained energy and power.
Moreover, battery degradation rate is another interesting system response to be considered.However, studying the degradation process, needs improvements to the current model by adding the impact of solid-electrolyte interface (SEI) growth and its effect on the capacity fade within the battery.
Finally, as mentioned earlier in this study, the optimised design of the cell was achieved through analysing 81 test cases.Given this number, it seems quite unlikely to get all the required data through experiments.On the other hand, to make a comparable test case, all of the materials (including the separator and electrolyte) must be similar to those of the commercial cell.Furthermore, the production procedure for the cell itself has to be similar.Understanding such detailed design information for a commercial cell is known to be difficult because of issues of confidentiality.However, manufacturing a new proto-type cell (manufactured by the University) with optimised design accompanied by an experimental evaluation is another step to be considered to ultimately validate the proposed optimised cell design and simulation framework.

Conclusions
In this study a 1D electrochemical-thermal model of an electrode pair of a lithium ion battery is developed in Comsol Multiphysics.The mathematical model is validated against the literature data for a 10 Ah LFP pouch cell operating under 1 C to 5 C electrical load at 25 • C ambient temperature.The validated model is used to conduct statistical analysis of the most influential parameters that dictate cell performance, i.e., particle size (r p ), electrode thickness (L pos ), volume fraction of the active material (ε s,pos ) and C-rate, and their interaction on the two main responses, namely; specific energy and specific power.This is to achieve an optimised window for energy and power within the defined range of design variables.The range of variation of the design variables for LFP lithium ion battery is determined based on data from literature (r p : 30-100 nm, L pos : 20-100 µm, ε s,pos : 0.3-0.7,C-rate: 1-5).A statistical model is developed by ANOVA of the numerical data in a full factorial design frame work.A full factorial design methodology is carried out to analyse the obtained results of the 1D electrochemical-thermal model and to determine the optimum energy and power by manipulating key design variables of the positive electrode.The summary of the statistical results are as follows: The significant factors for the specific energy are ranked as: Similarly, for the specific power it is defined as: In conclusion, the main effect and the interaction effect of all design variables on the energy and power, it is observed that the optimum energy can be achieved when (r p < 40 nm), (75 µm < L pos < 100 µm), (0.4 < ε s,pos < 0.6) and while the C-rate is below 4 C.The optimum power is achieved for a thin electrode (L pos < 30 µm), with high porosity and high C-rate (5 C).It is clear that the optimum energy and power cannot be achieved at the same time, hence the battery should be designed so that the power to energy ratio for a specific application is satisfactory.Finally, it should be mentioned that the developed model is not limited to the defined design variables and the responses.By changing the factors, or adding new design variables to the existing model, a new set of simulation can be run and used for further analysis.

Appendix A
The design parameters of the reference case are presented Table A1 below.The following parameters have been employed to develop the reference model.

Figure 1 .
Figure 1.Schematic of the 1D geometry of a single pair.

Figure 1 .
Figure 1.Schematic of the 1D geometry of a single pair.

Table 2 .Figure 2 .
Figure 2. Discharge curve of the 10 Ah pouch cell at different C-rates, at room temperature, = 25 with natural cooling condition.

Figure 2 .
Figure 2. Discharge curve of the 10 Ah pouch cell at different C-rates, at room temperature, T amb = 25 with natural cooling condition.

Figure 3 .
Figure 3. Normal probability plot of standardized residual for: (a) the specific energy; (b) the specific power.

Figure 3 .
Figure 3. Normal probability plot of standardized residual for: (a) the specific energy; (b) the specific power.

Figure 4 .
Figure 4. Main effect plot for the specific energy as a function of particle size ( ), electrode thickness ( ), volume fraction of the active material ( , ), C-rate, at different levels.

Figure 4 .
Figure 4. Main effect plot for the specific energy as a function of particle size (r p ), electrode thickness (L pos ), volume fraction of the active material (ε s,pos ), C-rate, at different levels.

Figure 5 .
Figure 5. Interaction effect of design parameters, particle size ( ), electrode thickness (), volume fraction of the active material ( , ) and (C-rate) at different levels on the specific energy.

Figure 5 .
Figure 5. Interaction effect of design parameters, particle size (r p ), electrode thickness (L pos ), volume fraction of the active material (ε s,pos ) and (C-rate) at different levels on the specific energy.

Figure 6 .
Figure 6.Main effect plot for the specific power, as a function of particle size ( ), electrode thickness ( ), volume fraction of the active material ( , ), C-rate, at different levels.

Figure 6 .
Figure 6.Main effect plot for the specific power, as a function of particle size (r p ), electrode thickness (L pos ), volume fraction of the active material (ε s,pos ), C-rate, at different levels.

Figure 7 .
Figure 7. Interaction effect of design parameters, particle size ( ), electrode thickness (), volume fraction of the active material ( , ) and C-rate) at different levels on the specific power.

Figure 7 .
Figure 7. Interaction effect of design parameters, particle size (r p ), electrode thickness (L pos ), volume fraction of the active material (ε s,pos ) and C-rate) at different levels on the specific power.

Figure 8 .
Figure 8. Contour plot of: (a) the specific energy as a function of ( and ); (b) the specific power as a function of ( and ); (c) the specific energy as a function of ( and , ); (d) the specific power as a function of ( and , ); (e) the specific energy as a function of ( and C-rate); (f) the specific power as a function of ( and C-rate).

Figure 8 .
Figure 8. Contour plot of: (a) the specific energy as a function of (r p and L pos ); (b) the specific power as a function of (r p and L pos ); (c) the specific energy as a function of (r p and ε s,pos ); (d) the specific power as a function of (r p and ε s,pos ); (e) the specific energy as a function of (r p and C-rate); (f) the specific power as a function of (r p and C-rate).

Figure 9 .
Figure 9. Contour plot of: (a) the specific energy as a function of ( , , ); (b) the specific power as a function of ( , C-rate); (c) the specific energy as a function of ( , , ); (d) the specific power as a function of ( , C-rate).

Figure 9 .
Figure 9. Contour plot of: (a) the specific energy as a function of (L pos , ε s,pos ); (b) the specific power as a function of (L pos , C-rate); (c) the specific energy as a function of (L pos , ε s,pos ); (d) the specific power as a function of (L pos , C-rate).
i 0 exchange current density (A•m −2 ) j Li reaction current density (A•m −2 ) k i reaction rate L thickness of the electrode (µm) m cell total weight of the cell (kg) N pairs number of electrode-pairs P cell specific power (W/kg) Q heat generation (W) R universal gas constant J• mol −1

Table 1 .
Different factors and their range, used in the factorial design-positive electrode of a 10 Ah lithium phosphate (LFP) pouch cell.

Table 3 .
The analysis of variance (ANOVA) results of the responses, R 1 (Specific energy) and R 2 (specific power) for the full factorial design.

Table A1 .
Electrochemical parameters applied for modelling of the 10 Ah LFP pouch cell.