Tutorial 1 - Your First BattMo Model

Introduction

Let’s make your first BattMo model! ?⚡?
In this tutorial, we will build, run, and visualize a simple pseudo-two-dimensional (P2D) simulation for a Li-ion battery cell. We will first introduce how BattMo handles model parameters, then run the simulation, dashboard the results, and explore the details of the simulation output.

Define Parameters

BattMo uses JSON to manage input parameters. This allows you to easily save, document, and share complete parameter sets from specific simulations. We have used long and explicit key names for good readability. If you are new to JSON, you can learn more about it here.
For this example, we provide a sample JSON file sample_input.json describing a a simple NMC-Graphite cell [1] with LiPF6 [2] as the electrolyte in a 1D geometry according to BattMo's JSON input specification.
This sample_input.json contains the following schemas, needed to run the simulation.
In this tutorial, we simulate the constant current discharge profile for the cell with a CRate of 1 and a SoC=1.
We load and parse the JSON input file into BattMo using the command:
jsonstruct = parseBattmoJson('Examples/JsonDataFiles/sample_input.json');
This transforms the parameter data as a MATLAB structure jsonstruct that is used to setup the simulation. We can explore the structure within the MATLAB Command Window by navigating the different levels of the structure.
disp(jsonstruct)
use_thermal: 0 include_current_collectors: 0 Geometry: [1×1 struct] NegativeElectrode: [1×1 struct] PositiveElectrode: [1×1 struct] Separator: [1×1 struct] Control: [1×1 struct] Electrolyte: [1×1 struct] ThermalModel: [1×1 struct] TimeStepping: [1×1 struct] Output: [1×1 struct] SOC: 0.9900 initT: 298.1500
We can see that the structure contains various attributes for the different schema. For example, if we want to know the thickness of the negative electrode coating, we can navigate the hierarchy to find it:
disp(jsonstruct.NegativeElectrode.Coating.thickness)
6.4000e-05
Let's also check the time step parameter, we see that the maximum time for the simulation is set at ~5040 seconds with 40 timesteps.
disp(jsonstruct.TimeStepping)
totalTime: 5040 numberOfTimeSteps: 40 useRampup: 1
Unless otherwise specified, BattMo uses SI base units for physical quantities.

Run Simulation

we can now use the runBatteryJson function to simulate the constant current current discharge (half cycle) of the cell. At each time step in the simulation, BattMo computes the Li-ion concentrations in the electrodes and the electrolyte as well as the electrode and electrolyte potentials.
output = runBatteryJson(jsonstruct);
Solving timestep 01/45: -> 3 Seconds, 937 Milliseconds Solving timestep 02/45: 3 Seconds, 937 Milliseconds -> 7 Seconds, 875 Milliseconds Solving timestep 03/45: 7 Seconds, 875 Milliseconds -> 15 Seconds, 750 Milliseconds Solving timestep 04/45: 15 Seconds, 750 Milliseconds -> 31 Seconds, 500 Milliseconds Solving timestep 05/45: 31 Seconds, 500 Milliseconds -> 63 Seconds Solving timestep 06/45: 63 Seconds -> 126 Seconds Solving timestep 07/45: 126 Seconds -> 252 Seconds Solving timestep 08/45: 252 Seconds -> 378 Seconds Solving timestep 09/45: 378 Seconds -> 504 Seconds Solving timestep 10/45: 504 Seconds -> 630 Seconds Solving timestep 11/45: 630 Seconds -> 756 Seconds Solving timestep 12/45: 756 Seconds -> 882 Seconds Solving timestep 13/45: 882 Seconds -> 1008 Seconds Solving timestep 14/45: 1008 Seconds -> 1134 Seconds Solving timestep 15/45: 1134 Seconds -> 1260 Seconds Solving timestep 16/45: 1260 Seconds -> 1386 Seconds Solving timestep 17/45: 1386 Seconds -> 1512 Seconds Solving timestep 18/45: 1512 Seconds -> 1638 Seconds Solving timestep 19/45: 1638 Seconds -> 1764 Seconds Solving timestep 20/45: 1764 Seconds -> 1890 Seconds Solving timestep 21/45: 1890 Seconds -> 2016 Seconds Solving timestep 22/45: 2016 Seconds -> 2142 Seconds Solving timestep 23/45: 2142 Seconds -> 2268 Seconds Solving timestep 24/45: 2268 Seconds -> 2394 Seconds Solving timestep 25/45: 2394 Seconds -> 2520 Seconds Solving timestep 26/45: 2520 Seconds -> 2646 Seconds Solving timestep 27/45: 2646 Seconds -> 2772 Seconds Solving timestep 28/45: 2772 Seconds -> 2898 Seconds Solving timestep 29/45: 2898 Seconds -> 3024 Seconds Solving timestep 30/45: 3024 Seconds -> 3150 Seconds Solving timestep 31/45: 3150 Seconds -> 3276 Seconds Solving timestep 32/45: 3276 Seconds -> 3402 Seconds Solving timestep 33/45: 3402 Seconds -> 3528 Seconds Solving timestep 34/45: 3528 Seconds -> 1 Hour, 54 Seconds *** Simulation complete. Solved 34 control steps in 2 Seconds, 358 Milliseconds (termination triggered by stopFunction) ***

Show the Dashboard

We can dashboard the main results using plotDashboard. The left 3 columns of the dashboard shows the profiles for the main state quantities (concentration and electric potential) in the negative electrode, electrolyte, and positive electrode. The rightmost column shows the calculated cell current and voltage. Here, for example at time step 40, we can see that the Li has fully de-intercalated from the negative electrode (graphite) and intercalated almost uniformly into the positive electrode (NMC). The electric potentials at the electrodes also show oxidation at the anode and reduction at the cathode. We also see the cell voltage drop from 4.1 V to 2.4 V at the end of the discharge cycle.
plotDashboard(output.model, output.states, 'step', 40);
In the following subsections, we will explore how to access and plot this data from the simulation output.

Explore the Output

The output structure returns among other thing the model and the states.
disp(output)
model: [1×1 Battery] inputparams: [1×1 BatteryInputParams] schedule: [1×1 struct] initstate: [1×1 struct] time: [58×1 double] E: [58×1 double] I: [58×1 double] states: {58×1 cell} energy: [57×1 double] energyDensity: [57×1 double] specificEnergy: []
The model contains information about the setup of the model and initial conditions, while states contains the results of the simulation at each timestep. Plotting the simulation results requires information about the grid (i.e. what is the position where the quantity is calculated?) and the state (i.e. what is the value of the quantity in that position at a given time?).

Explore the Grid

The grid (or mesh) is one of the most used properties of the model, which can be accessed with the command:
disp(output.model.grid)
cells: [1×1 struct] faces: [1×1 struct] nodes: [1×1 struct] griddim: 1 type: 'generic'
We can see that the grid is stored as a structure with information about the cells, faces, nodes, etc. The values of the state quantities (e.g. Li-ion concentration and electric potential) are calculated at the centroids of the cells. To plot the positions of the centroids, we can use the following commands:
% get the centroid locations for the full grid
x = output.model.grid.cells.centroids;
 
% plot them in a figure
figure();
grid_axes = axes();
plot(grid_axes, x, zeros(size(x)), '+', 'LineWidth', 2)
 
% remove the y-axis ticks
yticks([]);
 
% add plot annotations
title('Grid Centroids');
xlabel('Position / m');
This shows the overall grid that is used for the model. However, BattMo models use a modular hierarchy where the overall cell model is composed of smaller submodels for electrodes, electrolyte, and current collectors. Each of these submodels has its own grid.If you would like more information about the BattMo model hierarchy, please see BattMo Model Architecture.
For example, if we want to plot the grid associated with the different submodels in different colors, we can use the following commands:
% get the centroid locations for the negative electrode (x_ne), separator
% (x_sep), and positive electrode (x_pe)
x_ne = output.model.NegativeElectrode.grid.cells.centroids;
x_sep = output.model.Separator.grid.cells.centroids;
x_pe = output.model.PositiveElectrode.grid.cells.centroids;
 
% plot them together on one figure
figure()
hold on
plot(x_ne, zeros(size(x_ne)), '+', 'LineWidth', 2)
plot(x_sep, zeros(size(x_sep)), '+k', 'LineWidth', 2)
plot(x_pe, zeros(size(x_pe)), '+r', 'LineWidth', 2)
% remove the y-axis tick marks
yticks([]);
 
% add plot annotations
title('Grid Centroids by Component')
xlabel('Position / m')
legend('Negative Electrode', 'Separator', 'Positive Electrode')
Go ahead and add the centroids for the electrolyte to the above figure to see its grid distribution in the cell. The electrolyte overlaps with the electrodes throughout the cell.

Explore the States

The values of the state quantities at each time step are stored in the states cell array. Each entry in the array describes the state of the simulation at a given timestep.
For example, we can look at the state of the simulation at timestep 10 using the command:
timestep = 40;
disp(output.states{timestep})
Electrolyte: [1×1 struct] NegativeElectrode: [1×1 struct] PositiveElectrode: [1×1 struct] Control: [1×1 struct] time: 3.5692e+03 ThermalModel: [1×1 struct]
We see that the time of the state is 3569 seconds and there are other structures containing the states of the electrodes and electrolyte. We can look into the state of the electrolyte using the command:
disp(output.states{timestep}.Electrolyte)
c: [30×1 double] phi: [30×1 double]
which shows that there are two quantities there for concentration, c, and electric potential, phi.

Plot a Result

Let’s plot the concentration in the electrolyte at timestep 40. We can plot the results using basic MATLAB commands this way:
% get the centroid locations and the values for the electrolyte
% concentration at the given timestep
x = output.model.grid.cells.centroids;
c = output.states{timestep}.Electrolyte.c;
 
% plot the concentration values at the given grid centroid locations
figure()
plot(x, c, 'LineWidth', 3)
 
% add plot annotations
xlabel('Position / m')
ylabel('Concentration / mol \cdot m^{-3}')
BattMo also includes dedicated plotting functions that will come in handy when we start working with more complex systems (e.g. P4D grids). To demonstrate how it works, we can plot the electrolyte concentration profile at every timestep using the BattMo function plotCellData:
% create a figure and use hold on to show multiple lines in one figure
figure();
hold on
 
% define a colormap for the line colors
cm = cmocean('matter', length(output.states));
 
% iterate through each of the output states and plot the calculated
% electrolyte concentration at each timestep
for timestep = 1:length(output.states)
plotCellData(output.model.grid, output.states{timestep}.Electrolyte.c, 'LineWidth', 3, 'Color', cm(timestep,:))
end
 
% add plot annotations and turn hold off
xlabel('Position / m')
ylabel('Concentration / mol \cdot m^{-3}')
hold off
You did it! ? You have run and post-processed your first BattMo simulation! But there are still a lot of exciting features to discover. Let’s keep going and explore making changes to the model. ?

References

[1] Safari, M., et al. "Multimodal physics–based aging model for life prediction of Li–ion batteries." Journal of The Electrochemical Society 156.3 (2008): A145.
[2] Xu, Meng, et al. "A pseudo three–dimensional electrochemical–thermal model of a prismatic LiFePO4 battery during discharge process." Energy 80 (2015): 303–317.