Investigation of van der Waals interactions in two-dimensional MXenes via first-principles simulations

Abstract

This study systematically examines van der Waals (vdW) interactions in MXenes using density functional theory (DFT) with dispersion-corrected methods. We evaluate Grimme D2, Grimme D3, the exchange-hole dipole moment (XDM) method, and the nonlocal vdW-DF2 functional, pointing out their impact on interlayer distance, lattice parameters, and in-plane elastic stiffness for ten M₂CX_T compositions (M = Sc, Ti, V, Cr, Zr, Nb, Mo, Hf, Ta, and W and X = O, F, H, and OH). Interlayer spacings range from 4.8 to 9.7 Å, lattice constants range from 2.9 to 3.3 Å, and corrected in-plane stiffness values range from ∼100 to 360 N m⁻¹. Oxygen terminations tend to produce the highest stiffness, while hydrogen terminations provide denser lattices and relatively strong stiffness, making them promising for both mechanically resilient and flexible devices. Fluorine and hydroxyl groups provide increased spacings and softer responses, offering additional tunability. Metal-dependent trends describe bond contraction across a period and relativistic effects down a group. Benchmark comparisons suggest that D3 and XDM most closely reproduce experimental ranges (e.g., Ti₃C₂T_x XRD and nanoindentation), but D2 overbinds and vdW-DF2 underbinds. No correction is superior in all cases; method choice must be guided by the properties and chemistry of interest. This reproducible dataset establishes the structure–property tendencies between metals, terminations, and vdW schemes, providing design rules for the engineering of MXenes in flexible electronics, energy storage, and catalysis.

---

1. Introduction

Since the experimental isolation of graphene in 2004,¹ interest in studying these new classes of two-dimensional van der Waals (vdW) materials has grown rapidly. Because of their dimensional reduction, these 2D layered materials show properties absent in bulk 3D materials, with the known advantages of nanomaterials largely due to quantum confinement and enhanced surface effects. Moreover, the binding of these materials occurs through weak van der Waals interactions, which allows one to stack different materials of different properties to form a unique heterostructure with valuable hybrid properties in a small volume.² Based on their structural and electronic properties, as well as tunable surface terminations, two-dimensional MXenes have emerged as extremely promising candidates for various applications, including energy storage, membrane separation, catalysis, and nanoelectronics.^3–8 A key factor supporting this versatility is the interlayer interaction, which is dominated by van der Waals (vdW) forces. These forces control stacking, stability, and processability and directly govern the properties of practical interest. Interlayer binding controls the ion intercalation capacity and kinetics in batteries and supercapacitors,^9,10 determines the interlayer distance that allows size- and charge-selective sieving in membranes,^11,12 and manages exfoliation and restacking that define accessible surface area and long-term stability. Experimental studies have always revealed that sub-nanometer spacing changes can significantly modify charge storage and ion transport.^4,13 Simultaneously, chemical functionalization and intercalation introduce competing interactions—electrostatic, hydrogen bonding, and solvation—which make interpretation difficult.

The vdW-tuned interlayer distance and bonds have direct consequences for the mechanical strength and flexibility of MXene composites and films and have consequences for their use in flexible electronics, sensors, and structural reinforcements. Since the origin of these phenomena is a fine balance between short-range bonding and long-range dispersion, predictions by computation are extremely sensitive to the selected vdW correction. It is thus important to benchmark dispersion schemes (e.g., D2, D3, XDM, and vdW-DF2) against experiments for a reliable approach to simulate MXenes. Although experimental reference data are currently limited to a subset of MXenes (e.g., Ti₃C₂T_x and Nb₂C), these measurements provide useful external checks for computational predictions. Interestingly, theoretical and computational simulations initially predicted most of these materials before isolating them in the laboratory.^14–20 These visuals highlight the pressing need for reliable, benchmarked computational tools to forecast how the composition, surface terminations, and stacking order interact to tune the MXene interlayer's properties and facilitate practical device optimization. Additionally, the extensive range of possibilities of these materials makes affordable and accurate computational simulations a vital step during their development, alongside conventional experiments.

Currently, the first principle of simulation is indispensable for understanding, designing, and developing MXene-based materials. Advances in density functional theory (DFT)^13–15 and its implementation in robust computational packages^{9–11,21–23} have made the study of real-life materials practically feasible as routine procedures with increasing accuracy and efficiency. Specifically, for 2D materials, the intralayer is typically firmly bound via covalent and ionic forces, which are well modeled by semi-local DFT (e.g., functionals like PBE). Thus, the essential energetic characteristics of a 2D layered material relate to the weak vdW-related forces between layers that are ultimately responsible for stability, and these long-range interactions are not adequately described by semi-local DFT. In this study, we focus on systematically evaluating various dispersion correction methods within the DFT framework to better understand the role of vdW forces in 2D-layered MXenes across different terminations and explore periodic and group-dependent trends.

van der Waals forces arise from spontaneous dipole–induced dipole interactions; the origin of this interaction is attributed to zero-point quantum and thermal fluctuations. Although individually weak and scaled as 1/R⁶, their cumulative effect can significantly influence key properties, such as the interlayer lattice constant, elastic constant (C₃₃), and binding energy (U_b). For example, the elastic modulus of graphite is around 1000 GPa in the intralayer but around 40 GPa in the interlayer, underscoring the importance of accurately capturing these weak interactions.^24,25

Despite significant progress, a few fundamental challenges remain in DFT, particularly in accurately predicting band gaps, vdW interactions, and excited states.^26–28 As the material of interest in the current study is 2D MXene, where the interlayer binding is due to vdW forces, which are not accurately described by DFT, this study focuses on methods that account for these vdW interactions. In recent decades, numerous dispersion correction methods and vdW functionals have been developed.^29–39 Benchmark studies on standard datasets, e.g., S22⁴⁰ (consists of 22 noncovalent molecular complexes), S66⁴¹ (consists of 66 weakly interacting molecular complexes with biochemical relevance), and S66x8 (extended version of S66)⁴² have evaluated several approaches, including D2, D3, MBD, XDM, vdW(TS), and D3-ATM, in conjunction with the PBE–GGA exchange–correlation functional.⁴³

Based on previous work on different dispersion correction methods^29,33,36,43 and mindful of computational efficiency, this work selects the D2, D3, XDM, and vdW-DF2 methods to systematically assess their effects on the electronic, structural, and mechanical properties of MXenes. The study also includes these results with those obtained without any dispersion corrections. Our findings indicate that no single dispersion correction method uniformly outperforms the others across all investigated properties. In many cases, simple pairwise dispersion correction methods, such as D2 and D3, have shown better performance than high-level sophisticated non-local exchange correlational (E_XC) functionals, such as vdW-DF2. This observation emphasizes that the choice of dispersion correction should be tailored to the specific material system and properties under investigation, thus necessitating a systematic evaluation for accurate predictive modeling.

2. Computational section

DFT calculations for MXenes are performed as implemented in the Quantum Espresso package²² to systematically investigate the structural and electronic properties of MXenes. In our calculations, the generalized gradient approximation (GGA) was employed with the Perdew–Burke–Ernzerhof (PBE) exchange–correlation functional.⁴⁴ The projected augmented wave (PAW) method was used in conjunction with a plane-wave basis set. A kinetic energy cut-off of 75 Ry and an energy convergence threshold of 10⁻⁷ Ry ensured the numerical accuracy of the calculations. Integration over the Brillouin zone was performed using a 14 × 14 × 1 Monkhorst–Pack k-point mesh established through extensive convergence tests.

To avoid interactions along the out-of-plane direction, the c-axis was extended to 35 Å, while the a and b lattice parameters were maintained as in the original unit cell. Structural optimizations were performed until the forces on all atoms were reduced below 1 × 10⁻⁴ Ry Å⁻¹. For model construction and visualization, Vesta software is utilized.⁴⁵ To calculate the elastic constant values, we employed the Thermopw package.⁴⁶ We froze both the geometry (i.e. unit cell dimensions) and lattice points. The elastic constants were obtained in kbar and converted to two-dimensional units (N m⁻¹) following standard practice for 2D materials.

Two stacking configurations were considered when building the MXene models. With terminal groups—oxygen (O), hydrogen (H), hydroxyl (OH), or fluorine (F)—positioned at the FCC hollow site on the surface, single-layer MXenes were constructed utilizing an ABC stacking pattern. To precisely capture the interlayer interactions for bi-layer systems, an octahedral stacking configuration was chosen. Every model constructed is an M₂CX_T configuration, where X is the previously described termination and M is the early transition metal. M₄C₂X_T denotes the bilayer models used in the current work.

The computational models employed herein closely correspond to the MXene systems obtainable via experimentation. These models consist of mono- and bilayer architectures with varying terminations. Careful top-down etching of MAX phase precursors can produce few-layer and even bilayer MXenes. This process uses hydrofluoric acid (HF), LiF/HCl mixtures, or molten salt methods. Several studies using AFM, TEM, and XRD have recently shown that intercalating agents, such as tetramethylammonium hydroxide (TMAOH), dimethyl sulfoxide (DMSO), and various ionic liquids, enable successful delamination into mono- and few-layer sheets upon synthesis.^3,6 X-Ray diffraction (XRD) experiments have explicitly verified computational predictions for pristine and intercalated MXenes. These experiments recorded interlayer separations from about 1.0 to 3.2 nm.^3,12 The synthesis and characterization of flexible films and suspended bilayer MXene membranes for mechanical and electrical analyses highlight the experimental relevance of the stacking and termination configurations under study. The findings show that modern synthesis and processing technologies directly enable these idealized structures, illustrating that they are no longer theoretical. Therefore, the mechanical, electrical, and structural understanding of our models provide guides for experimental optimization and tweaking. We also emphasize that the bilayer structures used here are not abstract but reflect experimentally accessible few-layer stacks of MXenes. Selective etching of MAX phases routinely produces multilayer aggregates of Ti₂C, V₂C, and Mo₂C.^5,47 Our bilayer models thus form a realistic description of the interlayer environments in synthesized MXenes that allow for the systematic testing of vdW corrections in a physically relevant framework.

The choice of the FCC hollow site for the terminal groups was made deliberately to focus the study on the effects of van der Waals (vdW) interactions. Previous studies have demonstrated that the FCC hollow site generally offers favourable adsorption energies and structural stability. As in close-packed metal surfaces, the term “FCC hollow site” refers here to the threefold symmetric hollow position centered above three neighboring surface metal atoms. In our MXene slabs, terminations such as O, F, H, and OH preferentially occupy these FCC hollow positions, which ensures maximum coordination with the underlying transition-metal layer and stabilizes the surface geometry.⁴⁸ By fixing the terminal group adsorption at this well-established site, we can isolate and assess the influence of different DFT dispersion correction methods—D2, D3, XDM, and vdW-DF2—on the electronic structure, interlayer stability, and mechanical properties of MXenes. Conducting an exhaustive search for the global minimum energy configuration across all possible adsorption sites would significantly increase computational costs and introduce additional complexity. This could obscure the specific influence of van der Waals (vdW) interactions, which is the primary focus of this study.

A systematic series of MXene models was generated by combining ten transition metals (Sc, Ti, V, Cr, Zr, Nb, Mo, Hf, Ta, and W) with the four aforementioned terminal groups. This diverse set allowed for a comprehensive assessment of how both metal identity and surface termination influence the intrinsic properties of MXenes, providing robust insights relevant to their application in hydrogen evolution, oxygen reduction, CO₂ reduction, electronics, water splitting, water treatment, desalination, electromagnetic interference (EMI) shielding and flexible sensing devices (Fig. 1).

Fig. 1 Model representing the Ti₂CF₂ structure: (a) monolayer single unit cell (side view), (b) monolayer supercell (2 × 2, top view), and (c) bilayer supercell (2 × 2, top view). The Ti, C and F atoms are represented by light sky blue, mahogany brown, and pale periwinkle, respectively.

An accurate description of dispersion interactions in MXenes is crucial due to the weak van der Waals (vdW) forces present between the layers. These interactions play a significant role in determining the structural and electronic properties, and a reliable computational method is necessary to model them accurately. This section details the dispersion correction methods employed: D2, D3, XDM, and vdW-DF2. Each method is chosen owing to its strengths in accounting for the van der Waals forces pertinent to MXenes.

2.1 Grimme's D2 dispersion correction

This method was proposed by Grimme and is a widely used semi-empirical approach that supplements standard density functional theory (DFT) calculations with dispersion correction. This is done by introducing a pairwise additive term to the DFT energy to account for long-range dispersion forces, using predefined C₆ coefficients for each atom pair.³¹
where C₆,_ij: dispersion coefficient for atom pair i and j; R_ij: distance between atoms i and j; S₆: scaling factor dependent on functional; f_damp (R_ij): damping function that prevents singularities at small spacing.

Where R_ij is the distance between atoms i and j, and f_damp is a damping function to correct short-range behaviour. Grimme's D2 is simple and computationally efficient, making it suitable for large-scale systems, like MXenes, known for their outstanding chemical tunability and conductive properties. However, its reliance on fixed coefficients limits its flexibility in complex environments, potentially leading to inaccuracies in systems with diverse chemical characteristics.

2.2 Grimme's D3 dispersion correction

Improving on the D2 computational approach, the D3 correction introduces environment-dependent dispersion coefficients and a three-body interaction term to improve accuracy. This approach adjusts its coefficients based on the local coordination number, hence making it more refined and accurate.⁴⁹
where C₆,_ij and C_8,ij: dispersion coefficients; R_ij: atomic spacing; f_d,6(R_ij) and f_d,8(R_ij): damping functions; optional 3-body term corrections can be included in D3.

The inclusion of a three-body term allows D3 to capture many-body dispersion interactions that are critical for accurate modeling in layered systems, like MXenes. This method provides a balance between computational cost and accuracy, making it a reliable choice where pairwise interactions alone may be insufficient.

2.3 Exchange-hole dipole moment (XDM) method

XDM is a non-empirical approach that models dispersion forces based on electron density distribution. Unlike the semi-empirical D2 and D3 computational approaches, XDM calculates dispersion energies using properties derived from the electron density, specifically the exchange-hole dipole moment.³²
where f(R_ij): damping function; C_6,ij, C_8,ij and C₁₀: dispersion coefficients derived from dipole moments.

This computational approach does not rely on predefined parameters, making it a physically grounded technique for capturing dispersion in complex systems. Although XDM is computationally more demanding, it provides a high level of accuracy by directly linking dispersion corrections to the electronic structure of the system. This density-based computational approach is advantageous because non-covalent interactions are sensitive to the local environment.

2.4 Second generation van der Waals density functional (vdW-DF2)

The vdW-DF2 computational approach is part of the van der Waals density functional family, incorporating non-local correlation terms directly into the functional form. Unlike D2 and D3, vdW-DF2 does not require additional correction terms, as the functional is designed to account for van der Waals forces.³⁰

The non-local correlation term E_nl is added to the general gradient approximation (GGA) energy, capturing long-range dispersion interactions inherently:

E_DFT-vdW = E_GGA + E_nl,

where n(r) and n(r′): electron densities at points r and r′, respectively. ϕ(r,r′): Kernel function that depends on the separation distance and the electron density gradient.

This makes vdW-DF2 particularly well-suited for layered materials, like MXenes, known for their outstanding chemical tunability and conductive properties, where long-range interactions significantly influence structural and electronic properties. However, the computational cost associated with the evaluation of non-local terms is higher, which may limit its application to large-scale simulations.

Depending on the computational implementation and pseudopotential choices, the absolute numerical values of certain properties may vary slightly, but the overall trends reported in this study are robust and physically meaningful. Finally, to assess the accuracy of the dispersion correction schemes, we qualitatively compared our computed structural and elastic properties with the reference data available in the literature. For example, our values of interlayer distances fall within the ranges reported for Ti₃C₂T_x MXenes with XRD measurements,^12,47 while in-plane elastic constants were compared with nanoindentation and tensile measurements.^12,50 In addition, previous DFT studies employing semi-empirical (D2 and D3) and nonlocal (vdW-DF and vdW-DF2) schemes^51,52 were used as theoretical reference points. We emphasize that direct experimental benchmarks are not yet available for the specific MXenes studied here; however, these comparisons provide a reasonable external check and a framework for the error analysis discussed in the results section.

3. Results and discussions

The prediction of the structural and mechanical properties of MXenes relies critically on the proper treatment of the dispersion forces. These interactions are essential for determining the interlayer separation and overall mechanical stability in these layered systems. As mentioned earlier, the standard DFT approaches, such as those employing the PBE functional, often fail to capture these weak interactions adequately, leading to significant deviations in the calculated properties when dispersion corrections are neglected.

Studies on 2D materials have shown that while dispersion-corrected methods like PBE + D3(BJ) can offer reliable predictions for lattice constants, they may substantially overestimate binding energies, with reported mean absolute errors (MAEs) as high as 59%. Similarly, approaches generally considered more accurate, such as PBE + D3(0), have exhibited MAEs of around 40% in binding energy predictions.⁴³ These discrepancies highlight the risk of relying on a single property, such as lattice constants, for validating computational methods, as good agreement in one parameter does not guarantee overall predictive accuracy.

Consequently, the present study systematically evaluates multiple key properties—including interlayer spacing, mechanical parameters and lattice constants—both with and without dispersion corrections. This comprehensive approach is aimed at identifying the most suitable dispersion-corrected DFT methodology for modeling MXenes. By doing so, we ensure that the computational predictions are consistently reliable across different material properties, thereby providing a clearer insight into the overall performance and limitations of each dispersion correction method. Few missing values correspond to the structurally unstable models in the considered configurations.

3.1 Interlayer spacing

In this work, we calculated the inter-layer distance. To avoid interaction in the c direction, we extended the c axis to 35 Å and kept the original unit cell dimensions in a and b. These calculations are performed for M₂C (M = Sc, Ti, V, Cr, Zr, Nb, Mo, Hf, Ta, and W) MXenes with oxygen, fluorine, hydrogen and hydroxyl termination groups.

The interlayer spacing is a critical parameter because it influences both the electronic and mechanical behaviours of these materials, and in layered systems, they are predominantly bound by van der Waals forces. Standard DFT methods, such as those employing the PBE functional without dispersion corrections, tend to overestimate these spacings because they do not adequately account for weak, long-range interactions. In contrast, the application of dispersion corrections—D2, D3, XDM, and vdW-DF2—systematically reduces the calculated interlayer spacing, bringing them closer to the experimental and other theoretical values. Tables 1–4 summarizes the interlayer spacing for the various terminations. In the present work, interlayer spacing is defined by the distance between two carbon atoms of alternate layers.

Table 1 Interlayer spacing calculated using different dispersion correction methods for oxygen-terminated MXenes

System	D3 (Å)	D2 (Å)	XDM (Å)	vdW-DF2 (Å)	Without vdW (Å)	Mean standard deviation
Sc4C2O4	6.75	6.41	6.38	6.80	6.66	6.60 ± 0.14
Ti4C2O4	6.95	6.88	6.75	7.11	7.78	7.10 ± 0.29
V4C2O4	6.98	6.92	6.91	7.06	7.79	7.19 ± 0.31
Cr4C2O4	7.14	7.10	7.00	7.15	8.38	7.35 ± 0.12
Zr4C2O4	6.22	6.49	6.07	6.92	6.50	6.32 ± 0.16
Nb4C2O4	6.99	6.87	6.85	7.16	7.81	7.13 ± 0.44
Mo4C2O4	6.87	6.73	6.76	6.99	7.67	7.01 ± 0.27
Ta4C2O4	7.04	6.63	6.85	7.35	7.78	7.07 ± 0.38
W4C2O4	6.88	6.88	6.75	7.29	7.66	7.04 ± 0.31

---

Table 2 Interlayer spacing calculated using different dispersion correction methods for fluorine-terminated MXenes

System	D3 (Å)	D2 (Å)	XDM (Å)	vdW-DF2 (Å)	Without vdW (Å)	Mean standard deviation
Sc4C2F4	6.88	6.80	6.65	6.90	7.43	6.93 ± 0.21
Ti4C2F4	7.10	7.05	7.00	7.13	7.85	7.23 ± 0.25
V4C2F4	6.98	6.89	6.90	6.97	7.75	7.10 ± 0.26
Cr4C2F4	6.75	6.59	6.61	6.64	7.48	6.81 ± 0.08
Zr4C2F4	7.16	7.02	7.02	7.16	7.81	7.23 ± 0.23
Nb4C2F4	7.06	6.92	7.02	7.04	7.80	7.17 ± 0.24
Hf4C2F4	7.22	6.61	6.98	7.33	7.69	7.13 ± 0.34
Ta4C2F4	7.10	7.77	7.07	7.00	7.77	7.34 ± 0.25
W4C2F4	6.62	6.23	6.55	6.68	7.31	6.68 ± 0.28

---

Table 3 Interlayer spacing calculated using different dispersion correction methods for hydrogen-terminated MXenes

System	D3 (Å)	D2 (Å)	XDM (Å)	vdW-DF2 (Å)	Without vdW (Å)	Mean standard deviation
Ti4C2H4	5.07	5.08	5.04	5.22	5.12	5.10 ± 0.06
V4C2H4	4.83	4.84	4.81	4.98	4.87	4.87 ± 0.05
Cr4C2H4	4.80	4.80	4.78	4.91	4.82	4.82 ± 0.03
Zr4C2H4	5.51	5.50	5.48	5.67	5.56	5.54 ± 0.04
Nb4C2H4	5.25	5.24	5.22	5.39	5.29	5.28 ± 0.03
Ta4C2H4	5.25	5.25	5.22	5.37	5.29	5.28 ± 0.03
W4C2H4	5.54	5.74	5.49	5.51	5.56	5.58 ± 0.08

---

Table 4 Interlayer spacing calculated using different dispersion correction methods for hydroxyl-terminated MXenes

System	D3 (Å)	D2 (Å)	XDM (Å)	vdW-DF2 (Å)	Without vdW (Å)	Mean standard deviation
Sc4C2(OH)4	7.75	7.69	7.63	8.01	8.01	7.82 ± 0.13
Ti4C2(OH)4	7.49	7.47	7.39	7.83	7.81	7.60 ± 0.15
V4C2(OH)4	7.27	7.25	7.18	7.66	9.5	7.77 ± 0.7
Cr4C2(OH)4	7.00	6.92	6.86	8.02	8.07	7.21 ± 0.16
Zr4C2(OH)4	7.87	7.79	7.78	8.21	8.15	7.96 ± 0.12
Nb4C2(OH)4	7.63	7.52	7.55	7.95	9.72	8.07 ± 0.64
Mo4C2(OH)4	7.11	7.01	7.02	7.4	9.2	7.55 ± 0.66
Hf4C2(OH)4	7.88	7.55	7.68	8.84	9.86	8.24 ± 0.80
Ta4C2(OH)4	8.94	7.44	8.95	9.1	9.47	8.70 ± 0.67

---

Overall, the dispersion-corrected approaches reduce separation compared to the non-corrected calculations; the magnitude of correction varies among the methods. In general, for MXenes based on Sc, Ti, V, Zr, Nb, Hf, and Ta, the order of interlayer distance concerning termination is H < O < F < OH, where H-terminated MXenes have the lowest interlayer distance, though in few cases F and O have very close values. This ordering reflects the increasing steric bulk from H to OH. However, this trend is slightly different for half-filled metals, like Cr, Mo, and W; the M–F bond benefits from strong ionic interactions and less steric hindrance than the M–O bond, causing the ordering to shift to H < F < O < OH. For fluorine-terminated MXenes, a similar reduction is observed when dispersion corrections are applied. The mean standard deviations across different methods (approximately 0.21–0.28 Å) suggest a moderate sensitivity of the interlayer distance to the specific dispersion correction employed; this indicates the robustness of different dispersion-correction methods for predicting interlayer spacing with few exceptions. MXenes, consistent with hydrogen terminations, yield the shortest interlayer spacing. This trend is particularly pronounced for 3d metals (e.g., Ti and V), implying that hydrogen termination facilitates stronger interlayer coupling. MXenes with hydroxyl terminations display a large interlayer spacing compared to other terminations. This could be attributed to its substantial steric bulk.

As we go along the 3d series (Sc, Ti, V, and Cr), both the M–C and M–X bonds (where X is the termination O, F, H, and OH) contract. For example, in Sc₂CO₂, M–C and M–O bond lengths are 2.51 Å and 2.01 Å, respectively, whereas in Cr₂CO₂, they decrease to 2.02 Å and 1.91 Å, respectively. This contraction arises due to the increase in effective nuclear charge,^53,54 which strengthens the intralayer interactions and is consistent with the decrease in atomic and ionic radii as we go along the period. The bond distances for M–C and M–X are presented in Table S1 (SI). For O-terminated MXenes, this contraction in size results in a decrease in the effective interlayer interaction, resulting in increasing interlayer spacing as we go along the period. Similarly, for the half-filled Cr-based MXene, there is an increase in the interlayer spacing of 7.14 Å compared to 6.98 Å in V₄C₂O₄, even though the Cr–O bond is shorter. The d⁵ half-filled configuration introduces the occupation of antibonding orbitals (as predicted by the ligand field theory),⁵⁵ which weakens the overall bonding energy. This reduces attractive interactions across layers, leading to increased interlayer spacing despite a shorter bond. For the 4d series (Zr, Nb, and Mo), the interlayer spacing increases from Zr₄C₂O₄ to Nb₄C₂O₄ and then slightly decreases for Mo₄C₂O₄. The 4d orbitals are more diffuse compared to 3d orbitals, leading to effective orbital overlapping in the intralayer, which weakens the interlayer coupling and thus increases the interlayer spacing. They are also less prone to antibonding destabilization compared to 3d and have more effective screening. In particular, in Mo, stabilization from screening leads to a slight decrease in the spacing after Nb. In the 5d series (Hf, Ta, and W), despite the significant increase in atomic size from 4d to 5d, relativistic contraction⁵⁴ significantly shortens the 5d orbitals; combined with the additional 4f screening, the interlayer spacing decreases as we go along the 5d series due to strong interlayer coupling. For F-terminated MXenes, the interlayer spacing decreases as we go along the period for the 3d, 4d, and 5d series. The high electronegativity of fluorine leads to a strong ionic nature in the M–F bonds, and combined with the increase in bond polarity, which enhances the attraction between adjacent layers, the interlayer spacing decreases. For the H-terminated MXenes, hydrogen being the smallest termination, the layers are more compact, and due to minimal steric hindrance, they allow the layers to settle closer together. This results in a decrease in interlayer spacing for the 3d and 4d series as we go along the period. For the 5d series, the results showed an increase in the spacing. This could be attributed to the repulsion of electron clouds in adjacent layers due to the more diffuse electron density in the heavier 5d metals. Similarly, for OH-terminated MXenes, due to the bulky nature of OH and larger M–(OH) bond distances compared to other terminations, they yield large interlayer spacing. Despite the bulkiness of OH, the M–OH bond contracts and becomes stronger over time due to the increased effective nuclear charge. For the 3d and 4d series, a decrease in the interlayer spacing is observed.

In continuation with the above discussion, the trends down the group for the d² series (Ti, Zr, and Hf, where d² refers to the valence configuration) for O-terminated MXenes indicate that the interlayer spacing decreases. The effective screening in the Zr despite its increase in atomic radii compared to Ti causes a net decrease in interlayer distance. For the d³ series (V, Nb, and Ta), as we go down the group, there is an increase in M–O bond distances, leading to weaker intralayer interaction. Stronger screening in Nb and Ta stabilizes the interlayer distance despite the increasing size, leading to an initial increase and then saturation. For the d⁵ series (Cr, Mo, and W), the spacing decreases. The better screening and relativistic effects of Mo and W lead to a decrease in spacing. For the F-terminated MXenes, in the d² series (Ti, Zr, Hf), the spacing increases from Ti–Zr due to longer M–F bonds. However, from Zr to Hf there can be a slight decrease or increase in the spacing, depending on the nature of dispersion correction method. The d³ series (V, Nb, and Ta) follows a similar trend. In the case of the d⁵ series (Cr, Mo, and W), the spacing decreases due to stronger bonding and relativistic effects in W, overcoming the size increase from Cr. For d², d³, and d⁵ series with H and OH-terminations, the spacing generally increases as we go down the group, with few exceptions. The increase in size is mainly due to the elongation of the M–H and M–OH bonds.

Since standard PBE-DFT does not account for long-range van der Waals forces, it overestimates interlayer distances in MXenes, leading to too widely spaced layers. With the addition of fixed C₆ coefficients and a simple damping function, Grimme's D2 correction addresses this. It shows a general downward trend in spacing down each group (e.g., Ti → Zr → Hf) well but sometimes overdoses on bond contraction and misses on small details, such as partial antibonding in half-filled systems. Grimme's D3 also optimizes dispersion with three-body terms and environment-dependent C₆ values; it can overcorrect in very polar systems, giving unrealistically large spacings, but it does well to simulate the way screening and relativistic effects bring an increase → slight decrease trend (e.g., Ti → Zr → Hf). The XDM method gives the lowest and most chemically detailed spacings by correctly simulating dispersion, which is directly derived from electron density. By precisely accounting for local polarizability and coordination variations, the XDM approach, which derives dispersion directly from the electron density, yields the lowest and most chemically nuanced spacings despite being more computationally costly and sensitive to the base functional. Notably, for many MXenes, the vdW-DF2 functional yields interlayer spacings higher than even uncorrected PBE. This anomalous behavior has also been reported in graphite and MoS₂.²⁶ This is most likely a result of an unphysical weakening of interlayer binding caused by an overestimation of the exchange–correlation contributions in the nonlocal functional. Therefore, vdW-DF2 results should be interpreted cautiously as upper bounds rather than as accurate predictions. Our data exhibit consistent method-dependent trends: Grimme's D2 tends to produce the smallest interlayer spacings at equilibrium (tighter binding), while vdW-DF2 tends to produce the largest spacings and, in a few instances, even larger than the without dispersion-corrected PBE result. D3 and XDM tend to place themselves in the middle positions, but their ranking between themselves and D2 changes with metal type and termination (see Tables 1–4). The average standard deviations between methods (∼0.21–0.28 Å) reveal a moderate sensitivity of spacing to the applied dispersion correction.

From a qualitative perspective, typical XRD investigations of Ti₃C₂T_x MXenes provide interlayer spacings in the ∼9–12 Å range with variation depending on termination and intercalation conditions.^12,47 Furthermore, recent experiments on MXenes have obtained interlayer spacings substantially larger than most terminations: for example, Zhao et al. (Adv. Funct. Mater. 2020)⁵⁶ report ∼1.77 nm in freestanding Nb₄C₃T_x films. Pre-intercalation approaches can swell Ti₃C₂T_x up to ∼2.2 nm interlayer distance under certain alkylammonium cation treatments.⁵⁷ These reports demonstrate realistic limits to which interlayer spacing may be adjusted through termination and intercalation, placing the lower values of our set into context. In the absence of direct experimental benchmarks for the MXenes investigated here, these reports serve as plausibility checks. Taken collectively, the tabulated data indicate that D3 and XDM tend to give physically sensible spacings, with D2 and vdW-DF2 forming the lower and higher extremes, respectively. The plots of the above results are presented in the SI (Fig. S1–S4).

The data in Table 1 demonstrate systematic differences among the dispersion schemes in terms of how they describe interlayer binding: D2 consistently predicts the most compressed structures, vdW-DF2 yields the most expanded ones, and D3 and XDM fall between these two limits. These trends have their origin in the underlying parametrization of each functional and reflect their expected qualitative behaviour. We note that these results are not intended to imply the strict numerical superiority of any single dispersion correction across all MXene chemistries. Rather, Table 1 is best understood as a demonstration of the relative trends of each scheme, to which the structural and mechanical variations discussed throughout this work are largely attributable.

3.2 Mechanical property: C₁₁ elastic modulus

In the case of layer systems, the accurate modeling of interatomic interactions plays a prominent role in controlling the mechanical properties of MXenes, such as Young's modulus, tensile strength, and elastic constants. The manner in which the layers respond to mechanical deformation is mainly determined by dispersion corrections. The estimated strength and versatility of MXenes in practical use can be influenced by mechanically unrealistic predictions due to incorrect handling of van der Waals forces. This work compares the influence of various dispersion correction methods on mechanical predictions. This discussion stresses the importance of accurate dispersion modeling in predicting MXene mechanical stability.

In this section, we calculate the C₁₁ elastic constant, which corresponds to the in-plane elastic moduli. The crystal's response to an externally imposed stress or strain is determined by its elastic stiffness, thereby revealing details regarding the structural stability, bonding properties, and mechanical stability. The results are shown in Tables 5–8.

Table 5 Elastic constant C₁₁ calculated using different dispersion correction methods for oxygen-terminated monolayer MXenes

System	Calculated C11 (N m^−1)
	D3	D2	XDM	vdW-DF2	Without vdW
Sc2CO2	205.79	221.30	204.71	232.68	209.35
Ti2CO2	290.76	298.77	285.73	299.11	292.79
V2CO2	249.73	249.90	246.70	257.41	251.36
Cr2CO2	218.27	213.61	216.94	224.85	222.67
Zr2CO2	304.35	327.38	295.36	327.16	307.07
Nb2CO2	339.82	338.72	338.72	365.16	344.06
Mo2CO2	259.46	252.50	258.61	279.87	263.33
Hf2CO2	353.54	319.90	358.46	339.33	364.09
Ta2CO2	321.29	321.29	359.97	374.34	365.96
W2CO2	302.81	261.04	300.65	324.02	306.17

---

Table 6 Elastic constant C₁₁ calculated using different dispersion correction methods for fluorine-terminated monolayer MXenes

System	Calculated C11 (N m^−1)
	D3	D2	XDM	vdW-DF2	Without vdW
Sc2CF2	95.04	108.75	114.52	94.70	98.30
Ti2CF2	178.21	141.76	313.19	178.17	183.26
V2CF2	231.42	246.91	255.86	230.24	233.82
Cr2CF2	148.50	168.40	163.48	147.46	151.32
Zr2CF2	178.24	190.90	194.03	176.39	180.85
Nb2CF2	153.90	136.96	162.53	153.45	157.52
Mo2CF2	143.10	155.84	154.29	141.63	145.88
Hf2CF2	119.49	131.56	129.76	117.63	121.43
Ta2CF2	199.55	183.78	230.49	199.78	203.16
W2CF2	171.30	192.80	187.07	170.12	174.30

---

Table 7 Elastic constant C₁₁ calculated using different dispersion correction methods for hydrogen-terminated monolayer MXenes

System	Calculated C11 (N m^−1)
	D3	D2	XDM	vdW-DF2	Without vdW
Sc2CH2	161.61	174.18	159.13	170.13	163.77
Ti2CH2	212.83	218.37	208.93	221.40	221.96
V2CH2	187.40	191.75	187.87	199.20	190.57
Cr2CH2	138.16	140.41	136.85	148.01	137.80
Zr2CH2	208.48	230.79	207.42	220.10	215.56
Nb2CH2	245.73	242.99	243.11	256.91	247.86
Mo2CH2	261.77	244.20	260.79	275.94	267.61
Ta2CH2	257.96	216.19	253.24	275.02	263.19
W2CH2	286.12	253.52	286.06	304.65	303.43

---

Table 8 Elastic constant C₁₁ calculated using different dispersion correction methods for hydroxyl-terminated monolayer MXenes

System	Calculated C11 (N m^−1)
	D3	D2	XDM	vdW-DF2	Without vdW
Sc2C(OH)2	185.29	196.63	189.30	172.55	179.50
Ti2C(OH)2	232.63	228.84	244.09	208.07	221.67
V2C(OH)2	233.49	224.78	240.48	245.81	188.80
Zr2C(OH)2	224.49	240.79	232.65	236.77	216.87
Mo2C(OH)2	219.65	232.68	219.14	231.77	156.22
Ta2C(OH)2	245.81	190.98	275.21	263.08	224.87

---

Thermo_pw returns elastic constants scaled to the total supercell volume (slab + vacuum) in units of kbar. Subsequently, we scale raw Thermo_pw outputs to slab-normalized three-dimensional moduli and to two-dimensional in-plane stiffness in a reproducible manner.⁵¹ Here, C^kbar₁₁ denotes the Thermo_pw output (kbar), c denotes the total supercell height (Å; slab + vacuum; typical value used here c ∼35 Å), and s denotes the interlayer distance (Å) used as the effective thickness (in this work, we use the carbon–carbon interlayer distance for carbon-containing MXenes, which is consistent with ref. 51). The corrected 3D modulus (GPa) is then

To obtain physically meaningful two-dimensional elastic stiffness for MXenes, we convert the corrected modulus into N m⁻¹ using the following equation:

This thickness-normalized 2D unit (N m⁻¹) is independent of the vacuum size and is the preferred descriptor for 2D mechanical properties. All values reported in the main text use N m⁻¹.

It is inherently ambiguous to define the effective thickness of a 2D material, as has been widely debated for graphene and other layered crystals. Reported values range from ∼0.66 Å (estimated covalent radius) to 3.4 Å (interlayer spacing of graphite).^58,59 For MXenes, Kurtoglu et al.⁶⁰ used the metal–metal distance between the outermost M layers, while Zha et al.⁵¹ used bilayer optimization and extracted the interlayer spacing. In this study, we consistently define the effective thickness as the carbon–carbon interlayer distance. This definition preserves internal consistency and is compatible with previous computational investigations of terminated MXenes.

The maximum C₁₁ values are exhibited by oxygen-terminated MXenes in all metal groups, signifying enhanced rigidity and in-plane stability, except for a few exceptions. The overall order of increasing C₁₁ concerning termination is F ∼ OH < H < O. Along the period for the 3d (Sc, Ti, V, and Cr), 4d (Zr, Nb, and Mo) and 5d (Hf, Ta, and W) series with O terminations, the C₁₁ rises across early metals as M–O bond length contracts. Then, at half-filled metals (Cr, Mo, and W), it decreases due to antibonding d-orbital occupation and steric repulsion. A similar trend is also observed for F and OH terminations. For hydrogen-terminated MXenes, this trend is observed for the 3d series. In the 4d and 5d series, the results show a gradual increase in the C₁₁ without any dip. The rise in C₁₁ for Mo₂CH₂ and W₂CH₂ results from a weaker antibonding effect in the 4d and 5d metals, greater d-orbital overlap, and relativistic bond contraction. The low steric nature of hydrogen facilitates the preservation of strong, short M–H bonds even in heavier transition metals. For the O–, F=, and H-terminated MXenes down the group, there is a gradual rise in C₁₁ values for d², d³, and d⁵ series mainly due to the more diffused d-orbitals and enhanced screening. The results also provide a valuable design guide for applications requiring high in-plane rigidity (e.g., structural reinforcement and coatings); O-terminated MXenes are optimal. For systems where a softer in-plane response or ductility is preferred, F- or OH-terminated MXenes may be more suitable.

The selection of the van der Waals correction approach has a strong bearing on the predicted MXene elastic constants (C₁₁) by tuning the intralayer bonding strength. Grimme D3 often provides a balanced description and manages to reproduce the major stiffness trends across the metal series. For instance, it accurately predicts the rise in stiffness from Ti₂CO₂ (290.76 N m⁻¹) to Zr₂CO₂ (304.35 N m⁻¹) and Hf₂CO₂ (353.54 N m⁻¹), and the modest decline at W₂CO₂ (302.81 N m⁻¹). XDM typically improves upon this trend in highly polarizable or heavy-metal complexes (Mo- and W-based) due to its density-dependent dispersion but at increased computational expense. Grimme D2 tends to give systematically larger stiffnesses for rigid M–O bonds due to fixed dispersion coefficients, while vdW-DF2 can underpredict stiffness in several cases, reflecting its known tendency toward more diffuse short-range repulsion. The full plots are provided in SI, Fig. S5–S8. We expand on these observations below with a quantitative discussion of elastic-unit conversions, termination-dependent C₁₁ trends, and methodological sensitivity across vdW schemes.

Throughout the diverse set of structural and mechanical quantities presented in this work, D3 and XDM yield internally consistent predictions. D3 reproduces chemically reasonable spacing contractions, relativistic effects down the transition-metal series, and stable stiffness trends across terminations. Due to its density-dependent construction, XDM performs particularly well for highly polarizable or heavy-metal MXenes, such as Mo- and W-based compounds. In turn, D2 tends to overbind due to fixed dispersion coefficients, while vdW-DF2 systematically underbinds, resulting in diffuse interlayer separations. Our assessment of D3 and XDM therefore reflects their overall balanced performance across properties rather than a conclusion drawn from any single metric.

All C₁₁ values reported here are expressed directly in two-dimensional units (N m⁻¹), obtained from the raw Thermo_pw elastic outputs following the conversion procedure described in the Methods section. Throughout our data set, the in-plane stiffness values span roughly from ∼100 to 360 N m⁻¹. Examination of Tables 5–8 indicates that oxygen-terminated MXenes most commonly exhibit the largest in-plane stiffness (e.g., range for O-terminated systems in our set: ∼200–360 N m⁻¹), while hydrogen-terminated MXenes most commonly have the densest lattice parameters and frequently have the second-highest stiffness (H-terminated range in our set: ∼150–300 N m⁻¹). Thus, oxygen termination mostly correlates with the highest absolute C₁₁ in this dataset, while hydrogen termination features small lattice constants together with relatively high in-plane stiffness. The larger C₁₁ on average for O-termination can be explained by more stable, directional M–O bonding that enhances in-plane stiffness; hydrogen, as a sterically small atom, creates tight lattices and close M–H contacts that retain high stiffness but generally do not exceed M–O rigidity. These trends are robust over the vdW corrections we considered (see Tables 5–8).

Earlier computational and experimental works routinely reported the elastic constants of MXenes. For instance, Zha et al.⁵¹ provided C₁₁ values of 366 GPa for Ti₂CO₂ and 432 GPa for Nb₂CO₂ with PBE-D2. Experimental nanoindentation of multilayer Ti₃C₂T_x films produces moduli in the 330–480 GPa range.⁶¹ Corresponding measurements and computational reports for other MXene compositions also lie in this range.^62,63 Furthermore, a wider computational study by Sutherland et al.⁶⁴ contains C₁₁ values in Table S15 for different MXene compositions (M₃C₂T₂), and these reveal relative stiffness trends consistent with ours (increasing stiffness for O-terminated and heavier metals) although the absolute values depend on functional and thickness conventions. All the computed elastic constants in the present work are reported in two-dimensional units (N m⁻¹) alone, following the usual convention for 2D materials. Literature values in GPa are quoted only for contextual comparison and are not used for direct numerical benchmarking. Across the MXenes analyzed in this study, two-dimensional C₁₁ values range from ∼100 to 360 N m⁻¹ (depending on termination and chemistry), positioning and placing them within the mechanically robust range of known 2D materials and comparable to systems such as h-BN and graphene. Oxygen-terminated MXenes have the overall stiffest, and hydrogen-terminated MXenes are generally competitive and have dense lattice constants, making them a good choice for flexible and mechanically robust device applications. Our calculated C₁₁ values (∼100–360 N m⁻¹) fall within the mechanically relevant range reported for MXenes. For comparison, representative 2D benchmarks include graphene (∼342 N m⁻¹),⁵⁸ h-BN (∼240–300 N m⁻¹),⁶⁵ and MoS₂ (∼120–180 N m⁻¹).⁶⁶ More recent studies report stiffness values of ∼120–180 N m⁻¹ for Janus TMDs^67,68 and systematic ranges for group-VI TMD monolayers.⁶⁹ These comparisons highlight that O- and H-terminated systems are among the stiffest 2D materials. We note, however, that these findings are relative computational tendencies since experimental benchmarks for hydrogen-terminated MXenes are currently unavailable.

The calculated stiffnesses are sensitive to equilibrium geometry and to the selected dispersion correction: those that decrease interlayer spacing tend to increase the computed in-plane stiffness, and vice versa for schemes that increase spacing. Although the numerical values of C₁₁ depend on the chosen conversion convention (interlayer spacing in this work), the relative trends between terminations and between vdW models are consistent. Higher-level theory approaches are valuable for establishing precise absolute benchmarks on specific systems but are expensive computationally; our systematic dataset thus provides a reproducible baseline to guide targeted high-accuracy calculations and experiments.

3.3 Lattice constant

MXene lattice constants give a basic measure of structural parameters of MXenes, with an impact on properties like electronic band structure and density of states. In contrast to binding energies, lattice constants can be well-predicted using ordinary DFT functionals and dispersion-corrected methods. The accuracy in predicting lattice constants can, however, create a false sense of the overall validity of a method, especially if it does not predict energies and interlayer forces correctly. Although in this work, the focus is on the investigation of different dispersion approaches, it is important to note that the choice of E_XC also affects the lattice constant predictions. In the current study, as mentioned in Section 2, we employ the PBE functional, and several pioneering studies have been conducted to showcase the ability of PBE to accurately predict structural properties.^{3,37,44,70,71} Hence, the comparison of lattice constants based on various dispersion corrections provides a baseline study, providing indications of each approach's overall accuracy in representing MXene's structural architecture and locating possible deficiencies in energy calculations. For representative MXenes where experimental data are available, our values agree with those reported: Ti₂CO₂ = 3.01 Å, V₂CF₂ = 2.99 Å, Nb₂CO₂ = 3.14 Å, and Mo₂CF₂ = 3.28 Å.^5,72–75 Although experimental benchmarks are not directly found for most compositions examined here, these representative systems give valuable reference points, indicating that our values are consistent with those reported: Ti₂CO₂ = 3.01 Å, V₂CF₂ = 2.99 Å, Nb₂CO₂ = 3.14 Å, and Mo₂CF₂ = 3.28 Å.^5,72–75 Although experimental benchmarks are not directly found for most compositions examined here, these representative systems provide valuable reference points and indicate that our values are within the expected experimental window. The results for all metal series and terminations are displayed in Fig. 2–5.

Fig. 2 Lattice constant values of 9 different oxygen-terminated MXenes calculated with 4 different dispersion corrections and without correction.

Fig. 3 Lattice constant values of 10 different fluorine-terminated MXenes calculated with 4 different dispersion corrections and without correction.

Fig. 4 Lattice constant values of 10 different hydrogen-terminated MXenes calculated with 4 different dispersion corrections and without correction.

Fig. 5 Lattice constant values of 10 different hydroxyl-terminated MXenes calculated with 4 different dispersion corrections and without correction.

Across all the chosen MXene compositions and dispersion approaches, the termination group (O, F, H, and OH) shows a trend concerning the lattice constant value. The hydrogen-terminated MXenes have the smallest lattice constant among the chosen terminations. This is directly a result of the atomic size and M–H bond distance. The oxygen-terminated MXenes are marginally longer in comparison to H-terminated MXenes. Although the M–O bonds are larger, their strong covalent character compensates for this, keeping the cell compact. Fluorine-terminated MXenes have a larger lattice constant than H or O. This is due to the longer M–F bond lengths and moderate steric hindrance of the F atom. The largest lattice constant is observed in hydroxyl-terminated MXenes. Although in few cases, it is comparable to F-terminated MXenes, the large lattice constant is attributed to its bulky OH group and its high dipole moment pushing the neighboring layers apart. In summary, the observed order by increasing the lattice constant is H < O < F ≤ OH.

Across the period, for the 3d series (Sc, Ti, V, and Cr), with all four terminations, the lattice constant decreases regardless of the dispersion correction used. The gradual decrease is due to the contraction in atomic size, bond lengths and an increase in effective nuclear charge as we go along the period. This trend is also observed in the 4d (Zr, Nb, and Mo) and 5d series (Hf, Ta, and W) with H-termination. This monotonous decrease for H-terminated MXenes is due to their minimal steric bulk and low polarizability. For the 4d and 5d series with O, F, and OH termination, the lattice constant initially decreases; then, it increases when it comes to half-filled metals, like Mo and W. Although the size decreases as we go along the period due to the half-filled configuration in Mo and W, the antibonding orbital population begins to weaken intralayer bonding. Along with steric repulsion, the polarizability of O, F, and OH offsets orbital contraction and results in an increase in the lattice constant.

For the d² (Ti, Zr, and Hf), d³ (V, Nb, and Ta), and d⁵ (Cr, Mo, and W) series with O, F, H, and OH terminations, down the group, there is an initial increase and then a marginal decrease in the lattice constant. The initial increase is driven by an increase in M–X and M–C, i.e., larger atomic and ionic radii, so the lattice expands. However, upon coming to the 5d metals (Hf, Ta, and W), there is a strong relativistic contraction of the 5d orbitals, and combined with effective electronic screening, this leads to the shortening of these bonds, resulting in a reduction in the lattice constant.

Overall, D3 and XDM corrections provided lattice constants in the closest agreement with the reported values, verifying their balanced treatment of long- and short-range interactions. D2 tends to overestimate lattice constants for heavier metals owing to fixed dispersion coefficients, while vdW-DF2 sometimes overestimates²⁶ owing to its diffuse treatment of short-range repulsion. We note, nonetheless, that absolute values are modestly sensitive to functional choice. Relative trends between metals and terminations, however, are robust and reproducible. In situations where structural dimensions must be controlled accurately (i.e., in nanoelectronics), D3 and XDM yield the most reliable predictions.

4. Conclusions

In this study, we have shown that a proper treatment of van der Waals interactions is crucial for predicting the structural and mechanical properties of MXenes. Using four widely adopted dispersion-corrected methods (Grimme D2, D3, XDM, and vdW-DF2) within the PBE-DFT, we systematically examined ten M₂CX_T compositions (M = Sc, Ti, V, Cr, Zr, Nb, Mo, Hf, Ta, W; X = H, O, F, OH). Our findings reveal several termination- and metal-dependent design criteria. Interlayer spacings range from 4.8 to 9.7 Å and always adhere to the sequence H < O < F ≤ OH, where hydrogen-terminated heavy metals (e.g., W₂CH₂) result in ultracompact geometries and oxygen-terminated early-transition metals (e.g., Sc₂CO₂) result in expanded interlayers that are highly favorable for ion intercalation and storage. Lattice constants are strongly influenced by termination chemistry, varying from dense H-terminated structures to loose OH-terminated systems, with periodic and group-wise systematic trends resulting from bond contraction, occupation of antibonding orbitals, and relativistic effects.

The in-plane stiffness values calculated here (∼100–360 N m⁻¹) place many MXenes at the top end of previously known two-dimensional materials. Oxygen termination always shows the stiffest response, while hydrogen terminations give dense lattices with similarly strong mechanical responses, both desirable for mechanically robust and flexible-device use. Fluorine and hydroxyl terminations result in softer responses and greater spacings, expanding the tunability of the MXene design.

Among dispersion corrections, D3 and XDM most reliably reproduce experimental trends, D3 provides a stable compromise of accuracy and efficiency, and XDM captures polarizable, heavy-metal systems with high fidelity at higher computational expense. In contrast, D2 tends to overbind, and vdW-DF2 often underbinds. Although absolute moduli vary with the thickness convention adopted and the computational parameters used, the relative trends between terminations and between vdW schemes are robust. Our dataset thus furnishes an internally coherent and reproducible baseline against which future high-level theoretical benchmarks and experimental studies can be compared.

In summary, this work establishes systematic structure–property relationships in MXenes by linking dispersion corrections, termination chemistry, and transition-metal identity. The findings reported here provide actionable design principles: O-terminations optimize stiffness, H-terminations allow dense and competitive mechanical performance, and F/OH-terminations add flexibility and tunability. By clarifying how interlayer and in-plane properties are reliant on computational treatment and chemical environment, this study offers both methodological insight and practical routes to customizing MXenes in energy storage, catalysis, and flexible electronics.

Conflicts of interest

The authors declare no competing financial interests.

Data availability

With reference to the manuscript titled “Investigation of van der Waals Interactions in two-dimensional MXenes via first-principles simulations,” data for this article containing CIF files are available at GitHub https://github.com/PPISR/PPISR_Supplementary/blob/main/Theoretical_Sciences/All_cif.zip.

Supplementary information (SI) is available. See DOI: https://doi.org/10.1039/d5cp03872c.

Acknowledgements

The authors gratefully acknowledge the Admar Mutt Education Foundation (AMEF), Bengaluru, for funding and providing research facilities and Gitam University, Bengaluru, for access to high-performance computing clusters.

References

1. K. S. Novoselov, A. K. Geim, S. V. Morozov, D.-E. Jiang, Y. Zhang, S. V. Dubonos, I. V. Grigorieva and A. A. Firsov, Science, 2004, 306, 666–669.
2. A. K. Geim and I. V. Grigorieva, Nature, 2013, 499, 419–425.
3. M. Naguib, M. Kurtoglu, V. Presser, J. Lu, J. Niu, M. Heon, L. Hultman, Y. Gogotsi and M. W. Barsoum, Adv. Mater., 2011, 23, 4248–4253.
4. M. Naguib, V. N. Mochalin, M. W. Barsoum and Y. Gogotsi, Adv. Mater., 2014, 26, 992–1005.
5. M. Naguib, O. Mashtalir, J. Carle, V. Presser, J. Lu, L. Hultman, Y. Gogotsi and M. W. Barsoum, ACS Nano, 2012, 6, 1322–1331.
6. P. Simon and Y. Gogotsi, Nat. Mater., 2020, 19, 1151–1163.
7. X. Li, Z. Huang, C. E. Shuck, G. Liang, Y. Gogotsi and C. Zhi, Nat. Rev. Chem., 2022, 6, 389–404.
8. J. Zou, J. Wu, Y. Wang, F. Deng, J. Jiang, Y. Zhang, S. Liu, N. Li, H. Zhang and J. Yu, Chem. Soc. Rev., 2022, 51, 2972–2990.
9. M. Ghidiu, M. Naguib, C. Shi, O. Mashtalir, L. Pan, B. Zhang, J. Yang, Y. Gogotsi, S. J. Billinge and M. W. Barsoum, Chem. Commun., 2014, 50, 9517–9520.
10. M. R. Lukatskaya, O. Mashtalir, C. E. Ren, Y. Dall’Agnese, P. Rozier, P. L. Taberna, M. Naguib, P. Simon, M. W. Barsoum and Y. Gogotsi, Science, 2013, 341, 1502–1505.
11. L. Ding, Y. Wei, L. Li, T. Zhang, H. Wang, J. Xue, L.-X. Ding, S. Wang, J. Caro and Y. Gogotsi, Nat. Commun., 2018, 9, 155.
12. A. Lipatov, M. Alhabeb, M. R. Lukatskaya, A. Boson, Y. Gogotsi and A. Sinitskii, Adv. Electron. Mater., 2016, 2, 1600255.
13. M. R. Lukatskaya, S. Kota, Z. Lin, M.-Q. Zhao, N. Shpigel, M. D. Levi, J. Halim, P.-L. Taberna, M. W. Barsoum and P. Simon, MXenes, Jenny Stanford Publishing, 2023, pp. 723–743.
14. S. Cahangirov, M. Topsakal, E. Aktürk, H. Şahin and S. Ciraci, Phys. Rev. Lett., 2009, 102, 236804.
15. H. Liu, A. T. Neal, Z. Zhu, Z. Luo, X. Xu, D. Tománek and P. D. Ye, ACS Nano, 2014, 8, 4033–4041.
16. S. Zhang, Z. Yan, Y. Li, Z. Chen and H. Zeng, Angew. Chem., Int. Ed., 2015, 54, 3112–3115.
17. I. Boustani, Surf. Sci., 1997, 370, 355–363.
18. A. J. Mannix, X.-F. Zhou, B. Kiraly, J. D. Wood, D. Alducin, B. D. Myers, X. Liu, B. L. Fisher, U. Santiago and J. R. Guest, Science, 2015, 350, 1513–1516.
19. M. Krinninger, F. Kraushofer, N. B. Refvik, M. Blum and B. A. Lechner, ACS Appl. Mater. Interfaces, 2024, 16, 27481–27489.
20. J. T. Küchle, A. Baklanov, A. P. Seitsonen, P. T. Ryan, P. Feulner, P. Pendem, T.-L. Lee, M. Muntwiler, M. Schwarz and F. Haag, 2D Mater., 2022, 9, 045021.
21. M. Frisch, G. Trucks, H. Schlegel, G. Scuseria, M. Robb, J. Cheeseman, G. Scalmani, V. Barone, G. Petersson and H. Nakatsuji, Gaussian Inc., Wallingford CT, 2016, 1, 572.
22. P. Giannozzi, S. Baroni, N. Bonini, M. Calandra, R. Car, C. Cavazzoni, D. Ceresoli, G. L. Chiarotti, M. Cococcioni and I. Dabo, J. Phys.: Condens. Matter, 2009, 21, 395502.
23. F. Neese, F. Wennmohs, U. Becker and C. Riplinger, J. Chem. Phys., 2020, 152, 224108.
24. J. Klimeš, D. R. Bowler and A. Michaelides, Phys. Rev. B: Condens. Matter Mater. Phys., 2011, 83, 195131.
25. L. Girifalco and R. Lad, J. Chem. Phys., 1956, 25, 693–697.
26. K. Berland, V. R. Cooper, K. Lee, E. Schröder, T. Thonhauser, P. Hyldgaard and B. I. Lundqvist, Rep. Prog. Phys., 2015, 78, 066501.
27. M. J. Allen and D. J. Tozer, J. Chem. Phys., 2002, 117, 11113–11120.
28. Y. Zhang, W. Pan and W. Yang, J. Chem. Phys., 1997, 107, 7921–7925.
29. O. A. Vydrov and T. Van Voorhis, J. Chem. Phys., 2010, 133, 24.
30. K. Lee, É. D. Murray, L. Kong, B. I. Lundqvist and D. C. Langreth, Phys. Rev. B: Condens. Matter Mater. Phys., 2010, 82, 081101.
31. S. Grimme, J. Comput. Chem., 2006, 27, 1787–1799.
32. A. D. Becke and E. R. Johnson, J. Chem. Phys., 2005, 123, 154101.
33. A. Marini, P. García-González and A. Rubio, Phys. Rev. Lett., 2006, 96, 136404.
34. F. Furche, J. Chem. Phys., 2008, 129, 114105.
35. R. Hatcher, J. A. Kittl and C. Bowen, arXiv, 2019, preprint arXiv:1903.05550 DOI:10.48550/arXiv.1903.05550.
36. J. Klimeš and A. Michaelides, J. Chem. Phys., 2012, 137, 120901.
37. S. Grimme, A. Hansen, J. G. Brandenburg and C. Bannwarth, Chem. Rev., 2016, 116, 5105–5154.
38. M. Dion, H. Rydberg, E. Schröder, D. C. Langreth and B. I. Lundqvist, Phys. Rev. Lett., 2004, 92, 246401.
39. A. Tkatchenko and M. Scheffler, Phys. Rev. Lett., 2009, 102, 073005.
40. P. Jurečka, J. Šponer, J. Černý and P. Hobza, Phys. Chem. Chem. Phys., 2006, 8, 1985–1993.
41. J. Rezác, K. E. Riley and P. Hobza, J. Chem. Theory Comput., 2011, 7, 2427–2438.
42. B. Brauer, M. K. Kesharwani, S. Kozuch and J. M. Martin, Phys. Chem. Chem. Phys., 2016, 18, 20905–20925.
43. S. A. Tawfik, T. Gould, C. Stampfl and M. J. Ford, Phys. Rev. Mater., 2018, 2, 034005.
44. J. P. Perdew, K. Burke and M. Ernzerhof, Phys. Rev. Lett., 1996, 77, 3865.
45. K. Momma and F. Izumi, Appl. Crystallogr., 2008, 41, 653–658.
46. Thermo_pw is a driver of the QUANTUM ESPRESSO (QE) routines which provides an alternative organization of the QE workflow for the most common tasks. For more information see https://dalcorso.github.io/thermo_pw.
47. O. Mashtalir, M. Naguib, V. N. Mochalin, Y. Dall’Agnese, M. Heon, M. W. Barsoum and Y. Gogotsi, Nat. Commun., 2013, 4, 1716.
48. M. Khazaei, M. Arai, T. Sasaki, C. Y. Chung, N. S. Venkataramanan, M. Estili, Y. Sakka and Y. Kawazoe, Adv. Funct. Mater., 2013, 23, 2185–2192.
49. S. Grimme, J. Antony, S. Ehrlich and H. Krieg, J. Chem. Phys., 2010, 132, 154104.
50. T. Schultz, N. C. Frey, K. Hantanasirisakul, S. Park, S. J. May, V. B. Shenoy, Y. Gogotsi and N. Koch, Chem. Mater., 2019, 31, 6590–6597.
51. X.-H. Zha, K. Luo, Q. Li, Q. Huang, J. He, X. Wen and S. Du, Europhys. Lett., 2015, 111, 26007.
52. C. R. Rêgo, L. N. Oliveira, P. Tereshchuk and J. L. Da Silva, J. Phys.: Condens. Matter, 2015, 27, 415502.
53. C. J. Cramer, Essentials of computational chemistry: theories and models, John Wiley & Sons, 2013.
54. P. Pyykkö, Annu. Rev. Phys. Chem., 2012, 63, 45–64.
55. C. J. Ballhausen and H. B. Gray, Molecular orbital theory: an introductory lecture note and reprint volume, WA Benjamin New York, 1965.
56. S. Zhao, C. Chen, X. Zhao, X. Chu, F. Du, G. Chen, Y. Gogotsi, Y. Gao and Y. Dall’Agnese, Adv. Funct. Mater., 2020, 30, 2000815.
57. K. Liang, R. A. Matsumoto, W. Zhao, N. C. Osti, I. Popov, B. P. Thapaliya, S. Fleischmann, S. Misra, K. Prenger, M. Tyagi, E. Mamontov, V. Augustyn, R. R. Unocic, A. P. Sokolov, S. Dai, P. T. Cummings and M. Naguib, Adv. Funct. Mater., 2021, 31, 2104007.
58. C. Lee, X. Wei, J. W. Kysar and J. Hone, Science, 2008, 321, 385–388.
59. B. I. Yakobson and P. Avouris, in Carbon Nanotubes: Synthesis, Structure, Properties, and Applications, ed. M. S. Dresselhaus, G. Dresselhaus and P. Avouris, Springer Berlin Heidelberg, Berlin, Heidelberg, 2001, pp. 287–327 DOI:10.1007/3-540-39947-x_12.
60. M. Kurtoglu, M. Naguib, Y. Gogotsi and M. W. Barsoum, MRS Commun., 2012, 2, 133–137.
61. A. Lipatov, H. Lu, M. Alhabeb, B. Anasori, A. Gruverman, Y. Gogotsi and A. Sinitskii, Sci. Adv., 2018, 4, eaat0491.
62. D. Tan, C. Jiang, X. Cao, N. Sun, Q. Li, S. Bi and J. Song, RSC Adv., 2021, 11, 19169–19184.
63. Y. Ibrahim, A. Mohamed, A. M. Abdelgawad, K. Eid, A. M. Abdullah and A. Elzatahry, Nanomaterials, 2020, 10, 1916.
64. E. Sutherland, B. Traverso and N. A. Deskins, Mater. Adv., 2025, 6, 6787–6802.
65. A. Falin, Q. Cai, E. J. G. Santos, D. Scullion, D. Qian, R. Zhang, Z. Yang, S. Huang, K. Watanabe, T. Taniguchi, M. R. Barnett, Y. Chen, R. S. Ruoff and L. H. Li, Nat. Commun., 2017, 8, 15815.
66. S. Bertolazzi, J. Brivio and A. Kis, ACS Nano, 2011, 5, 9703–9709.
67. Q. Peng and S. De, Phys. Chem. Chem. Phys., 2013, 15, 19427–19437.
68. S. Imani Yengejeh, J. Liu, S. A. Kazemi, W. Wen and Y. Wang, ACS Omega, 2020, 5, 5994–6002.
69. R. C. Cooper, C. Lee, C. A. Marianetti, X. Wei, J. Hone and J. W. Kysar, Phys. Rev. B: Condens. Matter Mater. Phys., 2013, 87, 035423.
70. A. Dal Corso, J. Phys.: Condens. Matter, 2016, 28, 075401.
71. B. Anasori, M. R. Lukatskaya and Y. Gogotsi, 2D Metal Carbides and Nitrides (MXenes) for Energy Storage, 2023, pp. 677–722.
72. T. Sakhraoui and F. e Karlický, ACS Omega, 2022, 7, 42221–42232.
73. R. Ibragimova, P. Erhart, P. Rinke and H.-P. Komsa, J. Phys. Chem. Lett., 2021, 12, 2377–2384.
74. M. Khazaei, M. Arai, T. Sasaki, C.-Y. Chung, N. S. Venkataramanan, M. Estili, Y. Sakka and Y. Kawazoe, Adv. Funct. Mater., 2013, 23, 2185–2192.
75. A. N. Enyashin and A. L. Ivanovskii, J. Phys. Chem. C, 2013, 117, 13637–13643.

---