Line data Source code
1 : /*
2 : * Cell_ECM.hpp
3 : *
4 : * Created on: 17 Dec 2019
5 : * Author(s): Jorn Reniers, Volkan Kumtepeli
6 : */
7 :
8 : #pragma once
9 :
10 : #include "../Cell.hpp"
11 : #include "State_ECM.hpp"
12 : #include "../../utility/utility.hpp"
13 : #include "../../settings/settings.hpp"
14 :
15 : #include <cstring>
16 : #include <cassert>
17 : #include <iostream>
18 : #include <fstream>
19 : #include <string>
20 : #include <cstdlib>
21 : #include <cmath>
22 : #include <algorithm>
23 : #include <array>
24 :
25 : namespace slide {
26 :
27 : template <size_t N_RC = 1>
28 : class Cell_ECM : public Cell
29 : {
30 :
31 : protected:
32 : State_ECM<N_RC> st{ settings::T_ENV, 0.5 }; //!< States T, SOC, , I, Ir, ... ;
33 : //!< parameters:
34 :
35 : std::array<double, N_RC> Rp{}, inv_tau{}; // inv_tau = 1/(RC). All initialised zero.
36 : // double Rp{ 15.8e-3 }, Cp{ 38e3 }; //!< parallel resistance and capacitance
37 : XYdata_ff OCV; //!< SOC vs voltage curve.
38 : double Rdc{ 2e-3 }; //!< DC resistance [Ohm]
39 :
40 : public:
41 : Cell_ECM();
42 : Cell_ECM(double capin, double SOCin);
43 : Cell_ECM(double capin, double SOCin, double Rdc_, std::array<double, N_RC> Rp_, std::array<double, N_RC> inv_tau_);
44 :
45 : Cell_ECM(std::string IDi, double capin, double SOCin, double Rdc_, std::array<double, N_RC> Rp_, std::array<double, N_RC> inv_tau_)
46 : : Cell_ECM(capin, SOCin, Rdc_, Rp_, inv_tau_)
47 : {
48 : ID = std::move(IDi);
49 : }
50 :
51 :
52 1 : Cell_ECM(std::string IDi, double capin, double SOCin)
53 1 : : Cell_ECM(capin, SOCin)
54 : {
55 1 : ID = std::move(IDi);
56 1 : }
57 :
58 : Cell_ECM(std::string IDi) : Cell_ECM() { ID = std::move(IDi); }
59 :
60 2200934 : inline double I() const override { return st.I(); }
61 2 : inline double getIr() { return st.Ir(); } //!< current through the parallel resistance
62 158 : inline double SOC() override { return st.SOC(); }
63 662687 : inline double T() override { return st.T(); }
64 :
65 : //!< overwrite from Cell
66 0 : std::span<double> viewStates() override { return std::span<double>(st.begin(), st.end()); }
67 34 : void getStates(getStates_t s) override { s.insert(s.end(), st.begin(), st.end()); }
68 :
69 : auto &getStateObj() { return st; }
70 :
71 : double V() override; //!< crit is an optional argument
72 : Status setStates(setStates_t s, bool checkStates = true, bool print = true) override;
73 :
74 8 : double getRtot() override { return Rdc; } //!< Return the total resistance, V = OCV - I*Rtot
75 411782 : double getThotSpot() override { return T(); }
76 0 : double getThermalSurface() override { return 0; }; //!< Not implemented?
77 2888602 : double getOCV() override { return OCV.interp(st.SOC(), settings::printBool::printCrit); } // Linear interpolation #TODO add a OCV model.
78 :
79 : Status setSOC(double SOCnew, bool checkV = true, bool print = true) override;
80 : Status setCurrent(double Inew, bool checkV = true, bool print = true) override;
81 : Status setVoltage(double Vnew, bool checkI = true, bool print = true) override;
82 :
83 0 : inline void setT(double Tnew) override { st.T() = Tnew; }
84 :
85 : virtual bool validStates(bool print = true) override;
86 : void timeStep_CC(double dt, int steps = 1) override;
87 :
88 100 : ThroughputData getThroughputs() override { return { st.time(), st.Ah(), st.Wh() }; }
89 :
90 4 : Cell_ECM<N_RC> *copy() override { return new Cell_ECM<N_RC>(*this); }
91 : };
92 :
93 : // Implementation:
94 :
95 : /**
96 : * Default constructor for Cell_ECM class template.
97 : */
98 : template <size_t N_RC>
99 126 : inline Cell_ECM<N_RC>::Cell_ECM()
100 : {
101 126 : ID = "Cell_ECM<" + std::to_string(N_RC) + ">";
102 126 : capNom = 16;
103 : /// OCV curve, dummy linear curve with 3 points from 2.0V to 4.4V
104 126 : OCV.x = slide::linspace_fix(0.0, 1.0, 3);
105 126 : OCV.y = slide::linspace_fix(VMIN(), VMAX(), 3);
106 :
107 : if constexpr (N_RC >= 1) {
108 31 : constexpr double Cp0 = 38e3; // first parallel capacitance
109 31 : Rp[0] = 15.8e-3; // fist parallel (polarisation) resistance default value.
110 31 : inv_tau[0] = 1.0 / (Rp[0] * Cp0);
111 : }
112 :
113 : if constexpr (N_RC == 2) {
114 : Rp[1] = 2.5e-3; // second parallel (polarisation) resistance default value.
115 : inv_tau[1] = 1.0 / 100.0;
116 : }
117 :
118 126 : OCV.check_is_fixed();
119 126 : cellData.initialise(*this);
120 126 : }
121 :
122 : /**
123 : * Constructor for Cell_ECM class template with given capacity and state of charge.
124 : * @param capin Capacity input.
125 : * @param SOCin State of charge input.
126 : */
127 : template <size_t N_RC>
128 2 : inline Cell_ECM<N_RC>::Cell_ECM(double capin, double SOCin) : Cell_ECM()
129 : {
130 : /// check that the input argument is valid
131 2 : if (!free::check_SOC(SOCin)) throw 10;
132 :
133 2 : st.SOC() = SOCin;
134 2 : setCapacity(capin);
135 2 : }
136 :
137 : /**
138 : * Constructor for Cell_ECM class template with given capacity, state of charge, and other parameters.
139 : * @param capin Capacity input.
140 : * @param SOCin State of charge input.
141 : * @param Rdc_ DC resistance.
142 : * @param Rp_ Array of parallel resistance values.
143 : * @param inv_tau_ Array of inverse time constant values.
144 : */
145 : template <size_t N_RC>
146 : inline Cell_ECM<N_RC>::Cell_ECM(double capin, double SOCin, double Rdc_,
147 : std::array<double, N_RC> Rp_,
148 : std::array<double, N_RC> inv_tau_)
149 : : Cell_ECM(capin, SOCin)
150 : {
151 : Rdc = Rdc_;
152 : Rp = Rp_;
153 : inv_tau = inv_tau_;
154 : }
155 :
156 :
157 : /**
158 : * sets the current
159 : *
160 : * checkV true, the voltage is checked after setting the current
161 : * if it is outside the safety limits of the cell, error 3 is thrown and the old current is restored
162 : * if it is outside the valid limits of the cell, error 2 is thrown but the new current is kept
163 : * if inside allowed Vrange, it returns the voltage
164 : * false, the voltage is not checked (function returns 0, no errors are thrown)
165 : * if no value of checkV is given, it is set to true
166 : * print controls the printing of error messages
167 : * if true, error messages are printed (if the global printing variable is high enough)
168 : * if false, no messages are printed, but the errors are still thrown
169 : * if no value, the default is true
170 : *
171 : * returns the voltage if checkV = true, else it returns 0
172 : *
173 : * THROWS
174 : * 2 checkV is true && the voltage is outside the allowed range but still in the safety range
175 : * and current is in the correct direction. I.e. if charging and V > Vmax or discharging and V < Vmin
176 : * 3 checkV is true && the voltage is outside the safety limits, old current is restored
177 : * and current is in the correct direction. I.e. if charging and V > VMAX or discharging and V < VMAX
178 : * if currents are in the 'wrong' direction (e.g. charging but V < Vmin or V < VMIN) then don't throw errors
179 : * since this current is helping to rectify the situation
180 : */
181 : template <size_t N_RC>
182 278976 : inline Status Cell_ECM<N_RC>::setCurrent(double Inew, bool checkV, bool print)
183 : {
184 278976 : const double Iold = I();
185 278976 : st.I() = Inew;
186 :
187 278976 : const auto status = checkCurrent(checkV, print);
188 :
189 278976 : if (isStatusBad(status))
190 124 : st.I() = Iold;
191 :
192 278976 : return status;
193 : }
194 :
195 : /**
196 : * Sets the voltage of the cell, updates the current accordingly, and checks the current if specified.
197 : * @param Vnew New voltage value.
198 : * @param checkI If true, checks the current after setting the voltage (default is true).
199 : * @param print If true, prints error messages (default is true).
200 : * @return The status of the operation.
201 : */
202 : template <size_t N_RC>
203 0 : inline Status Cell_ECM<N_RC>::setVoltage(double Vnew, bool checkI, bool print)
204 : {
205 0 : const double Iold = st.I();
206 : // #TODO check if V is sensible here.
207 :
208 0 : const double ocv = getOCV();
209 0 : double v_now = ocv - Vnew;
210 :
211 0 : for (size_t i{}; i < N_RC; i++)
212 0 : v_now -= Rp[i] * st.Ir(i);
213 :
214 0 : const auto Inew = v_now / Rdc;
215 :
216 0 : st.I() = Inew;
217 :
218 0 : const auto status = checkCurrent(checkI, print);
219 :
220 0 : if (isStatusBad(status))
221 0 : st.I() = Iold;
222 :
223 0 : return status;
224 : }
225 :
226 : /**
227 : * Sets the state of charge (SOC) of the cell and checks the voltage if specified.
228 : * @note This function is mainly used for testing purposes.
229 : * @param SOCnew New SOC value (must be between 0 and 1).
230 : * @param checkV If true, the voltage is checked after setting the SOC (default is true).
231 : * - If the voltage is outside the safety limits, an error is thrown and the old SOC is restored.
232 : * - If the voltage is outside the valid limits, an error is thrown but the new SOC is kept.
233 : * - If the voltage is within the allowed range, the function returns the voltage.
234 : * @param print If true, error messages are printed based on the global printing variable (default is true).
235 : * If false, no messages are printed, but errors are still thrown.
236 : * @return The status of the operation.
237 : * @throws 10 If SOCnew is illegal (values must be between 0 and 1).
238 : */
239 : template <size_t N_RC>
240 14 : inline Status Cell_ECM<N_RC>::setSOC(double SOCnew, bool checkV, bool print) //!< Also not used except test functions.
241 : {
242 14 : if (!free::check_SOC(SOCnew))
243 0 : return Status::SOC_limits_violation;
244 :
245 14 : const double SOCold = st.SOC();
246 :
247 14 : st.SOC() = SOCnew;
248 :
249 14 : if (checkV) {
250 : double v;
251 12 : const auto status = checkVoltage(v, print); //!< get the voltage Does not throw anymore!
252 :
253 12 : if (isStatusBad(status))
254 0 : st.SOC() = SOCold; //!< Restore states here.
255 :
256 12 : return status;
257 : }
258 :
259 2 : return Status::Success;
260 : }
261 :
262 : /**
263 : * Calculates the cell voltage based on the current state of the cell.
264 : * @note Throws an error if the SOC is outside the allowed range.
265 : * @tparam N_RC The number of RC-elements in the ECM model.
266 : * @return The calculated cell voltage.
267 : * @throws 1 If the SOC is outside the allowed range (passed on from linear interpolation).
268 : */
269 : template <size_t N_RC>
270 2888602 : inline double Cell_ECM<N_RC>::V()
271 : {
272 2888602 : const bool verb = settings::printBool::printCrit; //!< print if the (global) verbose-setting is above the threshold
273 : try {
274 2888602 : const double ocv = getOCV();
275 2888602 : double v_now = ocv - Rdc * st.I();
276 :
277 4297817 : for (size_t i{}; i < N_RC; i++)
278 1409215 : v_now -= Rp[i] * st.Ir(i);
279 :
280 2888602 : return v_now;
281 0 : } catch (int e) {
282 : if (verb)
283 0 : std::cerr << "ERROR in Cell_ECM::getV when getting the OCV.\n";
284 0 : return 0;
285 : }
286 : }
287 :
288 : template <size_t N_RC>
289 10 : inline Status Cell_ECM<N_RC>::setStates(setStates_t s, bool checkV, bool print)
290 : {
291 : /*
292 : */
293 10 : const auto st_old = st; //!< Back-up values.
294 :
295 10 : std::copy(s.begin(), s.begin() + st.size(), st.begin()); //!< Copy states.
296 10 : s = s.last(s.size() - st.size()); //!< Remove first Nstates elements from span.
297 :
298 10 : const Status status = free::check_Cell_states(*this, checkV);
299 :
300 10 : if (isStatusBad(status))
301 0 : st = st_old; //!< Restore states here.
302 :
303 10 : return status;
304 : }
305 :
306 : template <size_t N_RC>
307 11 : inline bool Cell_ECM<N_RC>::validStates(bool print)
308 : {
309 : /*
310 : * note: does NOT check the voltage, only whether all fields are in the allowed range
311 : * throws
312 : * 10 invalid array (wrong length)
313 : */
314 :
315 11 : const bool verb = print && (settings::printBool::printCrit); //!< print if the (global) verbose-setting is above the threshold
316 :
317 : //!< check if each value is in the allowed range
318 :
319 11 : bool range = free::check_SOC(SOC());
320 :
321 11 : if (T() < Tmin() || T() > Tmax()) {
322 0 : if (verb)
323 0 : std::cerr << "ERROR in Cell_ECM::validState, T is outside of the range, "
324 0 : << Tmin() << " <= T <= " << Tmax()
325 0 : << ", value is " << T() << '\n';
326 0 : range = false;
327 : }
328 : //!< there is no range on the current (Ir or I)
329 :
330 11 : return range;
331 : }
332 :
333 : template <size_t N_RC>
334 273135 : inline void Cell_ECM<N_RC>::timeStep_CC(double dt, int nstep)
335 : {
336 : /*
337 : * take a time step of dt seconds while keeping the current constant
338 : */
339 273135 : if (dt < 0) {
340 : if constexpr (settings::printBool::printCrit)
341 : std::cerr << "ERROR in Cell_ECM::timeStep_CC, the time step dt must be "
342 0 : << "0 or positive, but has value " << dt << '\n';
343 0 : throw 10;
344 : }
345 :
346 273135 : const auto dth = dt / 3600.0;
347 : //!< take the specified number of time steps
348 2308373 : for (int t = 0; t < nstep; t++) {
349 : //!< Using forward Euler time integration.
350 2035238 : const auto dAh = st.I() * dth;
351 2035238 : st.SOC() -= dAh / Cap();
352 :
353 3046313 : for (size_t i{}; i < N_RC; i++) // dIr/dt = (I - Ir)/(RC)
354 1011075 : st.Ir(i) += dt * inv_tau[i] * (st.I() - st.Ir(i));
355 :
356 : //!< increase the cumulative variables of this cell
357 : if constexpr (settings::data::storeCumulativeData) {
358 2035238 : st.time() += dt;
359 2035238 : st.Ah() += std::abs(dAh);
360 2035238 : st.Wh() += std::abs(dAh * V());
361 : }
362 : }
363 273135 : }
364 :
365 : using Cell_Bucket = Cell_ECM<0>;
366 :
367 : } // namespace slide
|