Article pubs.acs.org/IECR
Embedding Monotonicity in the Construction of Polynomial OpenCircuit Voltage Model for Lithium-Ion Batteries: A Semi-infinite Programming Formulation Approach Yi-Jun He,† Jia-Ni Shen,† Ji-Fu Shen,† and Zi-Feng Ma*,†,‡ †
Shanghai Electrochemical Energy Devices Research Center, Department of Chemical Engineering, Shanghai Jiao Tong University, Shanghai 200240, China ‡ Sinopoly Battery Research Center, Shanghai 200241, China ABSTRACT: Accurate modeling of open-circuit-voltage (OCV) plays important roles both in state-of-charge (SOC) estimation and state-of-health (SOH) monitoring for lithium-ion batteries (LIBs). Monotonicity violation in OCV model would lead to inaccurate SOC estimation and ineffective of incremental capacity analysis (ICA) for SOH monitoring. In this study, first-order derivative of OCV, with respect to SOC is introduced to theoretically ensure the satisfaction of monotonicity and a nonlinear semi-infinite programming (NSIP) problem is constructed to parameter estimation. A global optimization approach via restriction of the right-hand side is used to efficiently and globally optimize the NSIP. Both fitting and ICA results demonstrate the effectiveness of the proposed method. Moreover, in comparison to the traditional polynomial and sigmoid models, the NSIP polynomial model is the best choice for performing further SOC estimation and SOH monitoring. The results thus indicate that a NSIP framework for embedding prior knowledge not only provides a promising approach to automatically capture OCV-SOC monotonicity constraint in LIBs, but also serves as a universal methodology for process modeling with the requirements of embedding derivative constraints.
1. INTRODUCTION Lithium-ion batteries (LIBs), as an important type of energy storage equipment, have been widely used in portable electronic devices, electric and hybrid vehicles, etc. The commercial batteries require optimal design battery management systems (BMS) for reliable, safe, and efficient operation. Two key functions of BMS are state of charge (SOC) estimation and state of health (SOH) monitoring.1−4 Without carefully development of SOC and SOH prediction models, BMS would fail to prevent batteries from operating outside its efficient and safe region. It is well-known that the determination of open-circuitvoltage (OCV) and SOC relationship is the basic principle of almost all methods for model-based SOC estimation approach.5−9 Although the OCV-SOC look-up table has been widely employed in BMS, it often requires higher computational cost and is inconvenient for analysis. Thus, analytical models are used to capture the relationship between OCV and SOC. However, the lithium-ion staging phenomena at the graphite anode lead to the exhibition of voltage plateaus, which displays a wide flat area on the OCV-SOC curve.9−11 That means a small OCV error can result in a large SOC error. In addition, OCV curve usually has some rich features especially for low and high values of the SOC, which requires a complex model with many parameters to adequately represent the OCV.12 Therefore, it is very important to increase the prediction accuracy of OCV model for meeting the requirement of the accuracy of SOC estimation. On the other hand, since OCV data often reflect battery capacity degradation,13 an accurate OCV-SOC model is also crucial to SOH monitoring. Based on OCV data, incremental capacity analysis (ICA) technique has been validated to be © 2015 American Chemical Society
effective for detecting the gradual changes in LIB cell behavior from life cycle test.4,9,10,14 ICA is through differentiating the battery charged capacity versus the terminal voltage to obtain the incremental capacity (IC) curve. Although the extracted peaks on IC curve have been successfully used to characterize the cell degradation trend, one major difficulty in performing ICA is the sensitivity to noise in battery voltage measurement.4,9,15 This can be explained since all the peaks on an IC curve lie within the flat region of the OCV curve, direct numerical derivative computation from the data could lead to inaccurate and undesirable results.4,9 Hence, analytical derivatives based on OCV model is usually calculated to mitigate the effect of measurement noise on ICA performance based on the numerical differentiation.9 However, it should be noted that ICA method can be ineffective, in case the monotonicity constraint is violated in the OCV model. Although it has been recognized that it is very important to ensure the monotonicity constraint in the OCV model,12 it lacks a systematic approach to tackle such problems in all developed parametric OCV models. It is our motivation that the monotonicity violation in OCV model can result not only in a poor prediction accuracy of model-based SOC estimation, but also in the ineffectiveness of ICA method for SOH monitoring. In the past decade, prior knowledge-based datadriven modeling approach has attracted significant attention in the chemical engineering community, such as boiling point prediction of crude oil, polypropylene melt index prediction, Received: Revised: Accepted: Published: 3167
November 6, 2014 February 13, 2015 March 5, 2015 March 5, 2015 DOI: 10.1021/ie5044049 Ind. Eng. Chem. Res. 2015, 54, 3167−3174
Article
Industrial & Engineering Chemistry Research
⎡ ⎤ ⎡ ⎤ 1 1 + OCV(y) = x0 + x1⎢ x ⎥ ⎢ ⎥ 2 ⎣ 1 + e x6(y − x10) ⎦ ⎣ 1 + e x7(y − x11) ⎦ ⎡ ⎤ ⎡ 1 ⎤ 1 + x3⎢ ⎥ + x4⎢⎣ ⎥ + x5y x y x ( y − 1) ⎣1 + e 8 ⎦ 1+e9 ⎦
and p-xylene oxidation modeling.16−21 In this study, we propose a systematic solution framework to satisfy the monotonicity relationship in the OCV model, with respect to SOC. The first-order derivative constraint of OCV with respect to SOC is introduced to parameter estimation problem and a global optimization via restriction of the right-hand22 is used to globally solve the constructed nonlinear semi-infinite programming (NSIP) model. To our knowledge, this work is the first attempt to explicitly include monotonicity constraint into the OCV model. The remainder of this paper is structured as follows. The OCV polynomial model and its corresponding nonlinear semiinfinite programming model are first described. Then, a global optimization algorithm for the NSIP model is introduced. The effectiveness of the proposed method is next validated by the experimental data set. Finally, conclusions are given.
(5)
The above sigmoid model includes 12 parameters, where x0−5 and x6−11 are linear and nonlinear parameters, respectively. For this complex nonlinear parameter estimation problem, the initial values of parameters will significantly affect the fitting accuracy of the model when local search algorithms such as quasi-newton and trust-region are adopted. Hence, it requires numerous trials to adjust the initial values for making satisfactory parameter estimation. It is worth noting that only global optimization methods are guaranteed to find the bestpossible fit in the cases without a perfect fit. 2.2. Nonlinear Semi-infinite Programming Model. To theoretically ensure the satisfaction of monotonicity of OCV over the entire SOC range of [0,1], the first-order derivative of polynomial OCV model, with respect to the SOC shown in eq 6, is explicitly embedded into the parameter estimation problem.
2. METHODOLOGY 2.1. Open-Circuit-Voltage Model. There are several types of OCV models available in the reference, including polynomial, exponential, logarithm, and sigmoid functions.4,9 Based on the study presented in ref 9, it was found that none of the five models4 is satisfied in terms of fitting accuracy and is suitable for further IC analysis except for the proposed sigmoid model. Moreover, Hu et al. stated that a very high-order polynomial is often needed for modeling the OCV curve, which usually has some rich features, especially for low and high values of the SOC.12 In this study, we focus on investigating whether a polynomial OCV model with different orders could capture the voltage plateaus and transitions on the OCV curve of LIBs. The main reason for selecting the polynomial model is that polynomial function can be easy to efficiently implement parameter estimation procedure by comparing to exponential, logarithm, and sigmoid functions. A general polynomial OCV model is shown as
∂OCV(y) = ∂y
∑ xiyn i=1
min x
x
1 K
k=1
x̂ = H z
(1)
(7)
∀ y ∈ [0, 1]
i=1
x ∈ [x L , x U]
where xL and xU are the lower and upper bound of x, respectively. In this study, the lower and upper bounds of each decision variable are set to −1000 and 1000, respectively. The above optimization problem, consisting of a finite number of decision variables x and an infinite number of constraints, is called a semi-infinite programming (SIP) problem. It has been widely recognized that SIP problems are very challenging to solve, even in linear SIP problems. Note that, for SIP problems, the feasibility of a point x̂ needs to be established, which in principle, calls for determining the global optimum of the lower level programming (LLP), which is shown as
(2)
n
∑ xî yn− 1 y ∈ [0,1] min
(3)
i=1
where z = [OCV1, OCV2, ..., OCVK]T, H = [h1, h2, ..., hK]T with hk = [1, yk, ..., ynk]T, and H+ is the Moore−Penrose generalized inverse of H shown as
H + = (HTH )−1HT
k=1
∑ xiyn− 1 ≥ 0
Since polynomial model is linear in parameters, fast parameter estimation based on the Moore−Penrose generalized inverse can be implemented as follows: +
K
∑ (OCV(yk ) − OCVk)2
n
K
∑ (OCV(yk ) − OCVk)2
1 K
subject to
where y is the SOC, x0−n are the parameters, and n is the order of polynomial function. Given an experimental dataset with SOC-OCV pairs {yk, OCVk}Kk=1, where K is the number of experimental data, a leastsquares problem is constructed to minimize the error between predictions and experimental data. This can be shown as min
(6)
i=1
It is noted that a monotonic increasing of OCV requires the first-order derivative greater than or equal to 0. A general formulation of parameter estimation problem with monotonic increasing constraint is thus described as
n
OCV(y) = x0 +
n
∑ xiyn− 1
(8)
If the above global optimal objective value is non-negative, the feasibility of the point x̂ can be guaranteed; otherwise, vice versa. It should be noticed that if more-complex functions such as exponential and sigmoid forms are introduced to construct the OCV model, the global solution of LLP subproblem becomes more difficult and time-consuming. Hence, in this study, we choose the popular polynomial functions to
(4)
For comparison, we also include the sigmoid OCV model revised from ref 9, 3168
DOI: 10.1021/ie5044049 Ind. Eng. Chem. Res. 2015, 54, 3167−3174
Article
Industrial & Engineering Chemistry Research
3. EXPERIMENTAL SECTION A commercial HEADWAY D38120 LiFePO4 battery cell with nominal capacity of 8 A h and voltage of 3.2 V was used to test the OCV curve at five temperatures (i.e., 263, 278, 293, 308, and 323 K). The OCV test was conducted with a NEWARE Model CT-3002−5 V500A battery test system and a HARDY Model HLT4005P programmable temperature test chamber. The data acquisition system has a logging frequency of 1 Hz, and the measurement precision of both current and voltage is 0.05%. Before testing OCV, the usable capacity of battery was calibrated by a standard reference capacity test method. The reference capacity test was carried out three consecutive times to verify that the cell capacity is stable on a small performance fluctuation. In this study, after raising or lowering the cell ambient temperature to the target value, a soak period of 5 h was permitted for thermal equalization before the OCV test. First, the cell was fully charged using a constant current of 1Crate (8 A) until the voltage reached to the cutoff voltage of 3.8 V. Then, the cell was discharged using a constant current of 1Crate to a specific SOC step by step and then rested 2 h for measuring the close-to-equilibrium OCV at each SOC step. Finally, repeat the previous two steps for each specified temperature for constructing an OCV-SOC-T table. The OCV curves at five temperatures are given in Figure 1. As shown in Figure 1, in the high-temperature range of 293−
approximate the OCV-SOC curve, which could reduce the computational cost to some extent. 2.3. Global Optimization for NSIP Model. The numerical solution methods of SIP problems can be roughly divided into four classes: discretization methods, exchange methods, reduction methods, and dual methods.23,24 However, most of the methods are proposed to tackle linear SIP problems, only a few studies investigated the global optimization of SIP problems with nonconvex LLP.22,25−31 In this study, an algorithm developed by Mitsos22 is used for the global solution of NSIP problem in OCV model. The algorithm uses three subproblems, i.e., the aforementioned LLP, the lower bounding problem (LBP), and the upper bounding problem (UBP). The LBP is a relaxation of the original NSIP problem shown in eq 7 and is described as min x
1 K
K
∑ (OCV(yk ) − OCVk)2 k=1
(9)
subject to n
∑ xiyn− 1 ≥ 0
∀ y ∈ Y LBP
i=1
x ∈ [x L , x U]
where YLBP ⊂ [0,1] is a finite set. A global solution of this nonconvex LBP can give a lower bounding of the original NSIP. The UBP is an approximation of the original NSIP via restriction of the right-hand side, which can be described as min x
1 K
K
∑ (OCV(yk ) − OCVk)2 k=1
(10)
subject to n
∑ xiyn− 1 ≥ ε g
∀ y ∈ Y UBP
i=1
x ∈ [x L , x U]
where YUBP ⊂ [0,1] is a finite set and εg > 0 is the restriction parameter. The restriction parameter is decreased by a factor of r > 1 iteration by iteration. The algorithm sequentially solve the LBP, LLP, and UBP problems globally until meeting the termination criterion f UBP − f LBP ≤ εf, where f UBP and f LBP are the objective values of UBP and LBP, respectively, and εf > 0 is the termination tolerance. Three subproblems are solved globally using BARON.32 Note that, although the brute force approach may be more efficient to solve the LLP with a single variable shown in eq 8 than using BARON, we directly utilize BARON because of its relatively high computational efficiency. It has been proved that the algorithm will converge finitely under the assumption of the existence of the εf−optimal SIP− slate point. The detailed description of the global optimization algorithm of NSIP and the sensitivity of algorithm parameters on computational performance are referred to Mitsos’s study.22 In this study, the initial restriction parameter and decrease factor are set equal to 1 and 2, respectively. Both the absolute and relative tolerances are set to 0.001, and the maximum iteration number is set to 50.
Figure 1. OCV curves, with respect to SOC, at different temperatures.
323 K, two obvious voltage plateaus of 3.33 V and 3.30 V are exhibited in the SOC ranges of [0.4,0.6] and [0.7,0.95], respectively, which are resulting from lithium-ion staging phenomena. It is also found that below the temperature of 278 K, the voltage plateau starts to disappear; while at the temperature of 263 K, an obvious continuous decrease of OCV with the decrease of SOC is observed. In addition, in low and high SOC ranges, i.e., [0,0.1] and [0.95,1], there exhibits sharp OCV changes. Such a sudden change of OCV for high and low values of SOC would make OCV-SOC modeling very difficult, which usually calls for high complex function to adequately capture OCV-SOC relationship. Note that although OCV depends both on temperature and SOC, fitting an OCV model, with respect to SOC, at a specific temperature is a common treatment method. Hence, this study will fit different OCV models for different temperatures. 3169
DOI: 10.1021/ie5044049 Ind. Eng. Chem. Res. 2015, 54, 3167−3174
Article
Industrial & Engineering Chemistry Research
4. RESULTS AND DISCUSSION In this section, five sets of SOC-OCV curves shown in Figure 1 are used to validate the effectiveness of the proposed method. The global optimization method of NSIP model is implemented and solved in GAMS 24.3.3 and the global solver BARON 14.0.3 is used to solve nonconvex NLP problems. All computations are carried out on a PC with 2.60 GHz processor and 8 GB of RAM. An absolute and relative tolerance of 0.01 and 0.1 gap is used in BARON for all NSIP calculations. For clarity, in this study, our proposed OCV model by embedding monotonicity constraint is called NSIP polynomial OCV model, while the other two models are called polynomial OCV model and sigmoid OCV model. In addition, the fitting accuracy criterion is characterized by root-mean-square (RMS) analysis. 4.1. Effect of Polynomial Order on the Performance of the OCV Model. It is necessary to investigate whether highorder polynomial models could fit the complex OCV curve with plateaus and transitions. Figure 2 shows the RMS of
Figure 3. RMS analysis of the NSIP polynomial OCV models with different orders at different temperatures.
(2) the minimum RMS value is 0.4 mV for fitting the OCV curve at 263 K; (3) the minimum RMS values are ∼2.5 mV and 2 mV for fitting the OCV curves at 278 and 293 K, respectively; and (4) the minimum RMS values are ∼5 mV at both 308 and 323 K. It is found that the RMS values at low temperature are significantly lower than that at high temperature. This can be explained since OCV curves at high temperature show a much wider flat region than that at low temperature, which would make it difficult for the polynomial model to capture such characteristics by adding the monotonicity constraint. Compared to Figure 2, it is found that the fitting ability of the NSIP polynomial OCV model is usually lower than that of the traditional polynomial OCV model with the same order. That is caused by the introduction of monotonicity constraint in NSIP polynomial OCV model. However, it is worth noting that the maximum RMS value of 5 mV can meet the engineering requirements. Moreover, the NSIP polynomial OCV model can theoretically and practically satisfy the monotonicity constraint, which is greatly helpful for the effective implementation of ICA. Figure 4 shows the computational time for constructing the NSIP polynomial OCV model at different temperatures. It is found that the CPU time almost gradually increases as the polynomial order increases and the maximum CPU time is
