Identifiability of generalised Randles circuit models
1
Structural identifiability of battery equivalent circuit models
arXiv:1505.00153v1 [math.OC] 1 May 2015
S.M.Mahdi Alavi, Adam Mahdi, Stephen J. Payne and David A. Howey
Abstract This paper studies structural identifiability of battery equivalent circuit models. It is shown that the general Randles model structures, both in the continuous and discrete time domains, are locally identifiable. With respect to the battery electrochemical impedance spectroscopy, some conditions are added to the model structures, making them globally identifiable. Finally, numerical simulations are provided to demonstrate the results.
Index Terms System identification, battery management systems, batteries, equivalent circuits.
I. Introduction Safety and reliability of battery systems are two critical factors in the development and commercialization of battery-powered (hybrid) electric vehicles and renewable energies, [29], [32]. Conventional battery management systems (BMSs) are usually designed based on the voltage, current and temperature measurements, [30], [47]. However, experimental evidence shows that more information about the battery internal states is required for a satisfactory level of safety and reliability, emphasizing the need for a proper deep monitoring system, [30], [47]. A wide range of the proposed monitoring systems are based on a model of the battery. There are two widely employed models in the literature: the electrochemical and equivalent circuit models, see [31] and [34] for surveys of these models. The electrochemical models are more accurate but complex; typically consisting of a number of coupled partial differential equations with boundary conditions. In contrast, the equivalent circuit models (ECM) provide less information about the electrochemical reactions, but are much simpler and computationally more efficient. The current paper will be concerned only with ECM and the problem of the parameter identifiability, which will be introduced later. Authors are with the Department of Engineering Science, University of Oxford, Oxford, OX1 3PJ, UK. Emails: {mahdi.alavi, adam.mahdi, david.howey}@eng.ox.ac.uk
2
R1 +
Rn Cw i(t)
...
8
R
R2
C1
C2
Cn
+v1(t)-
+v2(t)-
+vn(t)-
v(t) -
+vw(t)-
Fig. 1. The battery Randles equivalent circuit model.
Randles in [1] proposed an equivalent circuit model based on the battery impedance information that has been the basis for the design of many BMSs, [7]. A modified version of the Randles model consists of an ohmic resistor, R∞ , in series with a number of parallel resistors and capacitors, and with a capacitor Cw , as shown in Figure 1. The ohmic resistor R∞ represents the conduction of charge carriers through electrolyte and metallic conductors. The resistors and capacitors in the parallel pairs represent the charge transfer resistance and the double layer capacitance, respectively, [37], [38]. The number of parallel R’s and C’s depends on how many of these pairs are required such that the frequency response of the ECM fits with the battery impedance spectra within the frequency range of interests [35], [36], [39]. The capacitor Cw , also known as the Warburg term, accounts for the diffusion process, [37], [38]; it can also be thought of as the battery state of charge [20]. The identification of Randles models (with different numbers of parallel R’s and C’s) has recently been the subject of much research [19]–[26]. The objective of this paper is to study the identifiability of the battery Randles models, that is, whether the model parameters can be quantified from the input-output data. Typically the identifiability problem is divided into two broad areas: the practical and structural identifiability. The practical identifiability considers the practical aspects of the problem that come with real data such as noise and bias [16]. In studying the structural identifiability, on the other hand, one assumes perfect data (noise-free, rich excitation signal, etc.), and therefore it is, in fact, a data-independent concept. Unidentifiable parameters can be assigned an infinite number of values yet still lead to an identical input-output data. Thus, structural identfiability is a necessary condition for paractical identifiability and parameter estimation. A number of analytical approaches to structural identifiability have been proposed, including Laplace transform (transfer function) [18], [28], Taylor series expansion [42], [44], similarity transformations [2], [9], [14], [17], [27], [43], [45], [46], and differential algebra [11], [41]. For linear systems it has been shown [2], [14], [27], [46] that controllability and observability properties are closely related with the concept of structural
3
identifiability. For example, it was shown that a single-input single-output linear time-invariant systems are structurally identifiable if and only if its observer canonical form is controllable, [46]. However, controllable and observable systems can still be unidentifiable in the general case [2]. Recently, there has been a significant interest in the identifiability analysis of battery models [3]– [6], [12], [13], [33], [40]. In [33], the structural identifiability of a five-element ECM including two capacitors was discussed by comparing the number of unknown parameters of the transfer function and the circuit. In [3], the structural indentifiability of a more general nonlinear ECM was analysed based on the observability conditions. The practical identifiability of battery electrochemical models were discussed in [6], [12], [13], [40]. In particular, it was shown that some of the electrochemical model parameters are not identifiable given the typical charge-discharge cycles [12], [13], [40]. In [6] it was shown that the shape of the charge-discharge cycles plays a crucial role in the identifiability of battery parameters. The practical identifiability problem of the Randles ECMs have been studied in [4], [5], where the authors consider models that include up to two capacitors, and the analysis is based on the Fisher information matrix, [10], [15]. In [4], a bound of estimation errors was developed by using the Cram´er-Rao theorem. In [5] a method was proposed to optimally shape the battery cycles and improve the identifiability. This paper focuses on the structural identifiability of the battery Randles model (see Figure 1) in both the continuous and discrete-time domain. It is shown that the ECM given in Figure 1 is globally identifiable for n = 1, and locally identifiable for any finite n > 1. With respect to the battery electrochemical impedance spectroscopy, some real conditions are derived for which the ECM Figure 1 becomes globally identifiable within the defined domain. The rest of the paper is organized as follows. In section II, the model structures of the Randles circuit Figure 1, in both the continuous and discrete time domains are derived, and the identifiability problem is stated. Section III describes the main results of this paper with details. In section IV, several examples are provided to validate the main results. II. The model parameterization and problem statement The state-space and transfer-function parametrized models of the Randles circuit given in Figure 1 are introduced. Define the battery current as the system input u(t) = i(t) ∈ R, the battery terminal voltage as the system output y(t) = v(t) ∈ R, and the voltages across the internal capacitors as the states vector x(t) =
4
v1 (t) · · · vn (t) vw (t)
⊤
∈ Rn+1 , and the model parameters as θ=
R∞ R1 · · · Rn C 1 · · · C n C w
⊤
(1)
,
where θ typically belongs to some open set D ⊂ R2n+2 .
A. Continuous-time state-space parametrization By using the kirchhoff’s laws, a continuous-time state-space model structure of the circuit, parametrized by θ , can be written as: d dt x(t) = A(θθ ) x(t) + B(θθ ) u(t) y(t) = C(θθ ) x(t) + D(θθ )u(t)
(2)
where A(θθ ) ∈ R(n+1)×(n+1) , B(θθ) ∈ Rn+1 , C(θθ )⊤ ∈ Rn+1 and D(θ) ∈ R are matrices that depend on the parameter vector θ , and are given by −a (θθ ) 0 ··· 0 0 1 0 −a2 (θθ ) · · · 0 0 .. .. .. .. .. A(θθ ) = . . . . . 0 · · · −an (θθ ) 0 0 0 0 ··· 0 −aw
,
b (θθ ) 1 b2 (θθ ) . B(θθ) = .. , bn (θθ ) b (θθ ) w
⊤ C(θθ ) =
1 1 .. . 1 1
,
D(θθ ) = d,
(3)
with 1 ai (θθ ) = , for i = 1, · · · , n τi aw = 0 1 bi (θθ ) = , for i = 1, · · · , n, Continuous time : Ci 1 , bw (θθ ) = Cw d(θθ) = R∞ ,
(4)
and the circuit’s time constants
τi = RiCi .
(5)
5
B. Discrete-time state-space parametrization With sample time T s and by using the forward approximation, i.e. dx(t)/dt ≈ (x((k + 1)T s ) − x(kT s ))/T s , a discrete-time state-space model structure of the Randles circuit, parametrized by θ , can be written as: x((k + 1)T s ) = A(θθ ) x(kT s ) + B(θθ) u(kT s ) y(kT s ) = C(θθ ) x(kT s ) + D(θθ ) u(kT s ),
(6)
where the structure of model’s matrices remain unchanged (see (3)), but now he parametrization is: ! Ts , for i = 1, · · · , n ai (θθ ) = − 1 − τi aw = −1 Ts (7) Discrete time : bi (θθ ) = , for i = 1, · · · , n C i Ts bw (θθ ) = Cw d(θθ ) = R . ∞
Remark 1: With respect to the electrochemical impedance spectroscopy [37], it is reasonable to assume
that τ1 < τ2 < · · · < τn and C1 < C2 < · · · < Cn which results in: an < an−1 < · · · < a1 , and bn < bn−1 < · · · < b1
(8)
In addition, experimental evidence shows that bw < bn
(9)
C. Transfer function parametrization Since in the continuous- and discrete-time the structure of the state-space models is the same and differ only by the parameterization, the transfer functions, in both cases, is given by: T (p, θ ) =
bn (θθ ) bw (θθ ) b1 (θθ ) + ···+ + + D(θθ ), p + a1 (θθ ) p + an (θθ ) p + aw (θθ )
where p represents either the Laplace (i.e. p = s) or Z-transform (i.e. p = z) operator.
(10)
6
D. Problem statement Given the model structure (3), (6), or (10), determine whether the unknown parameter vector θ is structurally identifiable. III. Main results First, the concept of structural identifiability is defined rigorously. Then, the local and global structural identifiability problem is studied for a general family of battery systems both in the continuous- and discrete-time domain, which is one of the main results of the paper. Following [46] the definition of structural identifiability is introduced. For the sake of simplicity, the term ‘structural’ is often omitted in the rest of the paper. Definition 1: Let M be a model structure with the transfer function T (p, θ ), parametrized by θ , where θ belongs to an open subset DT ⊂ Rm , and consider the equation T (p, θ ) = T (p, θ ∗ ),
for almost all p,
(11)
where θ , θ ∗ ∈ DT . Then, the model structure M is said to be - globally identifiable if (11) has a unique solution in DT , - locally identifiable if (11) has finite number of solutions in DT , - unidentifiabile if (11) has infinite number of solutions in DT . Note that in case when the system is globally identifiable, the unique solutions is θ = θ ∗ . If the system is locally identifiable, on the other hand, there always exists a smaller domain DT ′ ⊂ DT , in which it can be thought of as globally identifiable. It is worth mentioning that unidentifiability is quite common, in which case one looks for an identifiable combinations, i.e. a block of parameters that is identifiable despite the unidentifiability of individual parameters. Remark 2: Instead of defining the identifiability using the transfer function (see Definition 1), one can use (see [9], [17]) the so called coefficient map defined as follows. Consider the normalized (dividing by of the nonzero coefficient) transfer function, that is T (p, θ ) =
c0 (θθ ) + c1 (θθ )p + · · · + ck1 (θθ )pk1 , d0 (θθ ) + d1 (θθ )p + · · · + dk2 −1 (θθ )pk2 −1 + pk2
(12)
and associate to it the following coefficient map c : Rm ⊃ DT → Rk1 +k2 +1 defined as CT (θθ ) = c0 (θθ ), . . . , ck1 (θθ ), d0 (θθ ), . . . , dk2 −1 (θθ ) .
(13)
7
Thus, the model structure M is globally identifiable if the coefficient map CT is one-to-one (injective); locally identifiable if C is many-to-one; and unidentifiable if it is infinitely many-to-one. The following lemma will be used in the proof of the identifiability of the battery equivalent circuits. Lemma 1: Let M(k) be a model structure parametrized by θ = (a1 , . . . , ak , b1 , . . . , bk , d), where a j ∈ R for j = 1, . . . , k and are pairwise different, with the following transfer function T (p, θ ) =
b1 b2 bk + + ···+ + d, p + a1 p + a2 p + ak
Then the following statements hold. (a) If k = 1, then the model structure M(k) is globally identifiable. (b) if k > 1, then then the model structure M(k) is locally identifiable. Proof 1: The proof will be divided into two cases k = 1 and k > 1. Let k = 1. Then, the identifiability equation (11) is given by b∗1 b1 +d = + d∗ p + a1 p + a∗1 has a unique solutions (a1 , b1 , d) = (a∗1 , b∗1 , d ∗), which proves the global identifiability. Let k > 1. Then, the identifiability equation (11) is given by b∗k b∗1 bk b1 + ···+ +d = + · · · + + d∗. p + a1 p + ak p + a∗1 p + a∗k
(14)
We claim that equation (14) admits only finite (more precisely, k! = 1 · 2 · . . . · k) number of solutions. To prove the claim note that (14) is the equality of two rational functions, which is satisfied provided that (p + a1 ) . . . (p + ak ) = (p + a∗1 ) . . . (p + a∗k ).
(15)
Clearly, since the n distinct roots uniquely characterize a polynomial of degree n, and there are k! permutations of k roots of (p + a1 ) . . . (p + ak ), equations (15) has k! solutions. Now, let us fix the permutation s : {1, . . . k} → {1, . . . k} and consider an assignment (a1 , . . . , ak ) = (a∗s1 , . . . , a∗sk ). Since by the assumption a j for j = 1, . . . , k are pairwise distinct, the functions 1/(p + a j ) (thought of as functions of the variable p) are linearly independent. Finally, since each side of equation (14) is a linear combination of linearly independent functions, we immediately obtain that b j = b∗s j for j = 1, . . . , k. This concludes the proof of the claim and the lemma.
By slightly modifying the proof of the previous lemma (in particular equation (15)) one obtains the
8
following simple important corollary, which will be used in the proof of the main theorem of this section. Corollary 1: The conclusions of Lemma 1 are also obtained by considering the model structure M(k) parametrized by θ = (a1 , . . . , ak , b1, . . . , bk , bw, d) with the transfer function T (p, θ ) =
k X j=1
bj bw + + d. p + a j p + aw
It will be useful to introduce the following definition of the reparametrization (see e.g. [9]). Definition 2: A reparametrization of the model structure M with the coefficient map CT is a map R : Rk ⊃ DT → Rm such that Im CT ◦ R = Im CT ,
(16)
where Im denotes the image of the map. Moreover, the reparametrization is identifiable if the map CT ◦R : Rk ⊃ D → Rk1 +k2 +1 is identifiable. The main result of this section and one of the main results of the paper is the following theorem which describes the identifiability of the Randles circuit. Theorem 1: Let McRC (n) (respectively MdRC (n)) denote the continuous-time (respectively discrete-time) state-space model structures (2) (respectively (6)) with the matrices (3) parametrized by (4) (respectively (7)), where n is the number of parallel RC elements connected in series (see Figure 1). Then the following conditions hold. (a) If n = 1, then the model structures McRC (n) and MdRC (n) are globally identifiable. (b) If n > 1, then the model structures McRC (n) and MdRC (n) are locally identifiable. (c) If n > 1 and additionally the conditions (8) and (9) are assumed, then the model structures McRC (n) and MdRC (n) are globally identified. Proof 2: Consider the model structure McRC (n). Let T c (p, θ ) denote the the corresponding transfer function given by (10), where the parameter vector θ is given by (1) and is parameterized by (4). One can write T c (p, θ ) as a rational function (12) of degree k1 = k2 = n + 1 and associate to it the coefficient map CcT : R2n+2 ⊃ D → R2n+3 defined as CcT (θθ ) = c0 (θθ ), . . . , cn+1 (θθ ), d0(θθ ), . . . , dn (θθ ) .
(17)
9
Now, consider the following transfer function T (p, θ a,b ) =
n X j=1
bj bw + + d, p + aj p
(18)
where θ a,b = (a1 , . . . , an, b1 , . . . , bn , bw, d) ∈ R2n+2 and the related coefficient map CT : R2n+2 → R2n+3 in the variable θ a,b (see (13) and Remark 2), where the components of the map are not given here explicitly. By Corollary 1 (see also Lemma 1) the model structure with the transfer function (18) is globally identifiable for n = 1 (the coefficient map CT is one-to-one); and locally identifiable for n > 1 (the coefficient map CT is many-to-one). Note that the coefficient map CcT can be written as the following composition CcT = CT ◦ Rc ,
(19)
where the map Rc : R2n+2 ⊃ D → R2n+2 is the reparametrization R(θθ ) = a1 , . . . , an , b1, . . . , bn , bw , d defined by (4), see also Definition 2. Note Rc is a one-to-one map with an inverse R∞ = d,
,
Cw =
1 , bw
and
Rj =
bj , aj
Cj =
1 , bj
for
j = 1, . . . , n.
(20)
Finally, for n = 1, the coefficient map CcT is one-to-one, since the composition of two injective functions is injective. For n > 1, the map CcT is many-to-one, since the composition of a one-to-one function with many-to-one is a many-to-one function. Therefore the model structure McRC (n) is globally identifiable for n = 1 and locally identifiable for n > 1, which concludes the proof of part (a). For the model structure MdRC (n) the proof is almost exactly the same as in the continuous case and the details will be omitted. The identifiability of the family of Randles circuits given in Figure 1 has been considered due to its importance to battery systems. However, after minor modification one obtains the following corollary. Corollary 2: Let McR (n) (respectively MdR (n)) denote the continuous-time (respectively discrete-time) state-space model structures constructed as before and corresponding to the electric circuit given in Figure 1 with the Warburg term Cw eliminated. Then the following conditions hold. (a) If n = 1, then the model structures McR (n) and MdR (n) are globally identifiable. (b) If n > 1, then the model structures McR (n) and MdR (n) are locally identifiable. (c) If n > 1 and additionally the conditions (8) is assumed, then the model structures McR (n) and MdR (n) are globally identified.
10
R1
Cw 8
R
R2
+
C1
C2 i(t)
v(t) Fig. 2. The case study: 6-element R − R||C − R||C − C circuit.
IV. Simulations The main result of section III is validated through numerical simulations. Consider a battery modelled by a 6-element Randles circuit as shown in Figure 2. Based on theorem 1, the model structure (2) (res. (6)) with the matrices (3) parametrized by (4) (res. (7)) is globally identifiable provided conditions (8) and (9) are satisfied. The exact validation of the result requires parameters estimation of the model by using the battery data, i.e., mainly the current and voltage, and a reliable identification methodology, for all θ = [R∞ R1 C1 R2 C2 Cw ]⊤ ∈ D that satisfy the conditions (8) and (9), and under any random initial condition. In general, the domain of D could be all positive values in R6 . These conditions make the exact validation practically impossible. Instead, a test procedure is applied that is similar to the validation methodology in [4] based on the Monte-Carlo simulations. In this methodology, the following operating point, with respect to the battery impedance spectroscopy literature [37], is arbitrarily selected. R∞ = 0.05Ω, R1 = 0.2Ω, C1 = 0.3F, R2 = 0.4Ω, C2 = 0.6F, Cw = 300F For identification, a pseudo-random binary sequence (PRBS) signal, as shown in Figure 3(a), is applied to the battery. The battery voltage, Figure 3(b), is computed and recorded by using the models in section II. The data are sampled with sampling rate of T s = 0.01s, which is sufficiently smaller than the circuit’s time constants. Remark 3: It should be noted that the determination of the excitation signal, its shape, length, amplitude, frequency, etc., are the subject of practical identifiability which is beyond the subject of this paper as extensively surveyed in section I.
The discrete and continuous time transfer functions of the circuit are identified by using the data and Matlab identification toolbox [8]. The estimation of the circuit parameters is then extracted from the
11
10
5 Battery current
Battery voltage
5
0
0
−5
−10 0
1
2
3
4 Time(s)
5
6
7
(a) Current
8
−5 0
1
2
3
4 Time(s)
5
6
7
8
(b) Voltage
Fig. 3. Battery data.
coefficients of the transfer functions by using the inverse mapping technique that is fully discussed in the following subsection. Four scenarios are considered in this study: 1 and 2) the discrete time estimation without and with the measurement noise, and 3 and 4) the continuous time estimation without and with the measurement noise. The measurement noise is a zero mean random number with standard deviation of σ = 10−5 . Each scenario is repeated for 100 times, every run with a random initial guess of the parameters. Table I shows the means and standard deviations of the estimations. Because of the space limitation, only the probability density functions (PDFs) of the discrete-time identification method with the noisy measurement are shown in Figure 4. The results show that the parameters θ = [R∞ R1 C1 R2 C2 Cw ]⊤ are globally structurally identifiable. However, the estimation error in the continuous-time method is larger, emphasizing the need for either the modification of the identification method or the data. As expected, the estimation error with no measurement noise is much smaller than that of with the noise. In addition, estimation of the Warburg term, Cw , is very challenging because of its integral feature in the model structure. A. Inverse mapping R − R||C − R||C − C circuit Inverse mapping determines the relationship between the circuit parameters and coefficients of the transfer functions. In this subsection, the inverse mapping of the 6-element Randles model is presented. Table II provides inverse mapping of more Randles models in both the discrete and continuous time domains. Consider the R − R||C − R||C − C circuit as shown in Figure 2. By using (10), a general transfer function
12
8000
8000 R PDF
C PDF
1
7000
Actual R1
6000
6000
5000
5000
4000
4000
3000
3000
2000
2000
1000
1000
0 0.1988 0.199 0.1992 0.1994 0.1996 0.1998
0.2
0.2002 0.2004
1
7000
Actual C1
0 0.2996 0.2998
0.3
0.3002 0.3004 0.3006 0.3008 0.301 0.3012
(a)
(b)
10000
1200 R2 PDF
C2 PDF
Actual R
2
8000
Actual C
1000
2
800 6000 600 4000 400 2000
200
0 0.395
0.396
0.397
0.398
0.399
0.4
0.401
0.402
0 0.594
0.596
0.598
(c)
0.6
0.602
0.604
0.606
0.608
0.61
(d)
5
2
x 10
0.01 C PDF
R∞ PDF Actual R∞
w
Actual Cw
0.008 1.5 0.006 1 0.004 0.5 0.002
0 0.0499
0.0499
0.05
0.05
(e)
0.05
0.05
0.05
0 0
100
200
300
400
500
600
700
(f)
Fig. 4. Probability distribution functions (PDFs) of 100 times estimations, with discrete-time identification method under noisy measurement for the actual parameters given in Table I.
13
TABLE I
Discrete and continuous time estimation results with the true values given in the first column under both the noise-free and noisy measurements.
Actual Value R∞ (Ω) 0.05 R1 (Ω) 0.2 C1 (F) 0.3 R2 (Ω) 0.4 C2 (F) 0.6 Cw (F) 300
Discrete time noise free with Mean St.d. Mean 0.05 2.092e-17 0.05 0.2 3.626e-16 0.1999 0.3 2.232e-16 0.3001 0.4 2.789e-16 0.3994 0.6 1.004e-15 0.6007 300 1.714e-13 235.0347
noise St.d. 1.095e-05 1.758e-04 1.727e-04 9.267e-04 0.0016 118.3770
noise Mean 0.0489 0.0974 0.1300 0.4777 0.4007 189.1490
Continuous time free with noise St.d. Mean St.d. 0.0009 0.0491 0.0009 0.0885 0.1162 0.0894 0.1467 0.1627 0.1467 0.0671 0.4635 0.0677 0.1719 0.4373 0.1736 95.6382 209.5112 96.5761
of the circuit in the both discrete and continuous time domains has a form of: f 3 q3 + f 2 q2 + f 1 q + f 0 H(s) = 3 q + g2 q2 + g1 q + g0 where q denotes the domain’s operator, i.e. q = s in the continuous time domain and q = z in the discrete time domain. The relationships between the coefficients and the circuit parameters are given by f3 = d f2 = b1 + b2 + bw + (a1 + a2 + aw )d f1 = (a2 + aw )b1 + (a1 + aw )b2 + (a1 + a2 )bw + (a1 aw + a1 a2 + a2 aw )d f 0 = a2 aw b1 + a1 aw b2 + a1 a2 bw + a1 a2 aw d
(21)
g2 = a1 + a2 + aw g1 = a1 a2 + a1 aw + a2 aw g0 = a1 a2 aw where the parameters ai ’s, bi ’s, aw , bw and d are those defined in (4) for the continuous and in (7) for the discrete time domains. Assume that the following transfer function was identified after the identification: fˆ3 q3 + fˆ2 q2 + fˆ1 q + fˆ0 ˆ H(s) = 3 q + gˆ 2 q2 + gˆ 1 q + gˆ 0 Because the circuit has an integrator, the identification method should be set up such that a pole of the denominator is fixed at q = 0 (respectively at q = 1) in the continuous (respectively in the discrete) time domains. The identification software typically provides this type of flexibilities to fix a number of poles
14
and zeros at certain values, [8]. Then, by using the first equation of (21), (4) and (7), the estimation of R∞ is given by: Rˆ ∞ = fˆ3 The roots of q3 + gˆ 2 q2 + gˆ 1 q + gˆ 0 = 0 are aˆ 1 and aˆ 2 , and the one that has been fixed at q = 0 or 1, depending on the domain of operation. By using the condition (8), aˆ 2 < aˆ 1 . From (4) and (7), the estimation of the circuit’s time constants are obtained as follows: 1 continuous time aˆ i τˆ i = Ts 1+ˆ discrete time ai i = 1, 2.
From (21) and by using the conditions (8) and (9), the estimations of b1 , b2 and bw are obtained by solving the following set of equations
where,
1 1 1 aˆ 2 + aˆ w aˆ 1 + aˆ w aˆ 1 + aˆ 2 aˆ 2 aˆ w aˆ 1 aˆ w aˆ 1 aˆ 2
fˆ2 − (ˆa1 + aˆ 2 + aˆ w ) fˆ3 X = fˆ1 − (ˆa1 aˆ w + aˆ 2 aˆ w + aˆ 1 aˆ 2 ) fˆ2 fˆ − aˆ aˆ aˆ fˆ 0
1 2 w 3
bˆ 1 = max {X(1), X(2), X(3)} bˆ w = min {X(1), X(2), X(3)} , and bˆ 2 is inverse of the remaining element of X. The estimations of the circuit parameters are consequently obtained as follows: 1 bˆ w continuous time Cˆ w = bTˆ s discrete time w 1 bˆ i continuous time ˆ Ci = Tbˆ s discrete time i
τˆ i Rˆ i = , i = 1, 2 Cˆ i
15
V. Conclusions and Future Works In this paper, it was shown that the continuous and discrete time battery model structures associated with the Randles equivalent circuit are locally identifiable. It was further shown that the model structures become globally identifiable by taking the conditions of electrochemical impedance spectroscopy models into account. The results are confirmed through Monte-carlo based numerical simulations. The simulations show that a) continuous-time identification methods lead to larger estimation errors compared to the discrete time identification methods, b) the estimation of the Warburg term is very challenging because of its integral feature. Optimal shaping of the excitation signal through practical identifiability analysis, for the general Randles model is the subject of future research. Acknowledgements A.M. acknowledges the support of the EPSRC project EP/K036157/1. References [1] J. Randles, Kinetics of rapid electrode reactions, Discuss. Faraday Soc., Vol. 1, pp. 11-19, 1947. [2] J. Distefano, On the relationships between structural identifiability and the controllability, observability properties, IEEE Transactions on Automatic Control, Vol. 22, pp. 652, 1977. [3] M. Rausch, S. Streif, C. Pankiewitz, R. Findeisen, Nonlinear observability and identifiability of single cells in battery packs, IEEE International Conference on Control Applications (CCA), pp. 401-406, 2013. [4] A. Sharma, H.K. Fathy, Fisher identifiability analysis for a periodically-excited equivalent-circuit lithium-ion battery model, American Control Conference (ACC), pp. 274–280, 2014. [5] M.J. Rothenberger, J. Anstrom, S. Brennan, H.K. Fathy, Maximizing Parameter Identifiability of an Equivalent-Circuit Battery Model Using Optimal Periodic Input Shaping, ASME 2014 Dynamic Systems and Control Conference, pp. 1–10, 2014. [6] A.M. D’Amato, J.C. Forman, T. Ersal, A.A. Ali, J.L. Stein, H. Peng, D.S. Bernstein, Noninvasive battery-health diagnostics using retrospective-cost identification of inaccessible subsystems, ASME 2012 5th Annual Dynamic Systems and Control Conference Joint with the JSME 2012 11th Motion and Vibration Conference, pp. 1–9, 2012. [7] C.D. Rahn, C.-Y. Wang, Battery Systems Engineering, John Wiley & Sons, 2013. [8] L. Ljung, System Identification Toolbox For Use with MATLAB, The MathWorks Inc., 1988. [9] N. Meshkat, S. Sullivant, Identifiable reparametrizations of linear compartment models, Journal of Symbolic Computation, Vol. 63, pp. 46-67, 2014. [10] T. S¨oderstr¨om, P. Stoica, System identification, Prentice-Hall, 1989. [11] S. Audoly, G. Bellu, L. D’Angi`o, M.P. Saccomani, C. Cobelli, C., Global identifiability of nonlinear models of biological systems, IEEE Transactions on Biomedical Engineering, Vol. 48, 55–65, 2001. [12] S.J. Moura, N.A. Chaturvedi, M. Krsti´c, Adaptive Partial Differential Equation Observer for Battery State-of-Charge/State-of-Health Estimation Via an Electrochemical Model, Journal of Dynamic Systems, Measurement, and Control, Vol. 136, pp. 011015, 2013.
16
[13] J.C. Forman, S.J. Moura, J.L. Stein, H.K. Fathy, Genetic identification and fisher identifiability analysis of the Doyle-Fuller-Newman model from experimental cycling of a LiFePO 4 cell, Journal of Power Sources, Vol. 210, pp. 263–275, 2012. [14] K. Glover, J. Willems, Parametrizations of linear dynamical systems: Canonical forms and identifiability, IEEE Transactions on Automatic Control, Vol. 19, pp. 640–646, 1974. [15] R.L. Peeters, B. Hanzon, Symbolic computation of Fisher information matrices for parametrized state-space systems, Automatica, pp. Vol. 35, 1059–1071, 1999. [16] A. Raue, C. Kreutz, T. Maiwald, J. Bachmann, M. Schilling, U. Klingm¨uller, J. Timmer, Structural and practical identifiability analysis of partially observed dynamical models by exploiting the profile likelihood, Bioinformatics, Vol. 25, pp. 1923–1929, 2009. [17] A. Mahdi, N. Meshkat, S. Sullivant, Structural identifiability of viscoelastic mechanical systems, PLoS ONE, Vol. 9, 2014. [18] C. Cobelli, J. DiStefano, Parameter and structural identifiability concepts and ambiguities: a critical review and analysis, The American journal of physiology, Vol. 239, R7–R24, 1980. [19] J. Lee, J. Lee, O. Nam, J. Kim, B.-H. Cho, H.-S. Yun, S.-S. Choi, K. Kim, J.H. Kim, S. Jun, Modeling and Real Time Estimation of Lumped Equivalent Circuit Model of a Lithium Ion Battery, 2006 12th International Power Electronics and Motion Control Conference, 2006. [20] H. Rahimi-Eichi, F. Baronti, M.-Y. Chow, Battery modeling,observer,open-circuit voltage,parameter identification,piecewise linearization,state-of-charge (SOC) estimation, IEEE Transactions on Industrial Electronics, Vol. 61, 2053–2061, 2014. [21] A. Rahmoun, H. Biechl, Modelling of Li-ion batteries using equivalent circuit diagrams, PRZEGLAD ELEKTROTECHNICZNY,Vol. 2, 152–156, 2012. [22] S. Jiang, A Parameter Identification Method for a Battery Equivalent Circuit Model, SAE Technical Paper, 2011. [23] J. Jang, J. Yoo, Equivalent circuit evaluation method of lithium polymer battery using bode plot and numerical analysis IEEE Transactions on Energy Conversion, Vol. 26, 290–298, 2011 [24] N. Moubayed, J. Kouta, A. El-Ali, H. Dernayka, R. Outbib, Parameter identification of the lead-acid battery model, 2008 33rd IEEE Photovolatic Specialists Conference, 2008. [25] R. Al-Nazer, V. Cattin, A new optimization algorithm for a Li-Ion battery equivalent electrical circuit identification, 9th International Conference of Modeling, Optimization and Simulation - MOSIM12, 2012. [26] B. Pattipati, C. Sankavaram, K.R. Pattipati, System identification and estimation framework for pivotal automotive battery management system characteristics, IEEE Transactions on Systems, Man and Cybernetics Part C: Applications and Reviews, Vol. 41, 869–884, 2011. [27] J.M. Van Den Hof, Structural identifiability of linear compartmental systems , IEEE Transactions on Automatic Control, Vol. 43, 800–818, 1998. [28] R. Bellman, K.J. Åstr¨om, On structural identifiability, Mathematical Biosciences, Vol. 7, 329–339, 1970. [29] H. Rahimi-Eichi, U. Ojha, F. Baronti, M. Chow, Battery Management System: An Overview of Its Application in the Smart Grid and Electric Vehicles, IEEE Industrial Electronics Magazine, Vol. 7, pp. 4–16, 2013. [30] N.A. Chaturvedi, R. Klein, J. Christensen, J. Ahmed, A. Kojic, Algorithms for advanced battery-management systems, IEEE Control Systems Magazine, Vol. 30, pp 49–68, 2010. [31] A. Seaman, T.-s. Dao, J. McPhee, A survey of mathematics-based equivalent-circuit and electrochemical battery models for hybrid and electric vehicle simulation, Journal of Power Sources, Vol. 256, pp. 410–423, 2014. [32] D.H. Doughty, A. Pesaran, Vehicle Battery Safety Roadmap Guidance, National Renewable Energy Laboratory, 2012. [33] M. Sitterly, L.Y. Wang, G. Yin, C. Wang, Enhanced identification of battery models for real-time battery management, IEEE Transactions on Sustainable Energy, Vol. 2, 300–308, 2011.
17
[34] X. Hu, S. Li, H. Peng, A comparative study of equivalent circuit models for Li-ion batteries, Journal of Power Sources, Vol. 198, 359–367, 2012 [35] Y. Hu, S. Yurkovich, Y. Guezennec, B.J. Yurkovich, Electro-thermal battery model identification for automotive applications, Journal of Power Sources, Vol. 196, pp. 449–457, 2011. [36] D. Andre, M. Meiler, K. Steiner, C.H. Wimmer, T. Soczka-Guth, D.U. Sauer, Characterization of high-power lithium-ion batteries by electrochemical impedance spectroscopy. I. Experimental investigation, Journal of Power Sources, vol. 196, pp. 5334–5341, 2011. [37] E. Barsoukov, J.R. Macdonald, Impedance Spectroscopy: Theory, Experiment, and Applications, John Wiley & Sons, 2005. [38] S. Buller, M. Thele, R.W.A.A. De Doncker, E. Karden, Impedance-based simulation models of supercapacitors and li-ion batteries for power electronic applications, IEEE Transactions on Industry Applications, Vol. 41, pp 742–747, 2005. [39] C.R. Birkl, D.A. Howey, Model identification and parameter estimation for LiFePO4 batteries, Hybrid and Electric Vehicles Conference 2013 (HEVC 2013), 2013. [40] A. Schmidt, M. Bitzer, A. Imre, L. Guzzella, Experiment-driven electrochemical modeling and systematic parameterization for a lithium-ion battery cell, Journal of Power Sources, Vol. 195, pp. 5071–5080, 2010. [41] L. Ljung, T. Glad, On global identifiability for arbitrary model parametrizations, Automatica, Vol. 30, pp 265–276, 1994. [42] H. Pohjanpalo, System identifiability based on the power series expansion of the solution, Mathematical Biosciences, Vol. 41, 1978. [43] F. Anstett, G. Bloch, G. Mill´erioux, L. Denis-Vidal, Identifiability of discrete-time nonlinear systems: The local state isomorphism approach, Automatica, Vol. 44, 2884–2889, 2008 [44] M.J. Chappell, K.R. Godfrey, S. Vajda, Global identifiability of the parameters of nonlinear systems with specified inputs: A comparison of methods, Mathematical Biosciences, Vol. 102, pp 41–73, 1990. [45] S. Vajda, K.R. Godfrey, H. Rabitz, Similarity transformation approach to identifiability analysis of nonlinear compartmental models, Mathematical biosciences, 93, pp 217-248, 1989 [46] L. Ljung, System Identification Theory for the User, PTR Prentice Hall Upper Saddle River NJ, 1987. [47] M.T. Lawder, B. Suthar, P.W.C. Northrop, S. De, C.M. Hoff, O. Leitermann, M.L. Crow, S. Santhanagopalan, V.R. Subramanian, Battery Energy Storage System (BESS) and Battery Management System (BMS) for Grid-scale applications, Proceedings of the IEEE, Vol. 22, pp. 1014–1030, 2014
18
TABLE II
Inverse mapping of few equivalent circuit models. R1
R1 R i(t)
i(t) C2
C1
+
R
v(t)
v(t)
v(t)
-
-
-
-
Estimated T.F.: ˆ p+ fˆ0 ˆ H(p) = f1p+ˆ g0
Estimated T.F.: ˆ 2 ˆ p+ fˆ0 ˆ H(p) = f2p2p+ˆ+g1f1p+ˆ g0 Set up to identify a pole at
τˆ 1 =
1 aˆ 1 Ts 1+ˆa1
Rˆ ∞ = fˆ2
τˆ 1 =
1 aˆ 1 Ts 1+ˆa1
(
cont. disc. time
roots([1
τˆ i =
A=
B=
"
"
1 aˆ w
Rˆ ∞ = fˆ3
gˆ 1 gˆ 0 ]) ⇒ aˆ 1 , aˆ 2
choose
aˆ 2 < aˆ 1
1 aˆ i Ts 1+ˆai
(
AX = B
bˆ 1 = fˆ0 − aˆ 1 fˆ1
Estimated T.F.: ˆ 3 ˆ 2 ˆ p+ fˆ0 ˆ H(p) = fp3 3p+ˆ+g2f2pp2 +ˆ+g1f1p+ˆ g0 Set up to identify a pole at p = 1(res. p = 0) in the
(gˆ 1 gˆ 0 ]) ⇒ aˆ 1 , aˆ w 0 cont. time aˆ w = ( −1 disc. time gˆ 1 cont. time aˆ 1 = −ˆg0 disc. time
disc.
i(t)
discrete (res. continuous) time.
roots([1
cont.
C2
p = 1(res. p = 0) in the
Rˆ ∞ = fˆ2
aˆ 1 = gˆ 0
C1
discrete (res. continuous) time.
Rˆ ∞ = fˆ1
(
Estimated T.F.: ˆ 2 ˆ p+ fˆ0 ˆ H(p) = fp2 2p+ˆ+g1f1p+ˆ g0
R2 Cw
+
v(t)
R1 8
C1
+
R2
R1 Cw
8
C1
+
R
8
i(t)
8
R
choose
cont. disc. time
AX = B
1 aˆ 1
#
fˆ1 − (ˆa1 + aˆ w ) fˆ2 fˆ0 − aˆ 1 aˆ w fˆ2
A=
#
B=
"
"
1 aˆ 2
1 aˆ 1
gˆ 2 (gˆ 1 gˆ 0 ]) ⇒ aˆ 2 , aˆ 1 , aˆ w 0 cont. time aˆ w = −1 disc. time
roots([1
τˆ i =
aˆ 2 < aˆ 1 1 aˆ i Ts 1+ˆai
(
cont. disc.
AX = B
fˆ1 − (ˆa1 + aˆ 2 ) fˆ2 fˆ0 − aˆ 1 aˆ 2 fˆ2
1 aˆ 1 + aˆ 2 aˆ 1 aˆ 2
fˆ2 − (ˆa1 + aˆ 2 + aˆ w ) fˆ3 B = fˆ1 − (ˆa1 aˆ w + aˆ 2 aˆ w + aˆ 1 aˆ 2 ) fˆ2 fˆ0 − aˆ 1 aˆ 2 aˆ w fˆ3
1 A = aˆ 2 + aˆ w aˆ 2 aˆ w
#
#
1 aˆ 1 + aˆ w aˆ 1 aˆ w
bˆ 1 = max{X(1), X(2)}
bˆ 1 = max{X(1), X(2)}
bˆ 1 = max{X(1), X(2), X(3)}
bˆ w = min{X(1), X(2)}
bˆ 2 = min{X(1), X(2)}
bˆ w = min{X(1), X(2), X(3)} bˆ 2 is the remaining element in X
Cˆ 1 =
1 bˆ 1 Ts bˆ 1
Rˆ 1 =
cont. disc.
τˆ 1 Cˆ 1
Cˆ 1 =
Cˆ w =
1 bˆ 1 Ts bˆ 1
cont. disc.
1 bˆ w Ts bˆ w
Rˆ 1 =
cont.
Cˆ i =
1 bˆ i Ts bˆ i
cont. disc.
disc. τˆ 1 Cˆ 1
Rˆ i =
τˆ i , Cˆ i
i = 1, 2
Cˆ i =
Cˆ w = Rˆ i =
1 bˆ i Ts bˆ i
τˆ i , Cˆ i
1 bˆ w Ts bˆ w
cont. disc.
cont. disc.
i = 1, 2
Structural identifiability of battery equivalent circuit models
arXiv:1505.00153v1 [math.OC] 1 May 2015
S.M.Mahdi Alavi, Adam Mahdi, Stephen J. Payne and David A. Howey
Abstract This paper studies structural identifiability of battery equivalent circuit models. It is shown that the general Randles model structures, both in the continuous and discrete time domains, are locally identifiable. With respect to the battery electrochemical impedance spectroscopy, some conditions are added to the model structures, making them globally identifiable. Finally, numerical simulations are provided to demonstrate the results.
Index Terms System identification, battery management systems, batteries, equivalent circuits.
I. Introduction Safety and reliability of battery systems are two critical factors in the development and commercialization of battery-powered (hybrid) electric vehicles and renewable energies, [29], [32]. Conventional battery management systems (BMSs) are usually designed based on the voltage, current and temperature measurements, [30], [47]. However, experimental evidence shows that more information about the battery internal states is required for a satisfactory level of safety and reliability, emphasizing the need for a proper deep monitoring system, [30], [47]. A wide range of the proposed monitoring systems are based on a model of the battery. There are two widely employed models in the literature: the electrochemical and equivalent circuit models, see [31] and [34] for surveys of these models. The electrochemical models are more accurate but complex; typically consisting of a number of coupled partial differential equations with boundary conditions. In contrast, the equivalent circuit models (ECM) provide less information about the electrochemical reactions, but are much simpler and computationally more efficient. The current paper will be concerned only with ECM and the problem of the parameter identifiability, which will be introduced later. Authors are with the Department of Engineering Science, University of Oxford, Oxford, OX1 3PJ, UK. Emails: {mahdi.alavi, adam.mahdi, david.howey}@eng.ox.ac.uk
2
R1 +
Rn Cw i(t)
...
8
R
R2
C1
C2
Cn
+v1(t)-
+v2(t)-
+vn(t)-
v(t) -
+vw(t)-
Fig. 1. The battery Randles equivalent circuit model.
Randles in [1] proposed an equivalent circuit model based on the battery impedance information that has been the basis for the design of many BMSs, [7]. A modified version of the Randles model consists of an ohmic resistor, R∞ , in series with a number of parallel resistors and capacitors, and with a capacitor Cw , as shown in Figure 1. The ohmic resistor R∞ represents the conduction of charge carriers through electrolyte and metallic conductors. The resistors and capacitors in the parallel pairs represent the charge transfer resistance and the double layer capacitance, respectively, [37], [38]. The number of parallel R’s and C’s depends on how many of these pairs are required such that the frequency response of the ECM fits with the battery impedance spectra within the frequency range of interests [35], [36], [39]. The capacitor Cw , also known as the Warburg term, accounts for the diffusion process, [37], [38]; it can also be thought of as the battery state of charge [20]. The identification of Randles models (with different numbers of parallel R’s and C’s) has recently been the subject of much research [19]–[26]. The objective of this paper is to study the identifiability of the battery Randles models, that is, whether the model parameters can be quantified from the input-output data. Typically the identifiability problem is divided into two broad areas: the practical and structural identifiability. The practical identifiability considers the practical aspects of the problem that come with real data such as noise and bias [16]. In studying the structural identifiability, on the other hand, one assumes perfect data (noise-free, rich excitation signal, etc.), and therefore it is, in fact, a data-independent concept. Unidentifiable parameters can be assigned an infinite number of values yet still lead to an identical input-output data. Thus, structural identfiability is a necessary condition for paractical identifiability and parameter estimation. A number of analytical approaches to structural identifiability have been proposed, including Laplace transform (transfer function) [18], [28], Taylor series expansion [42], [44], similarity transformations [2], [9], [14], [17], [27], [43], [45], [46], and differential algebra [11], [41]. For linear systems it has been shown [2], [14], [27], [46] that controllability and observability properties are closely related with the concept of structural
3
identifiability. For example, it was shown that a single-input single-output linear time-invariant systems are structurally identifiable if and only if its observer canonical form is controllable, [46]. However, controllable and observable systems can still be unidentifiable in the general case [2]. Recently, there has been a significant interest in the identifiability analysis of battery models [3]– [6], [12], [13], [33], [40]. In [33], the structural identifiability of a five-element ECM including two capacitors was discussed by comparing the number of unknown parameters of the transfer function and the circuit. In [3], the structural indentifiability of a more general nonlinear ECM was analysed based on the observability conditions. The practical identifiability of battery electrochemical models were discussed in [6], [12], [13], [40]. In particular, it was shown that some of the electrochemical model parameters are not identifiable given the typical charge-discharge cycles [12], [13], [40]. In [6] it was shown that the shape of the charge-discharge cycles plays a crucial role in the identifiability of battery parameters. The practical identifiability problem of the Randles ECMs have been studied in [4], [5], where the authors consider models that include up to two capacitors, and the analysis is based on the Fisher information matrix, [10], [15]. In [4], a bound of estimation errors was developed by using the Cram´er-Rao theorem. In [5] a method was proposed to optimally shape the battery cycles and improve the identifiability. This paper focuses on the structural identifiability of the battery Randles model (see Figure 1) in both the continuous and discrete-time domain. It is shown that the ECM given in Figure 1 is globally identifiable for n = 1, and locally identifiable for any finite n > 1. With respect to the battery electrochemical impedance spectroscopy, some real conditions are derived for which the ECM Figure 1 becomes globally identifiable within the defined domain. The rest of the paper is organized as follows. In section II, the model structures of the Randles circuit Figure 1, in both the continuous and discrete time domains are derived, and the identifiability problem is stated. Section III describes the main results of this paper with details. In section IV, several examples are provided to validate the main results. II. The model parameterization and problem statement The state-space and transfer-function parametrized models of the Randles circuit given in Figure 1 are introduced. Define the battery current as the system input u(t) = i(t) ∈ R, the battery terminal voltage as the system output y(t) = v(t) ∈ R, and the voltages across the internal capacitors as the states vector x(t) =
4
v1 (t) · · · vn (t) vw (t)
⊤
∈ Rn+1 , and the model parameters as θ=
R∞ R1 · · · Rn C 1 · · · C n C w
⊤
(1)
,
where θ typically belongs to some open set D ⊂ R2n+2 .
A. Continuous-time state-space parametrization By using the kirchhoff’s laws, a continuous-time state-space model structure of the circuit, parametrized by θ , can be written as: d dt x(t) = A(θθ ) x(t) + B(θθ ) u(t) y(t) = C(θθ ) x(t) + D(θθ )u(t)
(2)
where A(θθ ) ∈ R(n+1)×(n+1) , B(θθ) ∈ Rn+1 , C(θθ )⊤ ∈ Rn+1 and D(θ) ∈ R are matrices that depend on the parameter vector θ , and are given by −a (θθ ) 0 ··· 0 0 1 0 −a2 (θθ ) · · · 0 0 .. .. .. .. .. A(θθ ) = . . . . . 0 · · · −an (θθ ) 0 0 0 0 ··· 0 −aw
,
b (θθ ) 1 b2 (θθ ) . B(θθ) = .. , bn (θθ ) b (θθ ) w
⊤ C(θθ ) =
1 1 .. . 1 1
,
D(θθ ) = d,
(3)
with 1 ai (θθ ) = , for i = 1, · · · , n τi aw = 0 1 bi (θθ ) = , for i = 1, · · · , n, Continuous time : Ci 1 , bw (θθ ) = Cw d(θθ) = R∞ ,
(4)
and the circuit’s time constants
τi = RiCi .
(5)
5
B. Discrete-time state-space parametrization With sample time T s and by using the forward approximation, i.e. dx(t)/dt ≈ (x((k + 1)T s ) − x(kT s ))/T s , a discrete-time state-space model structure of the Randles circuit, parametrized by θ , can be written as: x((k + 1)T s ) = A(θθ ) x(kT s ) + B(θθ) u(kT s ) y(kT s ) = C(θθ ) x(kT s ) + D(θθ ) u(kT s ),
(6)
where the structure of model’s matrices remain unchanged (see (3)), but now he parametrization is: ! Ts , for i = 1, · · · , n ai (θθ ) = − 1 − τi aw = −1 Ts (7) Discrete time : bi (θθ ) = , for i = 1, · · · , n C i Ts bw (θθ ) = Cw d(θθ ) = R . ∞
Remark 1: With respect to the electrochemical impedance spectroscopy [37], it is reasonable to assume
that τ1 < τ2 < · · · < τn and C1 < C2 < · · · < Cn which results in: an < an−1 < · · · < a1 , and bn < bn−1 < · · · < b1
(8)
In addition, experimental evidence shows that bw < bn
(9)
C. Transfer function parametrization Since in the continuous- and discrete-time the structure of the state-space models is the same and differ only by the parameterization, the transfer functions, in both cases, is given by: T (p, θ ) =
bn (θθ ) bw (θθ ) b1 (θθ ) + ···+ + + D(θθ ), p + a1 (θθ ) p + an (θθ ) p + aw (θθ )
where p represents either the Laplace (i.e. p = s) or Z-transform (i.e. p = z) operator.
(10)
6
D. Problem statement Given the model structure (3), (6), or (10), determine whether the unknown parameter vector θ is structurally identifiable. III. Main results First, the concept of structural identifiability is defined rigorously. Then, the local and global structural identifiability problem is studied for a general family of battery systems both in the continuous- and discrete-time domain, which is one of the main results of the paper. Following [46] the definition of structural identifiability is introduced. For the sake of simplicity, the term ‘structural’ is often omitted in the rest of the paper. Definition 1: Let M be a model structure with the transfer function T (p, θ ), parametrized by θ , where θ belongs to an open subset DT ⊂ Rm , and consider the equation T (p, θ ) = T (p, θ ∗ ),
for almost all p,
(11)
where θ , θ ∗ ∈ DT . Then, the model structure M is said to be - globally identifiable if (11) has a unique solution in DT , - locally identifiable if (11) has finite number of solutions in DT , - unidentifiabile if (11) has infinite number of solutions in DT . Note that in case when the system is globally identifiable, the unique solutions is θ = θ ∗ . If the system is locally identifiable, on the other hand, there always exists a smaller domain DT ′ ⊂ DT , in which it can be thought of as globally identifiable. It is worth mentioning that unidentifiability is quite common, in which case one looks for an identifiable combinations, i.e. a block of parameters that is identifiable despite the unidentifiability of individual parameters. Remark 2: Instead of defining the identifiability using the transfer function (see Definition 1), one can use (see [9], [17]) the so called coefficient map defined as follows. Consider the normalized (dividing by of the nonzero coefficient) transfer function, that is T (p, θ ) =
c0 (θθ ) + c1 (θθ )p + · · · + ck1 (θθ )pk1 , d0 (θθ ) + d1 (θθ )p + · · · + dk2 −1 (θθ )pk2 −1 + pk2
(12)
and associate to it the following coefficient map c : Rm ⊃ DT → Rk1 +k2 +1 defined as CT (θθ ) = c0 (θθ ), . . . , ck1 (θθ ), d0 (θθ ), . . . , dk2 −1 (θθ ) .
(13)
7
Thus, the model structure M is globally identifiable if the coefficient map CT is one-to-one (injective); locally identifiable if C is many-to-one; and unidentifiable if it is infinitely many-to-one. The following lemma will be used in the proof of the identifiability of the battery equivalent circuits. Lemma 1: Let M(k) be a model structure parametrized by θ = (a1 , . . . , ak , b1 , . . . , bk , d), where a j ∈ R for j = 1, . . . , k and are pairwise different, with the following transfer function T (p, θ ) =
b1 b2 bk + + ···+ + d, p + a1 p + a2 p + ak
Then the following statements hold. (a) If k = 1, then the model structure M(k) is globally identifiable. (b) if k > 1, then then the model structure M(k) is locally identifiable. Proof 1: The proof will be divided into two cases k = 1 and k > 1. Let k = 1. Then, the identifiability equation (11) is given by b∗1 b1 +d = + d∗ p + a1 p + a∗1 has a unique solutions (a1 , b1 , d) = (a∗1 , b∗1 , d ∗), which proves the global identifiability. Let k > 1. Then, the identifiability equation (11) is given by b∗k b∗1 bk b1 + ···+ +d = + · · · + + d∗. p + a1 p + ak p + a∗1 p + a∗k
(14)
We claim that equation (14) admits only finite (more precisely, k! = 1 · 2 · . . . · k) number of solutions. To prove the claim note that (14) is the equality of two rational functions, which is satisfied provided that (p + a1 ) . . . (p + ak ) = (p + a∗1 ) . . . (p + a∗k ).
(15)
Clearly, since the n distinct roots uniquely characterize a polynomial of degree n, and there are k! permutations of k roots of (p + a1 ) . . . (p + ak ), equations (15) has k! solutions. Now, let us fix the permutation s : {1, . . . k} → {1, . . . k} and consider an assignment (a1 , . . . , ak ) = (a∗s1 , . . . , a∗sk ). Since by the assumption a j for j = 1, . . . , k are pairwise distinct, the functions 1/(p + a j ) (thought of as functions of the variable p) are linearly independent. Finally, since each side of equation (14) is a linear combination of linearly independent functions, we immediately obtain that b j = b∗s j for j = 1, . . . , k. This concludes the proof of the claim and the lemma.
By slightly modifying the proof of the previous lemma (in particular equation (15)) one obtains the
8
following simple important corollary, which will be used in the proof of the main theorem of this section. Corollary 1: The conclusions of Lemma 1 are also obtained by considering the model structure M(k) parametrized by θ = (a1 , . . . , ak , b1, . . . , bk , bw, d) with the transfer function T (p, θ ) =
k X j=1
bj bw + + d. p + a j p + aw
It will be useful to introduce the following definition of the reparametrization (see e.g. [9]). Definition 2: A reparametrization of the model structure M with the coefficient map CT is a map R : Rk ⊃ DT → Rm such that Im CT ◦ R = Im CT ,
(16)
where Im denotes the image of the map. Moreover, the reparametrization is identifiable if the map CT ◦R : Rk ⊃ D → Rk1 +k2 +1 is identifiable. The main result of this section and one of the main results of the paper is the following theorem which describes the identifiability of the Randles circuit. Theorem 1: Let McRC (n) (respectively MdRC (n)) denote the continuous-time (respectively discrete-time) state-space model structures (2) (respectively (6)) with the matrices (3) parametrized by (4) (respectively (7)), where n is the number of parallel RC elements connected in series (see Figure 1). Then the following conditions hold. (a) If n = 1, then the model structures McRC (n) and MdRC (n) are globally identifiable. (b) If n > 1, then the model structures McRC (n) and MdRC (n) are locally identifiable. (c) If n > 1 and additionally the conditions (8) and (9) are assumed, then the model structures McRC (n) and MdRC (n) are globally identified. Proof 2: Consider the model structure McRC (n). Let T c (p, θ ) denote the the corresponding transfer function given by (10), where the parameter vector θ is given by (1) and is parameterized by (4). One can write T c (p, θ ) as a rational function (12) of degree k1 = k2 = n + 1 and associate to it the coefficient map CcT : R2n+2 ⊃ D → R2n+3 defined as CcT (θθ ) = c0 (θθ ), . . . , cn+1 (θθ ), d0(θθ ), . . . , dn (θθ ) .
(17)
9
Now, consider the following transfer function T (p, θ a,b ) =
n X j=1
bj bw + + d, p + aj p
(18)
where θ a,b = (a1 , . . . , an, b1 , . . . , bn , bw, d) ∈ R2n+2 and the related coefficient map CT : R2n+2 → R2n+3 in the variable θ a,b (see (13) and Remark 2), where the components of the map are not given here explicitly. By Corollary 1 (see also Lemma 1) the model structure with the transfer function (18) is globally identifiable for n = 1 (the coefficient map CT is one-to-one); and locally identifiable for n > 1 (the coefficient map CT is many-to-one). Note that the coefficient map CcT can be written as the following composition CcT = CT ◦ Rc ,
(19)
where the map Rc : R2n+2 ⊃ D → R2n+2 is the reparametrization R(θθ ) = a1 , . . . , an , b1, . . . , bn , bw , d defined by (4), see also Definition 2. Note Rc is a one-to-one map with an inverse R∞ = d,
,
Cw =
1 , bw
and
Rj =
bj , aj
Cj =
1 , bj
for
j = 1, . . . , n.
(20)
Finally, for n = 1, the coefficient map CcT is one-to-one, since the composition of two injective functions is injective. For n > 1, the map CcT is many-to-one, since the composition of a one-to-one function with many-to-one is a many-to-one function. Therefore the model structure McRC (n) is globally identifiable for n = 1 and locally identifiable for n > 1, which concludes the proof of part (a). For the model structure MdRC (n) the proof is almost exactly the same as in the continuous case and the details will be omitted. The identifiability of the family of Randles circuits given in Figure 1 has been considered due to its importance to battery systems. However, after minor modification one obtains the following corollary. Corollary 2: Let McR (n) (respectively MdR (n)) denote the continuous-time (respectively discrete-time) state-space model structures constructed as before and corresponding to the electric circuit given in Figure 1 with the Warburg term Cw eliminated. Then the following conditions hold. (a) If n = 1, then the model structures McR (n) and MdR (n) are globally identifiable. (b) If n > 1, then the model structures McR (n) and MdR (n) are locally identifiable. (c) If n > 1 and additionally the conditions (8) is assumed, then the model structures McR (n) and MdR (n) are globally identified.
10
R1
Cw 8
R
R2
+
C1
C2 i(t)
v(t) Fig. 2. The case study: 6-element R − R||C − R||C − C circuit.
IV. Simulations The main result of section III is validated through numerical simulations. Consider a battery modelled by a 6-element Randles circuit as shown in Figure 2. Based on theorem 1, the model structure (2) (res. (6)) with the matrices (3) parametrized by (4) (res. (7)) is globally identifiable provided conditions (8) and (9) are satisfied. The exact validation of the result requires parameters estimation of the model by using the battery data, i.e., mainly the current and voltage, and a reliable identification methodology, for all θ = [R∞ R1 C1 R2 C2 Cw ]⊤ ∈ D that satisfy the conditions (8) and (9), and under any random initial condition. In general, the domain of D could be all positive values in R6 . These conditions make the exact validation practically impossible. Instead, a test procedure is applied that is similar to the validation methodology in [4] based on the Monte-Carlo simulations. In this methodology, the following operating point, with respect to the battery impedance spectroscopy literature [37], is arbitrarily selected. R∞ = 0.05Ω, R1 = 0.2Ω, C1 = 0.3F, R2 = 0.4Ω, C2 = 0.6F, Cw = 300F For identification, a pseudo-random binary sequence (PRBS) signal, as shown in Figure 3(a), is applied to the battery. The battery voltage, Figure 3(b), is computed and recorded by using the models in section II. The data are sampled with sampling rate of T s = 0.01s, which is sufficiently smaller than the circuit’s time constants. Remark 3: It should be noted that the determination of the excitation signal, its shape, length, amplitude, frequency, etc., are the subject of practical identifiability which is beyond the subject of this paper as extensively surveyed in section I.
The discrete and continuous time transfer functions of the circuit are identified by using the data and Matlab identification toolbox [8]. The estimation of the circuit parameters is then extracted from the
11
10
5 Battery current
Battery voltage
5
0
0
−5
−10 0
1
2
3
4 Time(s)
5
6
7
(a) Current
8
−5 0
1
2
3
4 Time(s)
5
6
7
8
(b) Voltage
Fig. 3. Battery data.
coefficients of the transfer functions by using the inverse mapping technique that is fully discussed in the following subsection. Four scenarios are considered in this study: 1 and 2) the discrete time estimation without and with the measurement noise, and 3 and 4) the continuous time estimation without and with the measurement noise. The measurement noise is a zero mean random number with standard deviation of σ = 10−5 . Each scenario is repeated for 100 times, every run with a random initial guess of the parameters. Table I shows the means and standard deviations of the estimations. Because of the space limitation, only the probability density functions (PDFs) of the discrete-time identification method with the noisy measurement are shown in Figure 4. The results show that the parameters θ = [R∞ R1 C1 R2 C2 Cw ]⊤ are globally structurally identifiable. However, the estimation error in the continuous-time method is larger, emphasizing the need for either the modification of the identification method or the data. As expected, the estimation error with no measurement noise is much smaller than that of with the noise. In addition, estimation of the Warburg term, Cw , is very challenging because of its integral feature in the model structure. A. Inverse mapping R − R||C − R||C − C circuit Inverse mapping determines the relationship between the circuit parameters and coefficients of the transfer functions. In this subsection, the inverse mapping of the 6-element Randles model is presented. Table II provides inverse mapping of more Randles models in both the discrete and continuous time domains. Consider the R − R||C − R||C − C circuit as shown in Figure 2. By using (10), a general transfer function
12
8000
8000 R PDF
C PDF
1
7000
Actual R1
6000
6000
5000
5000
4000
4000
3000
3000
2000
2000
1000
1000
0 0.1988 0.199 0.1992 0.1994 0.1996 0.1998
0.2
0.2002 0.2004
1
7000
Actual C1
0 0.2996 0.2998
0.3
0.3002 0.3004 0.3006 0.3008 0.301 0.3012
(a)
(b)
10000
1200 R2 PDF
C2 PDF
Actual R
2
8000
Actual C
1000
2
800 6000 600 4000 400 2000
200
0 0.395
0.396
0.397
0.398
0.399
0.4
0.401
0.402
0 0.594
0.596
0.598
(c)
0.6
0.602
0.604
0.606
0.608
0.61
(d)
5
2
x 10
0.01 C PDF
R∞ PDF Actual R∞
w
Actual Cw
0.008 1.5 0.006 1 0.004 0.5 0.002
0 0.0499
0.0499
0.05
0.05
(e)
0.05
0.05
0.05
0 0
100
200
300
400
500
600
700
(f)
Fig. 4. Probability distribution functions (PDFs) of 100 times estimations, with discrete-time identification method under noisy measurement for the actual parameters given in Table I.
13
TABLE I
Discrete and continuous time estimation results with the true values given in the first column under both the noise-free and noisy measurements.
Actual Value R∞ (Ω) 0.05 R1 (Ω) 0.2 C1 (F) 0.3 R2 (Ω) 0.4 C2 (F) 0.6 Cw (F) 300
Discrete time noise free with Mean St.d. Mean 0.05 2.092e-17 0.05 0.2 3.626e-16 0.1999 0.3 2.232e-16 0.3001 0.4 2.789e-16 0.3994 0.6 1.004e-15 0.6007 300 1.714e-13 235.0347
noise St.d. 1.095e-05 1.758e-04 1.727e-04 9.267e-04 0.0016 118.3770
noise Mean 0.0489 0.0974 0.1300 0.4777 0.4007 189.1490
Continuous time free with noise St.d. Mean St.d. 0.0009 0.0491 0.0009 0.0885 0.1162 0.0894 0.1467 0.1627 0.1467 0.0671 0.4635 0.0677 0.1719 0.4373 0.1736 95.6382 209.5112 96.5761
of the circuit in the both discrete and continuous time domains has a form of: f 3 q3 + f 2 q2 + f 1 q + f 0 H(s) = 3 q + g2 q2 + g1 q + g0 where q denotes the domain’s operator, i.e. q = s in the continuous time domain and q = z in the discrete time domain. The relationships between the coefficients and the circuit parameters are given by f3 = d f2 = b1 + b2 + bw + (a1 + a2 + aw )d f1 = (a2 + aw )b1 + (a1 + aw )b2 + (a1 + a2 )bw + (a1 aw + a1 a2 + a2 aw )d f 0 = a2 aw b1 + a1 aw b2 + a1 a2 bw + a1 a2 aw d
(21)
g2 = a1 + a2 + aw g1 = a1 a2 + a1 aw + a2 aw g0 = a1 a2 aw where the parameters ai ’s, bi ’s, aw , bw and d are those defined in (4) for the continuous and in (7) for the discrete time domains. Assume that the following transfer function was identified after the identification: fˆ3 q3 + fˆ2 q2 + fˆ1 q + fˆ0 ˆ H(s) = 3 q + gˆ 2 q2 + gˆ 1 q + gˆ 0 Because the circuit has an integrator, the identification method should be set up such that a pole of the denominator is fixed at q = 0 (respectively at q = 1) in the continuous (respectively in the discrete) time domains. The identification software typically provides this type of flexibilities to fix a number of poles
14
and zeros at certain values, [8]. Then, by using the first equation of (21), (4) and (7), the estimation of R∞ is given by: Rˆ ∞ = fˆ3 The roots of q3 + gˆ 2 q2 + gˆ 1 q + gˆ 0 = 0 are aˆ 1 and aˆ 2 , and the one that has been fixed at q = 0 or 1, depending on the domain of operation. By using the condition (8), aˆ 2 < aˆ 1 . From (4) and (7), the estimation of the circuit’s time constants are obtained as follows: 1 continuous time aˆ i τˆ i = Ts 1+ˆ discrete time ai i = 1, 2.
From (21) and by using the conditions (8) and (9), the estimations of b1 , b2 and bw are obtained by solving the following set of equations
where,
1 1 1 aˆ 2 + aˆ w aˆ 1 + aˆ w aˆ 1 + aˆ 2 aˆ 2 aˆ w aˆ 1 aˆ w aˆ 1 aˆ 2
fˆ2 − (ˆa1 + aˆ 2 + aˆ w ) fˆ3 X = fˆ1 − (ˆa1 aˆ w + aˆ 2 aˆ w + aˆ 1 aˆ 2 ) fˆ2 fˆ − aˆ aˆ aˆ fˆ 0
1 2 w 3
bˆ 1 = max {X(1), X(2), X(3)} bˆ w = min {X(1), X(2), X(3)} , and bˆ 2 is inverse of the remaining element of X. The estimations of the circuit parameters are consequently obtained as follows: 1 bˆ w continuous time Cˆ w = bTˆ s discrete time w 1 bˆ i continuous time ˆ Ci = Tbˆ s discrete time i
τˆ i Rˆ i = , i = 1, 2 Cˆ i
15
V. Conclusions and Future Works In this paper, it was shown that the continuous and discrete time battery model structures associated with the Randles equivalent circuit are locally identifiable. It was further shown that the model structures become globally identifiable by taking the conditions of electrochemical impedance spectroscopy models into account. The results are confirmed through Monte-carlo based numerical simulations. The simulations show that a) continuous-time identification methods lead to larger estimation errors compared to the discrete time identification methods, b) the estimation of the Warburg term is very challenging because of its integral feature. Optimal shaping of the excitation signal through practical identifiability analysis, for the general Randles model is the subject of future research. Acknowledgements A.M. acknowledges the support of the EPSRC project EP/K036157/1. References [1] J. Randles, Kinetics of rapid electrode reactions, Discuss. Faraday Soc., Vol. 1, pp. 11-19, 1947. [2] J. Distefano, On the relationships between structural identifiability and the controllability, observability properties, IEEE Transactions on Automatic Control, Vol. 22, pp. 652, 1977. [3] M. Rausch, S. Streif, C. Pankiewitz, R. Findeisen, Nonlinear observability and identifiability of single cells in battery packs, IEEE International Conference on Control Applications (CCA), pp. 401-406, 2013. [4] A. Sharma, H.K. Fathy, Fisher identifiability analysis for a periodically-excited equivalent-circuit lithium-ion battery model, American Control Conference (ACC), pp. 274–280, 2014. [5] M.J. Rothenberger, J. Anstrom, S. Brennan, H.K. Fathy, Maximizing Parameter Identifiability of an Equivalent-Circuit Battery Model Using Optimal Periodic Input Shaping, ASME 2014 Dynamic Systems and Control Conference, pp. 1–10, 2014. [6] A.M. D’Amato, J.C. Forman, T. Ersal, A.A. Ali, J.L. Stein, H. Peng, D.S. Bernstein, Noninvasive battery-health diagnostics using retrospective-cost identification of inaccessible subsystems, ASME 2012 5th Annual Dynamic Systems and Control Conference Joint with the JSME 2012 11th Motion and Vibration Conference, pp. 1–9, 2012. [7] C.D. Rahn, C.-Y. Wang, Battery Systems Engineering, John Wiley & Sons, 2013. [8] L. Ljung, System Identification Toolbox For Use with MATLAB, The MathWorks Inc., 1988. [9] N. Meshkat, S. Sullivant, Identifiable reparametrizations of linear compartment models, Journal of Symbolic Computation, Vol. 63, pp. 46-67, 2014. [10] T. S¨oderstr¨om, P. Stoica, System identification, Prentice-Hall, 1989. [11] S. Audoly, G. Bellu, L. D’Angi`o, M.P. Saccomani, C. Cobelli, C., Global identifiability of nonlinear models of biological systems, IEEE Transactions on Biomedical Engineering, Vol. 48, 55–65, 2001. [12] S.J. Moura, N.A. Chaturvedi, M. Krsti´c, Adaptive Partial Differential Equation Observer for Battery State-of-Charge/State-of-Health Estimation Via an Electrochemical Model, Journal of Dynamic Systems, Measurement, and Control, Vol. 136, pp. 011015, 2013.
16
[13] J.C. Forman, S.J. Moura, J.L. Stein, H.K. Fathy, Genetic identification and fisher identifiability analysis of the Doyle-Fuller-Newman model from experimental cycling of a LiFePO 4 cell, Journal of Power Sources, Vol. 210, pp. 263–275, 2012. [14] K. Glover, J. Willems, Parametrizations of linear dynamical systems: Canonical forms and identifiability, IEEE Transactions on Automatic Control, Vol. 19, pp. 640–646, 1974. [15] R.L. Peeters, B. Hanzon, Symbolic computation of Fisher information matrices for parametrized state-space systems, Automatica, pp. Vol. 35, 1059–1071, 1999. [16] A. Raue, C. Kreutz, T. Maiwald, J. Bachmann, M. Schilling, U. Klingm¨uller, J. Timmer, Structural and practical identifiability analysis of partially observed dynamical models by exploiting the profile likelihood, Bioinformatics, Vol. 25, pp. 1923–1929, 2009. [17] A. Mahdi, N. Meshkat, S. Sullivant, Structural identifiability of viscoelastic mechanical systems, PLoS ONE, Vol. 9, 2014. [18] C. Cobelli, J. DiStefano, Parameter and structural identifiability concepts and ambiguities: a critical review and analysis, The American journal of physiology, Vol. 239, R7–R24, 1980. [19] J. Lee, J. Lee, O. Nam, J. Kim, B.-H. Cho, H.-S. Yun, S.-S. Choi, K. Kim, J.H. Kim, S. Jun, Modeling and Real Time Estimation of Lumped Equivalent Circuit Model of a Lithium Ion Battery, 2006 12th International Power Electronics and Motion Control Conference, 2006. [20] H. Rahimi-Eichi, F. Baronti, M.-Y. Chow, Battery modeling,observer,open-circuit voltage,parameter identification,piecewise linearization,state-of-charge (SOC) estimation, IEEE Transactions on Industrial Electronics, Vol. 61, 2053–2061, 2014. [21] A. Rahmoun, H. Biechl, Modelling of Li-ion batteries using equivalent circuit diagrams, PRZEGLAD ELEKTROTECHNICZNY,Vol. 2, 152–156, 2012. [22] S. Jiang, A Parameter Identification Method for a Battery Equivalent Circuit Model, SAE Technical Paper, 2011. [23] J. Jang, J. Yoo, Equivalent circuit evaluation method of lithium polymer battery using bode plot and numerical analysis IEEE Transactions on Energy Conversion, Vol. 26, 290–298, 2011 [24] N. Moubayed, J. Kouta, A. El-Ali, H. Dernayka, R. Outbib, Parameter identification of the lead-acid battery model, 2008 33rd IEEE Photovolatic Specialists Conference, 2008. [25] R. Al-Nazer, V. Cattin, A new optimization algorithm for a Li-Ion battery equivalent electrical circuit identification, 9th International Conference of Modeling, Optimization and Simulation - MOSIM12, 2012. [26] B. Pattipati, C. Sankavaram, K.R. Pattipati, System identification and estimation framework for pivotal automotive battery management system characteristics, IEEE Transactions on Systems, Man and Cybernetics Part C: Applications and Reviews, Vol. 41, 869–884, 2011. [27] J.M. Van Den Hof, Structural identifiability of linear compartmental systems , IEEE Transactions on Automatic Control, Vol. 43, 800–818, 1998. [28] R. Bellman, K.J. Åstr¨om, On structural identifiability, Mathematical Biosciences, Vol. 7, 329–339, 1970. [29] H. Rahimi-Eichi, U. Ojha, F. Baronti, M. Chow, Battery Management System: An Overview of Its Application in the Smart Grid and Electric Vehicles, IEEE Industrial Electronics Magazine, Vol. 7, pp. 4–16, 2013. [30] N.A. Chaturvedi, R. Klein, J. Christensen, J. Ahmed, A. Kojic, Algorithms for advanced battery-management systems, IEEE Control Systems Magazine, Vol. 30, pp 49–68, 2010. [31] A. Seaman, T.-s. Dao, J. McPhee, A survey of mathematics-based equivalent-circuit and electrochemical battery models for hybrid and electric vehicle simulation, Journal of Power Sources, Vol. 256, pp. 410–423, 2014. [32] D.H. Doughty, A. Pesaran, Vehicle Battery Safety Roadmap Guidance, National Renewable Energy Laboratory, 2012. [33] M. Sitterly, L.Y. Wang, G. Yin, C. Wang, Enhanced identification of battery models for real-time battery management, IEEE Transactions on Sustainable Energy, Vol. 2, 300–308, 2011.
17
[34] X. Hu, S. Li, H. Peng, A comparative study of equivalent circuit models for Li-ion batteries, Journal of Power Sources, Vol. 198, 359–367, 2012 [35] Y. Hu, S. Yurkovich, Y. Guezennec, B.J. Yurkovich, Electro-thermal battery model identification for automotive applications, Journal of Power Sources, Vol. 196, pp. 449–457, 2011. [36] D. Andre, M. Meiler, K. Steiner, C.H. Wimmer, T. Soczka-Guth, D.U. Sauer, Characterization of high-power lithium-ion batteries by electrochemical impedance spectroscopy. I. Experimental investigation, Journal of Power Sources, vol. 196, pp. 5334–5341, 2011. [37] E. Barsoukov, J.R. Macdonald, Impedance Spectroscopy: Theory, Experiment, and Applications, John Wiley & Sons, 2005. [38] S. Buller, M. Thele, R.W.A.A. De Doncker, E. Karden, Impedance-based simulation models of supercapacitors and li-ion batteries for power electronic applications, IEEE Transactions on Industry Applications, Vol. 41, pp 742–747, 2005. [39] C.R. Birkl, D.A. Howey, Model identification and parameter estimation for LiFePO4 batteries, Hybrid and Electric Vehicles Conference 2013 (HEVC 2013), 2013. [40] A. Schmidt, M. Bitzer, A. Imre, L. Guzzella, Experiment-driven electrochemical modeling and systematic parameterization for a lithium-ion battery cell, Journal of Power Sources, Vol. 195, pp. 5071–5080, 2010. [41] L. Ljung, T. Glad, On global identifiability for arbitrary model parametrizations, Automatica, Vol. 30, pp 265–276, 1994. [42] H. Pohjanpalo, System identifiability based on the power series expansion of the solution, Mathematical Biosciences, Vol. 41, 1978. [43] F. Anstett, G. Bloch, G. Mill´erioux, L. Denis-Vidal, Identifiability of discrete-time nonlinear systems: The local state isomorphism approach, Automatica, Vol. 44, 2884–2889, 2008 [44] M.J. Chappell, K.R. Godfrey, S. Vajda, Global identifiability of the parameters of nonlinear systems with specified inputs: A comparison of methods, Mathematical Biosciences, Vol. 102, pp 41–73, 1990. [45] S. Vajda, K.R. Godfrey, H. Rabitz, Similarity transformation approach to identifiability analysis of nonlinear compartmental models, Mathematical biosciences, 93, pp 217-248, 1989 [46] L. Ljung, System Identification Theory for the User, PTR Prentice Hall Upper Saddle River NJ, 1987. [47] M.T. Lawder, B. Suthar, P.W.C. Northrop, S. De, C.M. Hoff, O. Leitermann, M.L. Crow, S. Santhanagopalan, V.R. Subramanian, Battery Energy Storage System (BESS) and Battery Management System (BMS) for Grid-scale applications, Proceedings of the IEEE, Vol. 22, pp. 1014–1030, 2014
18
TABLE II
Inverse mapping of few equivalent circuit models. R1
R1 R i(t)
i(t) C2
C1
+
R
v(t)
v(t)
v(t)
-
-
-
-
Estimated T.F.: ˆ p+ fˆ0 ˆ H(p) = f1p+ˆ g0
Estimated T.F.: ˆ 2 ˆ p+ fˆ0 ˆ H(p) = f2p2p+ˆ+g1f1p+ˆ g0 Set up to identify a pole at
τˆ 1 =
1 aˆ 1 Ts 1+ˆa1
Rˆ ∞ = fˆ2
τˆ 1 =
1 aˆ 1 Ts 1+ˆa1
(
cont. disc. time
roots([1
τˆ i =
A=
B=
"
"
1 aˆ w
Rˆ ∞ = fˆ3
gˆ 1 gˆ 0 ]) ⇒ aˆ 1 , aˆ 2
choose
aˆ 2 < aˆ 1
1 aˆ i Ts 1+ˆai
(
AX = B
bˆ 1 = fˆ0 − aˆ 1 fˆ1
Estimated T.F.: ˆ 3 ˆ 2 ˆ p+ fˆ0 ˆ H(p) = fp3 3p+ˆ+g2f2pp2 +ˆ+g1f1p+ˆ g0 Set up to identify a pole at p = 1(res. p = 0) in the
(gˆ 1 gˆ 0 ]) ⇒ aˆ 1 , aˆ w 0 cont. time aˆ w = ( −1 disc. time gˆ 1 cont. time aˆ 1 = −ˆg0 disc. time
disc.
i(t)
discrete (res. continuous) time.
roots([1
cont.
C2
p = 1(res. p = 0) in the
Rˆ ∞ = fˆ2
aˆ 1 = gˆ 0
C1
discrete (res. continuous) time.
Rˆ ∞ = fˆ1
(
Estimated T.F.: ˆ 2 ˆ p+ fˆ0 ˆ H(p) = fp2 2p+ˆ+g1f1p+ˆ g0
R2 Cw
+
v(t)
R1 8
C1
+
R2
R1 Cw
8
C1
+
R
8
i(t)
8
R
choose
cont. disc. time
AX = B
1 aˆ 1
#
fˆ1 − (ˆa1 + aˆ w ) fˆ2 fˆ0 − aˆ 1 aˆ w fˆ2
A=
#
B=
"
"
1 aˆ 2
1 aˆ 1
gˆ 2 (gˆ 1 gˆ 0 ]) ⇒ aˆ 2 , aˆ 1 , aˆ w 0 cont. time aˆ w = −1 disc. time
roots([1
τˆ i =
aˆ 2 < aˆ 1 1 aˆ i Ts 1+ˆai
(
cont. disc.
AX = B
fˆ1 − (ˆa1 + aˆ 2 ) fˆ2 fˆ0 − aˆ 1 aˆ 2 fˆ2
1 aˆ 1 + aˆ 2 aˆ 1 aˆ 2
fˆ2 − (ˆa1 + aˆ 2 + aˆ w ) fˆ3 B = fˆ1 − (ˆa1 aˆ w + aˆ 2 aˆ w + aˆ 1 aˆ 2 ) fˆ2 fˆ0 − aˆ 1 aˆ 2 aˆ w fˆ3
1 A = aˆ 2 + aˆ w aˆ 2 aˆ w
#
#
1 aˆ 1 + aˆ w aˆ 1 aˆ w
bˆ 1 = max{X(1), X(2)}
bˆ 1 = max{X(1), X(2)}
bˆ 1 = max{X(1), X(2), X(3)}
bˆ w = min{X(1), X(2)}
bˆ 2 = min{X(1), X(2)}
bˆ w = min{X(1), X(2), X(3)} bˆ 2 is the remaining element in X
Cˆ 1 =
1 bˆ 1 Ts bˆ 1
Rˆ 1 =
cont. disc.
τˆ 1 Cˆ 1
Cˆ 1 =
Cˆ w =
1 bˆ 1 Ts bˆ 1
cont. disc.
1 bˆ w Ts bˆ w
Rˆ 1 =
cont.
Cˆ i =
1 bˆ i Ts bˆ i
cont. disc.
disc. τˆ 1 Cˆ 1
Rˆ i =
τˆ i , Cˆ i
i = 1, 2
Cˆ i =
Cˆ w = Rˆ i =
1 bˆ i Ts bˆ i
τˆ i , Cˆ i
1 bˆ w Ts bˆ w
cont. disc.
cont. disc.
i = 1, 2
