Merge branch 'ageing_issue'

bertapedret · Jan 19, 2020 · 3ccce9a · 3ccce9a
2 parents f32e528 + 6c9014c
commit 3ccce9a
Show file tree

Hide file tree

Showing 3 changed files with 28 additions and 23 deletions.
diff --git a/README.md b/README.md
@@ -1,7 +1,7 @@
 
 # LIONSIMBA - Lithium-ION SIMulation BAttery Toolbox
 
-<img style="float:right" height="330px" src="http://i68.tinypic.com/2nuhq9h.jpg" alt="Sublime's custom image"/>
+<img style="float:right" height="330px" src="https://i.ibb.co/1QhBgF3/logo.png" alt="Sublime's custom image"/>
 
  A Matlab framework based on a finite volume model suitable for Li-ion battery design, simulation, and control
 
@@ -14,7 +14,7 @@ Connect to the official web page to get the latest news
 [http://sisdin.unipv.it/labsisdin/lionsimba.php](http://sisdin.unipv.it/labsisdin/lionsimba.php)
 
 -----------------------------------------------------------------
-## Installation video on different platforms
+## Installation video on different platforms (Octave installation follows the same steps)
 
 Windows: https://www.youtube.com/watch?v=jmVh6F44C2I&t=19s
 
@@ -81,7 +81,7 @@ You can get LIONSIMBA in two ways:
 
 ### 1 - Download the latest version in zip format
 
- Download the latest zip package from [HERE](https://github.com/lionsimbatoolbox/LIONSIMBA/archive/master.zip) 
+ Download the latest zip package from [HERE](https://github.com/lionsimbatoolbox/LIONSIMBA/archive/master.zip)
 
 ### 2 - Clone the repository
  ```sh
@@ -100,6 +100,11 @@ Feel free to fork the project and modify at your best convenienve. The framework
 
 ## Changelog
 
+### Last Update 01/19/2020 - V 2.1 Released (Now supports Octave)
+**Major changes**
++ Fixed bug in the analytical initialisation of the model equations that was not allowing simulations with ageing (issue#6, thanks to mariapaygani for spotting it)
++ Implemented functions to support execution in Octave
+
 ### Last Update 04/02/2018 - V 2.0 Released
 **Major changes**
 + Constant and variable profile power input mode added

diff --git a/battery_model_files/simulator_tools/initialise_model.m b/battery_model_files/simulator_tools/initialise_model.m
@@ -45,7 +45,7 @@
 ce_init = param.ce_init*[ones(param.Np,1);ones(param.Ns,1);ones(param.Nn,1)];
 Phie_init = zeros(param.Np + param.Ns + param.Nn, 1); % designated as the ground potential for this system
 
-solverFlux  = zeros(param.Nn, 1);
+solverFlux  = zeros(param.Np+param.Nn, 1);
 film = zeros(param.Nn, 1);
 
 jflux_init = ionicFlux(ce_init, cs_star, Phis_init, Phie_init, T, solverFlux, film, param,sign_input_density,I_density);

diff --git a/startSimulation.m b/startSimulation.m
@@ -58,12 +58,12 @@
 %   LIONSIMBA: A Matlab framework based on a finite volume model suitable for Li-ion battery design, simulation, and control
 %   Copyright (C) 2016-2018 :Marcello Torchio, Lalo Magni, Davide Raimondo,
 %                            University of Pavia, 27100, Pavia, Italy
-%                            Bhushan Gopaluni, Univ. of British Columbia, 
+%                            Bhushan Gopaluni, Univ. of British Columbia,
 %                            Vancouver, BC V6T 1Z3, Canada
-%                            Richard D. Braatz, 
-%                            Massachusetts Institute of Technology, 
+%                            Richard D. Braatz,
+%                            Massachusetts Institute of Technology,
 %                            Cambridge, Massachusetts 02142, USA
-%   
+%
 %   Main code contributors to LIONSIMBA 2.0:
 %                           Ian Campbell, Krishnakumar Gopalakrishnan,
 %                           Imperial college London, London, UK
@@ -78,7 +78,7 @@
 end
 
 % Version of LIONSIMBA
-version       = '2.0';
+version       = '2.1';
 
 if(isempty(startParameters))
     % Load battery's parameters if not provided by the user
@@ -165,7 +165,7 @@
     if(param{i}.daeFormulation~=1)
         error(['Make sure that the daeFormulation flag is set to 1 for each cell parameters structure. Cell ', num2str(i),' does not respect this constraint.'])
     end
-    
+
     % When Fick's law of diffusion is used, at least 10 discretization
     % points are required. Raise an error if this condition is not met.
     if((param{i}.Nr_p<10 || param{i}.Nr_n<10) && param{i}.SolidPhaseDiffusion==3)
@@ -207,11 +207,11 @@
         param{i}.I_density = 0; % this is a dummy value
         param{i}.P_density = 0;
     end
-    
+
     % Preallocate the differentiation matrices used for the solid phase
     % diffusion in case of Fick's law.
     param{i} = solidPhaseDiffusionDifferentiationMatrices(param{i});
-    
+
     % Get the number of unknowns for the differential states
     [~, ~, ~, ~, ~, n_diff(i)] =   differentialInitialConditions(param{i});
     % Get the number of unknowns for the algebraic states
@@ -220,16 +220,16 @@
     % Store the number of differential and algebraic variables for each cell.
     param{i}.ndiff = n_diff(i);
     param{i}.nalg  = n_alg(i);
-    
+
     % The x_index variable will be used in the battery model file
     % for indexing purposes.
     param{i}.x_index    = (start_x_index:n_diff(i)+n_alg(i)+start_x_index-1);
-    
+
     param{i}.xp_index   = (start_xp_index:n_diff(i)+start_xp_index-1);
-    
+
     % Update the starting x_index value for the (possible) next cell
     start_x_index       = n_diff(i)+n_alg(i)+start_x_index;
-    
+
     start_xp_index      = n_diff(i)+n_alg(i)+start_xp_index;
 end
 
@@ -241,7 +241,7 @@
         % Solve the algebraic equations to find a set of semi-consistent initial
         % conditions for the algebraic equations. This will help the DAE solver as a warm startup.
         init_point = initialise_model(param{i});
-        
+
         % Build the initial values array for the integrator
         Yt0 = [ce_init;cs_average_init;T_init;film_init;Q_init;init_point];
         Y0  = [Y0;Yt0];
@@ -321,7 +321,7 @@
 
     % Get the model equations written in an implicit form in a symbolic way.
     [dx_tot, ~, ~] = batteryModel(0,xsym,xpsym,ida_user_data);
-    
+
     % Evaluate the Jacobian matrix. (Please refer to the Sundials guide for
     % further information about the Jacobian structure).
     J = jacobian(dx_tot,xsym) + cj*jacobian(dx_tot,xpsym);
@@ -346,7 +346,7 @@
 elseif(param{1}.UseJacobian==1 && ~isempty(param{1}.JacobianFunction))
     % If the Jacobian needs to be used and it has also been provided in the
     % parameters structure, use it directly.
-    
+
     disp('Analytical function of the Jacobian matrix provided by the user.')
     % If the Jacobian has been provided from the user, use it.
     JacFun = param{1}.JacobianFunction;
@@ -410,7 +410,7 @@
 sim_time = 0;  % Simulation time (i.e. wall) for reporting purposes only. Not used in controlling solver or time-loop
 % Loop until the integration time reaches tf.
 while(t<tf)
-    
+
     % Check if simulation stop conditions are met or not
     exit_reason = checkSimulationStopConditions(n_cells, Phis_t, cs_bar_t, T_t, param);
     % If a stop condition is reached, stop the simulation
@@ -423,7 +423,7 @@
     % refer to IDA manual for more information about syntax and its usage.
     tic
     [status, t, y]   = IDASolve(tf,'OneStep');
-    
+
     % If status > 0 it means that roots have been found during the
     % resolution of the equations
     if(status>0)
@@ -434,7 +434,7 @@
             % using the last known values of the states
             t = t*(1+1e-5);
             IDAReInit(t,y,0*y,options);
-            
+
             % Find consistent initial conditions starting from the new
             % point and keep on integrating
             [~, y, yp] = IDACalcIC(t+10,'FindAlgebraic');
@@ -475,7 +475,7 @@
             elseif (param{1}.OperatingMode==2) || (param{1}.OperatingMode==5)
                 fprintf(['Applied power density \t\t',num2str(param{1}.getPowerDensity(t,t0,tf,param{1}.extraData)),' W/m^2\n']);
             end
-				
+
 			switch param{1}.edge_values
 				% Linear interpolation
 				case 2