The change of Na storage mechanism from soft carbon to hard carbon

Abstract

Soft carbon and hard carbon materials are two of the most promising materials for sodium ion battery anodes. However, the relation between the Na storage mechanism and the microstructure of the carbon materials is not well understood. In this work, we extend our previous phenomenological model and develop an efficient mean field theory for the calculation of the OCV of carbon materials ranging from soft carbon materials to closed pore and open pore hard carbon materials. Our model utilizes a very simple model for the energetic sites within the microstructure of disordered carbon materials and aims for a unifying picture of the Na storage mechanism in carbon materials ranging from soft carbon to hard carbon. In addition, we also study the impact of the current density on the utilization of energetic sites in the microstructure of hard carbon materials.

---

1 Introduction

Sodium ion batteries (SIB) are potential alternatives to lithium ion batteries (LIB) because Na is a much more earth abundant element than Li and thus they are cheaper.^1,2 Unfortunately, graphite, the traditional anode for LIB,^3,4 is not suitable for SIB because of its low capacity^5,6 and high volume expansion.⁷ Other carbon material anodes suitable for SIB are expanded graphite,⁸ carbon nanotubes^9,10 and graphene oxide.^11,12 Soft carbon (SC), with a high surface area,^13,14 and hard carbon (HC), with a large layer distance and high porosity,^14–16 are promising SIB anode materials because they have low cost precursors¹³ and high capacity.^14,16 Moreover, the structures of SC and HC can be varied by changing the ratio of precursors¹⁷ or the synthesis temperature¹⁸ and they can be combined to create soft/hard carbon composites to improve the electrochemical performance.^19,20

The Na storage mechanisms in SC and HC materials are different due to their different microstructures. For SC, experimental studies suggest that due to the lack of closed pores Na can adsorb only on the surface of the graphite layers and may intercalate in between two graphite layers.¹⁷ Na ions can also adsorb on the defect points of defect-rich SC.^21,22 Density-functional theory calculations suggest that intercalation is energetically favourable only for layer distances larger than 0.35 nm.²³ Layer separation of this size is mostly absent in SC and also in the systems investigated in ref. 17. Therefore intercalation seems to be an unlikely storage mechanism in SC. Models for Na storage in SC have to adopt to the qualitatively different behaviour of SC to be consistent with the DFT calculations.

Additionally, the Na storage mechanism in HC materials is still under debate. Different combinations of mechanisms have been discussed: adsorption–intercalation,²⁴ adsorption–insertion,²⁵ intercalation–insertion,²⁶ adsorption–intercalation–insertion²⁷ and adsorption–intercalation–pore filling–deposition.²⁸ Moreover, the pore structure of HC materials, closed pore and open pore, also affects the Na storage mechanism. The open pore structure affects the slope capacity and the closed pore structure affects the plateau capacity.²⁹ A possible explanation for these different mechanisms is that each mechanism corresponds to a specific microstructure of HC. In other words, there is a strong relation between the microstructure of the HC material and the Na storage mechanism.

There have been several theoretical investigations of the relation between the material microstructure and the Na storage mechanism at the atomic and microstructural scale. At the atomic level, DFT was used to study Na interaction with graphite layers. Wasalathilake et al. calculated the adsorption energy of Na on a graphite surface with heteroatoms,³⁰ Cai et al. studied the intercalation energy of Na in between two graphite layers,²³ and Olsson et al. investigated the insertion energy of Na into a cylindrical pore (carbon nanotube).³¹ However, a complete understanding of Na insertion in carbon materials requires considering the interplay of different mechanisms. Na can occupy at least 4 energetically different sites in non-graphitizing carbon materials depending on the position of the Na relative to the graphite layers. The energetically different states correspond to Na inside a pore, in between two graphite layers, Na adsorbed on the graphite surface and Na adsorbed on defect sites.^32,33 Due to the simulation costs, DFT cannot cover all 4 mechanisms simultaneously.

At the microstructure scale, the correlations between the microstructure, synthesis parameters and performance are studied using machine learning (ML) algorithms. Ji et al. investigated the impact that the synthesis parameters have on the microstructure,³⁴ while Liu et al. analyzed data sets from the literature to identify performance and capacity indicators.³⁵ Oral et al. expanded the model by considering electrolytes and experimental conditions.³⁶ Also Li et al. demonstrated the correlation between synthesis conditions, microstructure and capacity.³⁷ However, one of the limitations of ML is the need for a lot of datasets.^34–37 More importantly, physical insights cannot be provided by ML techniques due to the lack of any mechanistic explanation of the identified correlations.

To the best of our knowledge there are no theoretical studies that quantify the relation between the microstructure of hard carbon materials and the Na storage mechanism in a given structure. In our paper we propose a simple mechanistic model for the relation between the microstructure parameter and Na storage mechanism. The purpose of our work is to propose a simple unified description of the OCV of disordered carbon materials which allows us to establish a correlation between microstructural changes and the OCV. In our previous work, we formulated a model for Na insertion in HC.³⁸ Here, we will modify this model to include both SC and HC. We also develop a simple mean field theory to quantify the insertion of Na in SC and HC materials. Also, the impact of non-zero currents on the distribution of ions across the internal energy states¹⁷ will be investigated.

This paper is arranged as follows. A new simple mean field theory for the calculation of chemical potentials is described in detail in Section 2. The model is parameterized using existing literature data in Section 3. The relation between the Na storage mechanism and microstructure with OCV and voltage curves at finite currents are also discussed in this section. Finally, conclusions with an outlook on future applications are given in Section 4.

2 Model

2.1 Unit cell of the carbon materials

In our previous work,³⁸ we proposed a simple model for Na insertion in HC. Here, we extend the model to include Na insertion in both SC and HC materials. A simplified cartoon of the unit cell for the carbon materials is shown in Fig. 1.

Fig. 1 Unit cell of the carbon materials with different types of Na: pore Na (blue), layer Na (red), adsorption Na (black), and defect Na (orange).

Characteristic parameters for the microstructure, the average length of a graphite layer L_a, the layer thickness L_c, layer distance d₀₀₂, and density ρ obtained from experiments are used to describe the unit cell of the carbon material. In our simplified model, the unit cell of the carbon materials has a pore confined by a triangular prism with an isosceles triangle base with a base length L_a and a side length of . Each side of the triangle has n layers with layer distance d₀₀₂ and a thickness of L_c = (n − 1)d₀₀₂. The side length is varied until the total density ρ is obtained. The change from SC to closed pore or open pore HC is determined by the density ρ, Fig. 2.

Fig. 2 The change from soft carbon (high density) to closed pore and open pore hard carbon (low density).

2.2 Mean field theory for the chemical potential

In our previous work,³⁸ we performed Monte Carlo NVT simulations to study Na insertion in HC. However, the simulation becomes very time consuming on a large scale. For example, to simulate Na storage in carbon materials with layer length L_a = 14 nm, a number of layers of n = 9, and a maximum capacity of 187 mAh g⁻¹,¹⁷ we need to model more than 10⁴ Na particles. To reduce simulation cost, a mean field theory for the calculation of the chemical potential will be developed. As in any mean field theory, the mean field theory for the chemical potential assumes that insertion of another ion takes place in the presence of all the other ions in the system. In an averaged sense this also includes local interactions. Therefore, the chemical potential is not restricted to the ideal dilute limit. Mean field theory is for the excess chemical potential which includes non-ideal behaviour (see Appendix B). On the other hand, the averaging excludes by definition the impact of clustering e.g. being close to defects. To go beyond mean field theory, the free energy functional would need to include the interplay between exclusion interaction due to the finite volume of the Na ions and attractive interactions. This could be part of future investigations.

2.2.1 4-Site model of Na insertion. As shown in Fig. 1, depending on the relative position between the Na and graphite in the carbon materials, Na can occupy at least 4 energetically different sites: ions in pores, ions intercalated in graphite layers, and adsorbed ions on the surfaces of the graphitic sheets and on defect sites. For simplicity, we call the different positions “pore Na”, “layer Na”, “adsorption Na” and “defect Na”, respectively.

The total energy in the system is:

where the subscript i = 1–4 labels the energetically different sites. N_i and ε_i are the number of Na ions in type i and the energy per Na ion of type i.

The total number of Na ions in the carbon material is:

At equilibrium, the chemical potential of each type of Na in the carbon material is equal to some value μ₀ which depends on the number of inserted Na ions. Please note that μ₀ is the chemical potential of the Na ion plus the chemical potential of the electron which is inserted in the graphite together with the Na ion.

μ_i(N_i) = μ₀, i = 1–4 (3)

By varying μ₀, the state of charge (SOC) is varied smoothly from 0 to 1. In this work, the chemical potential μ₀ is in the range μ₀ = −3 to 0 eV.

The chemical potential of each type is:

μ_i = μ^id_i + μ^ex_i (4)

where μ^id_i is the ideal chemical potential of type i. The expression for μ^id_i is shown in Appendix A.

μ ^ex _i is the excess chemical potential of type i.

From eqn (1), the excess chemical potential μ^ex_i is:

The expression of ε_i for each type will be given in detail in Section 2.2.2.

We use eqn (2) and (3) to obtain the capacity at each value of μ₀. The OCV voltage φ is calculated as:

eφ = ε⁰ − μ_NaC = ε⁰ − (μ^id_NaC + μ^ex_NaC) (6)

where ε⁰ (eV) is the reference potential of the hard carbon vs. Na/Na⁺.

μ ^id_NaC and μ^ex_NaC (eV) are the ideal and average excess chemical potential of Na in the carbon materials. The expression and detailed derivation of μ^id_NaC and μ^ex_NaC are shown in Appendices A and B, respectively.

2.2.2 Energy of each type of Na.
a. Defect Na. In ref. 38 we assumed all defect sites to have the same energy ε₄ and a limiting number N^max₄ of defects. The number of accessible defect sites N₄ can be obtained as:
b. Adsorption Na. The energy of adsorption Na can be approximated as:²⁶where ε⁰₃ is the energy at which Na starts to adsorb on the surface of graphite.³⁹ For NaClO₄ in experiments,¹⁷ we have ε⁰₃ = −2.52 eV.³⁹N^max₃ is the maximum number of adsorption sites for Na. a₃ is the exponent number for fitting with experiments. Higher values of a₃ lead to a larger occupation of adsorption sites at high voltage, Fig. 3.

Fig. 3 The OCV results for the SOC of Na adsorption with different values of a₃.

c. Layer Na. Since the exact intercalation energy of layer Na is not known, the impact of interaction between the Na ions on the intercalation energy is assumed to be a polynomial function of layer Na:

ε₂ = ε₂₀ + ε₂₁x₂ + ε₂₂x₂² + ⋯ + o(x₂ⁿ) (9)

where ε₂₀ is the dilute limit of intercalation energy of Na in between two graphite layers. This can be considered as the interaction between Na and the graphite layers. ε_2n is the energy scale for interaction between the layer Na ions. is the normalized equation for layer Na.³⁸

In At the low density limit (less than one Na per 24 carbon units, x₂ < 1), the interaction energy with the Na ions can be assumed to be linear:

ε₂ = ε₂₀ + ε₂₁x₂ + o(x₂). (10)

For higher concentrations of layer Na, we include a higher degree polynomial function, eqn (9).

d. Pore Na. In experiments, the pores in carbon materials can be open or closed.²⁹ In our model, Fig. 4, we choose the ratio as the descriptor for the type of pore. For small values the pore is closed. As the value increases pore opens up and like an adsorption site, since the surface of the layers become accessible to the Na. The unit cell of the carbon materials turns into that of a soft carbon material which only has adsorption sites and layer sites, Fig. 2. The closed pore affects the plateau capacity since the Na ions have to intercalate. Open pores affect the sloping capacity²⁹ as Na gets more and more access to the surface of the layers where they can adsorb.

Fig. 4 The change from closed pore Na to open pore Na and the approach to Na adsorption depending on the ratio .

To account for the impact of the pore size, we modify eqn (8) with the pore parameter a_p:

The pore parameter a_p is defined as:For when the pore parameter a_p → 1 and eqn (11) changes back to eqn (8).

The scale r_p of the pore radius is estimated from an idealized picture of a pore in the carbon materials, where it is assumed that the pore is formed as a perfect triangle (in 2D) with the width of the graphitic sheets as boundaries, as shown in Fig. 5:

whereand

Fig. 5 The carbon material pore.

2.3 Input energy parameter in mean field theory

For convenience, the input parameters are shown in Table 1.

Table 1 Input energy parameters in mean field theory

Type of Na	Parameter	Equation	Value
	μ 0	(3)	−3 to 0 eV
Adsorption Na	ε ^03	(8)	−2.52 eV^39
	a 3	(8)	Parameterizing to fit with experiments
	N ^max3	(8)
Layer Na	ε 20	(9)
	ε 21	(9)
Pore Na	N ^max1	(10)

---

As shown in Table 1, we have total of 7 input parameters for our phenomenological model. In these parameters, we fix the value of ε⁰₃ to −2.52 eV from DFT calculations,³⁹ and vary μ₀ in a range from −3 to 0 eV. Other parameters are fitted to experimental values, giving insights into the relation between the material microstructure and Na storage mechanism (see Section 3.2.2). All parameters are fitted to experiment. There are no derived parameters in this model.

3 Results and discussion

3.1 Experiment data

In experiments,¹⁷ a series materials ranging from SC to HC are synthesized by changing the molar ratio of terephthaloyl chloride (TPC) and perylene (PER), . The molar ratio is varied from x = 0.5 (SC) to x = 10 (HC). The microstructures of the carbon materials are shown in Table 2.

Table 2 The microstructure properties of the carbon materials obtained from experiments¹⁷

x	L a, nm	L c, nm	d 002, nm	n	ρ, cm^3 g^−1	S BET, m^2 g^−1
0.5	14.03	2.655	0.345	8.7	2.29	8.2
2	12.02	1.768	0.352	6.0	1.65	9.9
4	9.81	1.382	0.376	4.6	1.50	3.6
6	8.82	1.392	0.379	4.7	1.34	9.0
10	8.78	1.305	0.382	4.4	1.35	6.4

---

3.2 Microstructure effect

3.2.1 Comparison of the model with experiment. In Fig. 6, the OCV results obtained from experiment are compared to the results obtained from our model in 5 cases when x between is 0.5–10 and the error between experiment and our model Δϕ = |ϕ_exp − ϕ_mod| is also shown. The capacities obtained from our model are 182.0 mAh g⁻¹, 273.8 mAh g⁻¹, 314.6 mAh g⁻¹, 401.7 mAh g⁻¹ and 401.7 mAh g⁻¹; compared with 181.7 mAh g⁻¹, 272.7 mAh g⁻¹, 307.0 mAh g⁻¹, 396.4 mAh g⁻¹ and 399.3 mAh g⁻¹ obtained from experiment, respectively. The highest error is at low capacity when x = 0.5 and x = 2, where the greatest contribution comes from adsorption Na.³⁸ This could be because our assumption for adsorption Na interactions, eqn (8), is oversimplified. The overall excellent agreement between our model based on specific material parameters and the experimental data,¹⁷Table 2, motivated us to explore the impact that the microstructures of SC and HC materials have on the change of the Na storage mechanism.

Fig. 6 (a) The OCV results when x = 0 obtained from experiment¹⁷ (blue) and from our model (red). (b) The OCV results when x = 2 obtained from experiment¹⁷ (blue) and from our model (red). (c) The OCV results when x = 4 obtained from experiment¹⁷ (blue) and from our model (red). (d) The OCV results when x = 6 obtained from experiment¹⁷ (blue) and from our model (red). (e) The OCV results when x = 10 obtained from experiment¹⁷ (blue) and from our model (red). (f) Comparison of the voltage error of the carbon materials obtained from experiment¹⁷ and from our model.

3.2.2 Parameterizing microstructure effects on the Na storage mechanism.
3.2.2.1 Adsorption Na.
a. The exponent number a₃. As discussed in Section 2.2.2, a higher value of a₃ indicates an increase of Adsorption Na at high voltage. We obtain a₃ = 0.5 as the best fit to the experimental data, Fig. 7. Therefore, the value of a₃ = 0.5 is fixed to 0.5 for all cases.

Fig. 7 The OCV results when x = 0.5 obtained from experiment¹⁷ (blue) and from our model (red) with a₃ = 0.5.

b. The maximum capacity of the adsorption Na. As shown in Fig. 8, the normalized maximum capacity for adsorption Na has a linear relationship with the specific area S_BET for S_BET ≤ 10 m² g⁻¹ (see Table 2). However, the specific area of the carbon material can fluctuate in the range from 0.1 m² g⁻¹ to 1000 m² g⁻¹.³⁵ A considerable extrapolation of the linear relation beyond 10 m² g⁻¹ cannot be valid. Deviations from the linear relation have to be caused by limitations in the accessibility of Na to the inner surface measured in BET experiments. A better understanding of the relation between the BET surface area and capacity due to adsorption is required to find proper generalizations of the observed linear relation.

Fig. 8 The maximum adsorption Na capacity of the carbon materials.

3.2.2.2 Layer Na.
a. The maximum capacity of the layer Na. DFT shows that Na cannot intercalate in between two graphite layers when the layer distance is too narrow, d₀₀₂ < 0.35 nm.²³ Therefore, there is no layer Na in the SC, x = 0.5. For the other carbon materials, the maximum layer Na capacity first increases up to x = 4 and then decreases (Fig. 9). This could be explained by a change in the microstructure of the carbon materials. When the carbon material changes from SC to HC, the layer distance d₀₀₂ increases (see Table 2). The larger the layer distance d₀₀₂, the more the layer Na intercalates between the two graphite layers. Therefore, the maximum layer Na capacity increases from x = 2 to 4. However, when the layer distance d₀₀₂ is large enough, the number of intercalated Na ions depends strongly on the layer length L_a and the number of graphite layers n (see Fig. 10). The layer Na increases when L_a and n increase. However, as can be seen in Table 2, the layer length L_a and the number of graphite layers n decreases from x = 0.5 to 10. Therefore, the maximum capacity of layer Na begins to decrease above x = 4 despite better accessibility to the layer sites.

Fig. 9 The maximum layer Na capacity Q^max₂ of the carbon materials.

Fig. 10 The dependence of layer Na on the layer length L_a and number of layers n at a long enough layer distance d₀₀₂. The orange circles are the extra layer Na when the layer length L_a (top) and number of layers n (bottom) increase.

b. The energy of interaction between the Na layer and graphite layers in the carbon materials, ε₂₀. In the DFT calculations,²³ the initial intercalation energy of layer Na ε₂₀ depends strongly on the layer distance d₀₀₂. Based on the DFT results, it diverges when the layer d₀₀₂ is close to the layer distance of graphite d_graphite. To mimic this behaviour, we approximate the value ε₂₀ by:(see Fig. 11) where d_graphite = 0.335 nm is the layer distance in graphite. This is also the minimum distance between 2 graphite layers required for intercalation of Na.³⁸a_l > 0 is the scaling energy. ε^∞_l is the insertion energy for layer Na when d₀₀₂ → ∞. ε₂₀ diverges when d₀₀₂ = d_graphite. Na cannot be intercalated at smaller values.

Fig. 11 The relationship between and ε₂₀ and .

c. The energy of interaction between intercalated layer Na, ε₂₁. To obtain a simple fit formula for the energy of interaction between intercalated layer Na, ε₂₁, we note a couple of properties:

• ε₂₁ should be a function of the total length (n − 1)L_a, Fig. 10.

• ε₂₁ is inversely proportional to the total length (n − 1)L_a. The longer the total length (n − 1)L_a, the more layer Na N₂ can intercalate in between two graphite layers.

• ε₂₁ should always be positive, which means the layer interaction ε₂ increase when the layer Na N₂ increases, eqn (9).

• ε₂₁(L_a → 0)→ + ∞, which means that at very short layer lengths, L_a → 0 no stable intercalation between two graphite layers is possible.

From the listed properties of ε₂₁ and based on the experimental data in ref. 17, we could fit a linear relation between ε₂₁ and the inverse total length 1/[(n − 1)L_a], Fig. 12.

Fig. 12 The relationship between ε₂₁ and 1/[(n − 1)L_a].

3.2.2.3 Pore Na. From the parameters of the carbon materials in Table 2 we can obtain a parabolic relation between the pore capacity and the pore radius r_p. At small r_p, Q^max₁ increases with r_p, however, at high r_p, the value of Q^max₁ decreases because Na starts to form a quasi-metallic state which prevents further insertion in the pores of the carbon materials (Fig. 13).⁴⁰ This also explains why Na does not fill in mesopores at very large r_p.⁴¹

Fig. 13 The relationship between the maximum pore Na capacity Q^max₁ and pore radius r_p.

3.2.3 The change of the Na storage mechanism on changing from soft carbon to hard carbon. With the results obtained above we can now summarize the nature of the transition from a soft carbon to a hard carbon material.
a. Soft carbon (x = 0.5). SC has narrowly spaced layers when compared to graphite, d₀₀₂ ≈ d_graphite, which means that Na cannot intercalate in between two graphite layers. Moreover, Na can also not insert into the pores because of the low porosity or high density ρ, respectively. Therefore, Na can only adsorb on the surface of the graphite layers, Fig. 14, facilitating fast insertion kinetics.

Fig. 14 The maximum pore capacity of each type of Na in the carbon materials.

b. Transition state (x = 2). The ‘transition state’ (TS) in between SC and HC is characterized by a larger pore volume and layer distance, allowing the onset of layer insertion and pore occupation (Fig. 14). The major insertion mechanism is still adsorption on the surface of the graphite layers. Due to a high specific area S_BET, (see Table 2), adsorption Na still contributes most to the total capacity.
c. Hard carbon (x = 4–10). Compared with the SC and TS, most Na in this state occupies the pores and does not adsorb on the surface of the graphite layers anymore. The reduction of the adsorption sites is reflected in the reduction of both the length of the graphite layers L_a and the specific area S_BET. Here, we can see the continuous change of the Na storage mechanism from SC (adsorption Na only), via the TS (mostly adsorption Na, less layer Na and pore Na) to HC (most pore Na, less layer Na, adsorption Na and defect Na), Fig. 14.

At low density, ρ ≈ 1.35 for x ≥ 6, the presence of defect Na also explains the small shoulder at φ = 1.8 V observed in the experimental OCV when x = 6 and x = 10. In general, the defect Na absorbs on the defect sites like heteroatoms and graphite defects.³⁷ The edge defects of non-parallel graphite increase with decreasing ρ. Therefore, it can be expected that Na ions adsorb on the edge defect of non-parallel graphite layers.

The change in the Na storage mechanism when changing from SC to HC is summarized again in Fig. 15.

Fig. 15 The change in the Na storage mechanism when changing from a soft carbon to a hard carbon material.

3.3 The site occupation at non-zero current

In the following chapter we discuss the impact of finite current on the site distribution of the Na ions. To do this we heuristically adjust the energy parameters in Table 1 to fit the voltage curves in the experiment. Under current conditions the capacity decreases as usual compared to OCV conditions. Part of the capacity decrease is due to overpotentials caused by the reaction kinetics of insertion and transport. In addition to causing overpotentials, transport limitations also induce a change in the relative occupation of the energetic sites within the carbon material, since there are different energy barriers for the different occupation sites.⁴² The change of shape of the voltage curve is therefore an indicator for the transport limitations between the energetically different sites. For future modelling of transport in disordered carbon materials it is therefore important to investigate the characteristic change in site occupation during operation. By analysing the shape of the voltage curve¹⁷ under different current loads we try to first obtain an indication how the relative occupation changes under current conditions.

3.3.1 Comparison with experiment. In ref. 16 the sample x = 6 was chosen to study the impact of different C-rates. The voltage curves obtained from experiments and our model by modifying the internal distribution of Na occupation are shown in Fig. 16. The excellent agreement with our model based on the specific material structure despite the neglection of overpotentials motivated us to explore the microstructure effect and the change of Na storage mechanism from SC to HC under load conditions. At high rates, the agreement decreases due to the kinetic effect i.e. non-equilibrium properties. This will be explored in future works.

Fig. 16 (a) The OCV results obtained from experiment¹⁷ (blue) and our model (red) at a scan rate of 0.025 A g⁻¹ when x = 6. (b) The OCV results obtained from experiment¹⁷ (blue) and our model (red) at a scan rate of 0.1 A g⁻¹ when x = 6. (c) The OCV results obtained from experiment¹⁷ (blue) and our model (red) at a scan rate of 0.2 A g⁻¹ when x = 6. (d) The OCV results obtained from experiment¹⁷ (blue) and our model (red) at a scan rate of 0.5 A g⁻¹ when x = 6. (e) The OCV results obtained from experiment¹⁷ (blue) and our model (red) at a scan rate of 1 A g⁻¹ when x = 6. (f) The error between experiment¹⁷ and our model when x = 6 at different scan rates.

3.3.2 The impact of finite currents on the Na storage mechanism. In Fig. 17 it can be seen that layer Na quickly approaches a limiting capacity while other types decrease with increasing C-rate. Specifically, the capacity of pore Na decreases considerably with increasing current. This could be explained by a strong limitation of transport between pores due to a high tortuosity in the random microstructure which prevents Na filling into the pores at very high C-rates. The only efficient insertion mechanism at large specific currents is adsorption. Intercalation and transport within the layers still seem to be possible even at large currents. The constant capacity contribution of the layer sites, especially at high currents, is most probably an artefact due to the neglection of intercalation energy barriers. We will investigate the impact of intercalation overpotentials in a future paper.

Fig. 17 The maximum capacity of each type of Na at different scan rates when x = 6.

4 Conclusions

In this paper, we modify our previous phenomenological model and develop a unified mean field theory for the OCV of carbon materials ranging from soft carbon to hard carbon materials. Our model OCV as well as maximum capacity values are in good agreement with experiment. The change in the occupation of different energetic sites for Na storage in the SC and HC carbon materials allows us not only to illustrate the effect of microstructure on the Na storage mechanism but also to clarify the change in Na storage that occurs upon changing from SC to HC.

To the best of our knowledge, for the first time, a unified fit for the OCV of disordered carbons including a simple parameterization of different microstructures is obtained. The change in the OCV with changing microstructure is interpreted as a result of the changing distribution of the Na ions across physically different occupation sites. The maximum capacity of the adsorption Na follows a linear relationship with the specific area when S_BET ≤ 10 m² g⁻¹. The initial intercalation energy of layer Na depends on the layer distance d₀₀₂, while the interaction energy of Na inserted within the layers depends on the total length (n − 1)L_a of the graphitic domains. The maximum pore Na capacity increases with pore radius r_p at small values of r_p. However, for large values of r_p, the maximum pore Na capacity decreases due to quasi-metallic Na aggregation, which blocks the pores against the further insertion of Na.

The change of Na storage mechanism on changing from SC to HC is also demonstrated. SC has a narrow layer distance and no porosity; therefore, Na can only adsorb on the surface of the graphite layers. Upon changing to HC the amount of Na inserting into the graphitic layers and pores increases due to the increasing layer distance in HC and more accessible pores. At lower density, Na can also adsorb on the edge defect sites of HC.

The impact of electric currents on the Na storage mechanism was also studied. By modifying the internal distribution of Na occupation, quasi-OCV curves were fitted to the voltage curves at different currents. Even in the ideal situation where overpotentials are neglected, a shift in the site distribution towards adsorption becomes obvious. This observation is consistent with the observed high rate capability of hard carbon, which is facilitated by the adsorption of Na ions on internal graphitic surfaces.

In this work, the purpose is to obtain an analytic model for the OCV, which can be easily parameterized by DFT and/or experiment. In future publications it is planned to use this OCV model within a non-equilibrium framework to study the behaviour of the anode under current load conditions, as e.g., in the theory of strongly interacting electrolytes by Schammer et al.⁴⁵ This framework is capable of including all mentioned dynamic effects via appropriate modelling of the Onsager coefficients. Moreover, the OCV model will be integrated into full cell models, including kinetic properties, for sodium ion batteries with carbon material anodes.

Author contributions

Huy Sy Nguyen contributed to the conceptualization, data curation, formal analysis, investigation, methodology, resources, writing – original draft, and writing – reviewing and editing. Arnulf Latz contributed to project administration, supervision, formal analysis, validation, visualization, and writing – review and editing.

Conflicts of interest

There are no conflicts to declare.

Data availability

The dataset used to prepare this manuscript can be found at Zenodo (https://doi.org/10.5281/zenodo.15711067).

Appendices

Appendix A: Ideal chemical potential of Na in carbon materials

The free energy with layers that have different maximal occupation:Consider only the ideal part and use Stirling's formula:⁴³

ln N! = N ln N − N + o(ln N). (A2)

Therefore we have:The ideal chemical potential of each type of Na is:where: and .

The ideal contribution to the chemical potential of the anode is:

Appendix B: Excess chemical potential of Na in carbon materials

Considering an NVT system for Na insertion in HC as shown in Fig. 18.

Fig. 18 NVT for Na insertion in hard carbon. The blue, red, black, and orange areas indicate available sites in the pores, within graphite layers, on the surface of graphitic sheets, and on defects, respectively.

For an NVT system, the excess chemical potential is calculated via Widom's method as:⁴⁴

where U_j is the change of energy in the system for each trial j in n trials of particle insertion.

In mean field theory, the change of energy in each site i is equal to the excess chemical potential of this type:

U_j = μ^ex_i, j ∈ V_i. (B3)

where V_i, i = 1–4 indicate the pore, layer, adsorption, and defect sites, respectively.

For a given distribution of Na ions across the different sites eqn (B2) may be approximated as

where ξ_i is the probability of a trial particle j inserting at the same site in V_i with .

In mean field theory, the insertion probability ξ_i is approximated as the probability of choosing a site in V_i. This probability is given by the number fraction of site i with a maximal number of sites N^max_i

Substituting eqn (B4) and (B5) into eqn (B2), we have:

Acknowledgements

This work contributes to research performed at CELEST (Center for Electrochemical Energy Storage Ulm – Karlsruhe) and was funded by the German Research Foundation (DFG) under project ID 390874152 (POLiS Cluster of Excellence). This work was also supported by computer resources from the University of Ulm and German Aerospace Center (DLR).

References

1. I. Moeez, B. Ali, K. W. Nam and K. Y. Chung, Chem. Eng. J., 2025, 503, 158150,  DOI:10.1016/j.cej.2024.158150.
2. D. Xia, W. Huang, C. Shi, A. Promi, D. Huo, C. Sun, S. Hwang, G. Kwon, H. Huang and F. Lin, ACS Sustainable Chem. Eng., 2025, 13(1), 241–250,  DOI:10.1021/acssuschemeng.4c06900.
3. L. He, S. Wei, X. Zhang, S. Wang, Y. Xia, Z. Ni, C. Li, W. Dong, D. Shen and S. Yang, J. Energy Storage, 2025, 111, 115426.
4. K. Xiong, T. Qi and X. Zhang, Electroanalysis, 2025, 37, e202400318.
5. R. Karunarathna, H. R. Arachchige, S. Karunarathne, W. P. S. L. Wijesinghe, C. Sandaruwan, M. M. M. P. G. Mantilaka, Y. Y. Kannangara and A. M. Abdelkader, Small Sci., 2025, 5, 2400405,  DOI:10.1002/smsc.202400405.
6. X. Feng, F. Wu, Y. Fu, Y. Li, Y. Gong, X. Ma, P. Zhang, C. Wu and Y. Bai, Small, 2025, 21, 2409120,  DOI:10.1002/smll.202409120.
7. A. V. Desai, V. R. Seymour, R. Ettlinger, A. Pramanik, A. G. Manche, D. N. Rainer, P. S. Wheatley, J. M. Griffin, R. E. Morris and A. R. Amstrong, Chem. Commun., 2023, 59, 1321–1324.
8. K. Huang, S. Tian, H. Xu, C. Fang, L. Wu, H. Wang, F. Liu, W. He and X. Zhang, Chem. Eng. J., 2025, 505, 159326.
9. Y. Jin, H. Yoo, H. Seong, J. H. Moon, G. Kim, T. Jung, Y. Myung, W. Lee, S. Kim and J. Choi, Chem. Eng. J., 2025, 503, 158146,  DOI:10.1016/j.cej.2024.158146.
10. Z. Lin, J. Wu, Z. Wu, M. Yu, H. Jia, K. Zhou, X. Huang and S. Ying, J. Energy Storage, 2025, 106, 114924,  DOI:10.1016/j.est.2024.114924.
11. Z. Zhan, R. Li, Y. Mei, D. Liu, C. Zheng and S. Huang, J. Alloys Compd., 2025, 1017, 178960,  DOI:10.1016/j.jallcom.2025.178960.
12. M. M. Almarzoge, M. Gencten and G. Ozsin, ChemElectroChem, 2025, 12, e202400564,  DOI:10.1002/celc.202400564.
13. J. Gu, Y. You, M. Liu, L. Huang, Z. Sun, J. Liu, L. Sha, M. Chen, S. Li, S. Wu, Q. Zhang and L. Zhang, Chem. Eng. J., 2025, 505, 159331,  DOI:10.1016/j.cej.2025.159331.
14. T. Qi, K. Xiong and X. Zhang, J. Power Sources, 2025, 626, 235721,  DOI:10.1016/j.jpowsour.2024.235721.
15. Y. Liu, J. Yin, R. Wu, H. Zhang, R. Zhang, R. Huo, J. Zhao, K. Y. Zhang, J. Yin, X. L. Wu and H. Zhu, Energy Storage Mater., 2025, 75, 104008,  DOI:10.1016/j.ensm.2025.104008.
16. C. Jiang, Y. Li, C. Wan, Q. Wu, S. Yang, K. Chen, Y. You, Y. Tian, J. Huang and H. Xie, Chem. Eng. J., 2025, 506, 160124,  DOI:10.1016/j.cej.2025.160124.
17. D. Cheng, X. Zhou, H. Hu, Z. Li, J. Chen, L. Miao, X. Ye and H. Zhang, Carbon, 2021, 182, 758–769.
18. H. D. Asfaw, C. W. Tai, M. Valvo and R. Younesi, Mater. Today Energy, 2020, 18, 100505.
19. Y. Yin, Y. Tan, Y. Lu, Y. Wang, J. Yang, Y. Li and B. Huang, J. Energy Storage, 2025, 113, 115649,  DOI:10.1016/j.est.2025.115649.
20. D. Xiang, W. Zhu, Z. Sun, L. Li, K. Sun, H. Xu, Y. Chen and H. Yang, Diamond Related Mater., 2025, 152, 111894,  DOI:10.1016/j.diamond.2024.111894.
21. X. Yao, Y. Ke, W. Ren, X. Wang, F. Xiong, W. Yang, M. Qin, Q. Li and L. Mai, Adv. Energy Mater., 2018, 8, 1803260,  DOI:10.1002/aenm.201900094.
22. F. Xie, Y. Niu, Q. Zhang, Z. Guo, Z. Hu, Q. Zhou, Z. Xu, Y. Li, R. Yan, Y. Lu, M. M. Titirici and Y. S. Hu, Angew. Chem., Int. Ed., 2022, 61, e202116394,  DOI:10.1002/anie.202116394.
23. C. Cai, Y. Chen, P. Hu, T. Zhu, X. Li, Q. Yu, L. Zhou, X. Yang and L. Mai, Small, 2022, 18, 2105303.
24. S. Alvin, C. Chandra and J. Kim, Chem. Eng. J., 2021, 411, 128490,  DOI:10.1016/j.cej.2021.128490.
25. N. Sun, Z. Guan, Y. Liu, Y. Cao, Q. Zhu, H. Liu, Z. Wang, P. Zhang and B. Xu, Adv. Energy Mater., 2019, 9, 1901351,  DOI:10.1002/aenm.201901351.
26. M. P. Mercer, S. Affleck, E. M. G. Arriazu, A. A. Zülke, P. A. Maughan, S. Trivedi, M. Fichtner, A. R. Munnangi, E. P. M. Leiva and H. E. Hoster, ChemPhysChem, 2022, 23, e202100748,  DOI:10.1002/cphc.202100748.
27. Y. Jin, S. Sun, M. Ou, Y. Liu, C. Fan, X. Sun, J. Peng, Y. Li, Y. Qiu, P. Wei, Z. Deng, Y. Xu, J. Han and Y. Huang, ACS Appl. Energy Mater., 2018, 1(5), 2295–2305.
28. X. X. He, J. H. Zhao, W. H. Lai, R. Li, Z. Yang, C. Xu, Y. Dai, Y. Gao, X. H. Liu, L. Li, G. Xu, Y. Qiao, S. L. Chou and M. Wu, ACS Appl. Mater. Interfaces, 2021, 13, 44358–44368,  DOI:10.1021/acsami.1c12171.
29. G. Huang, H. Zhang, F. Gao, D. Zhang, Z. Zhang, Y. Liu, Z. Shang, C. Gao, L. Luo, M. Terrones and Y. Wang, Carbon, 2024, 228, 119354,  DOI:10.1016/j.carbon.2024.119354.
30. K. C. Wasalathilake, G. A. Ayoko and C. Yan, Carbon, 2018, 140, 276–285,  DOI:10.1016/j.carbon.2018.08.071.
31. E. Olsson, J. Cottom, H. Au, Z. Guo, A. C. S. Jensen, H. Alptekin, A. J. Drew, M. M. Titirici and Q. Cai, Adv. Funct. Mater., 2020, 30, 1908209,  DOI:10.1002/adfm.201908209.
32. M. A. Reddy, M. Helen, A. Groβ, M. Fichtner and H. Euchner, ACS Energy Lett., 2018, 3, 2851–2857,  DOI:10.1021/acsenergylett.8b01761.
33. H. Euchner, B. P. Vinayan, M. A. Reddy, M. Fichtner and A. Groβ, J. Mater. Chem. A, 2020, 8, 14205–14213,  10.1039/d0ta04868b.
34. T. Ji, X. Liu, D. Sheng, Y. Li, H. Ruan, H. Guo, Z. X. Shen and L. Lai, Energy Storage Mater., 2024, 71, 103563,  DOI:10.1016/j.ensm.2024.103563.
35. X. Liu, T. Wang, T. Ji, H. Wang, H. Liu, J. Li and D. Chao, J. Mater. Chem. A, 2022, 10, 8031–8046,  10.1039/d1ta10588d.
36. B. Oral, B. Tekin, D. Eroglu and R. Yildirim, J. Power Sources, 2022, 549, 232126.
37. J. Li, Y. Jin, K. Sun, A. Wang, G. Zhang, L. Zhou, W. Yang, M. Fan, J. Jiang, Y. Wen and S. Wang, J. Power Sources, 2025, 631, 236323,  DOI:10.1016/j.jpowsour.2025.236323.
38. H. S. Nguyen and A. Latz, Phys. Chem. Chem. Phys., 2023, 25(41), 28196–28204,  10.1039/D3CP03186A.
39. H. Alptekin, H. Au, E. Olsson, J. Cottom, A. C. S. Jensen, T. F. Headen, Q. Cai, A. J. Drew, M. C. Ribadeneyra and M. M. Titirici, Adv. Mater. Interfaces, 2022, 9, 2101267,  DOI:10.1002/admi.202101276.
40. H. Au, H. Alptekin, A. C. S. Jensen, E. Olsson, C. A. O’Keefe, T. Smith, M. C. Ribadeneyra, T. F. Headen, C. P. Grey, Q. Cai, A. J. Drew and M. M. Titirici, Energy Environ. Sci., 2020, 13, 3469–3479,  10.1039/d0ee01363c.
41. M. Yuan, B. Cao, H. Liu, C. Meng, J. Wu, S. Zhang, A. Li, X. Chen and H. Song, Chem. Mater., 2022, 34, 3489–3500,  DOI:10.1021/acs.chemmater.2c00405.
42. A. C. S. Jensen, E. Olsson, H. Au, H. Alptekin, Z. Yang, S. Cottrell, K. Yokoyama, Q. Cai, M. M. Titirici and A. J. Drew, J. Mater. Chem. A, 2020, 8, 743–749,  10.1039/C9TA10113F.
43. G. J. O. Jameson, Mathematical Gazette, 2015, 99(544), 68–74,  DOI:10.1017/mag.2014.9.
44. D. M. Heyes, Chem. Phys., 1992, 159, 149–167.
45. M. Schammer, B. Horstmann and A. Latz, J. Electrochem. Soc., 2021, 168(2), 026511,  DOI:10.1149/1945-7111/abdddf.

---