Chemical bonding at the metal–organic framework/metal oxide interface: simulated epitaxial growth of MOF-5 on rutile TiO₂

Abstract

Thin-film deposition of metal–organic frameworks (MOFs) is now possible, but little is known regarding the microscopic nature of hybrid hetero-interfaces. We first assess optimal substrate combinations for coherent epitaxy of MOFs based on a lattice matching procedure. We then perform a detailed quantum mechanical/molecular mechanical investigation of the growth of (011) MOF-5 on (110) rutile TiO₂. The lowest energy interface configuration involves a bidentate connection between two TiO₆ polyhedra with deprotonation of terephthalic acid to a bridging oxide site. The epitaxy of MOF-5 on the surface of TiO₂ was modelled with a forcefield parameterised to quantum chemical binding energies and bond lengths. The microscopic interface structure and chemical bonding characteristics are expected to be relevant to other hybrid framework-oxide combinations.

---

1 Introduction

Metal–organic frameworks (MOFs) are formed of inorganic metal clusters and organic ligands, self assembled into extended porous networks. Progress in the understanding of crystallisation mechanisms¹ and methods for controlling pore geometry and morphology^2,3 has allowed the fabrication of thin film MOFs for applications including molecular sensors, smart membranes and catalytic coatings.^4,5

There are two distinct approaches for fabricating MOF heterointerfaces.⁶ The first is the deposition of a preformed MOF onto a surface, e.g. using a one-pot solvothermal reaction. The second is templated growth where the surface of a substrate can be functionalised, for example, with –thio, –COOH or –CF₃ monolayers.^7–9 The latter is considered advantageous for uniform thin-film growth since it offers selective and directional growth with a reduced number of defects at interfacing sites.

There are numerous examples of MOF interface formation in the recent literature. For example, Ameloot et al. demonstrated controlled thickness and crystal size of thin film growth of HKUST-1 (a Cu²⁺ paddle-wheel structure with attractive physical properties¹⁰) on Cu metal electrodes using electrochemical deposition.¹¹ In accordance with classical nucleation theory, they report decreasing MOF crystal size with an increased applied voltage due to an increased rate of crystal growth. In 2010, Yoo et al. reported the heterostructured layered growth of the isorecticular IRMOF-3 and MOF-5 structures on the surface of Al₂O₃.¹² Using a core/shell approach, where one MOF is seeded on top of the other, a layered structure formed of Al₂O₃ substrate/IRMOF-3/MOF-5 layers was grown. Growth was demonstrated regardless of which MOF formed the core and shell. This approach allows for the construction of complex heterostructured devices with differing porosity and functional properties throughout the material. More recently, Fischer et al. reported a heterostructure between HKUST-1 and [Cu₂ndc₂dabco] (ndc: 1,4-naphthalene dicarboxylate, dabco: 1,4-diazabicyclo(2.2.2) octane), which were grown upon a pyridyl-functionalized Au substrate.¹³ The layered systems displayed enhanced adsorption affinity for small organic molecules, such as methanol, at standard conditions. A significant development has been the application of vapour phase deposition techniques: Ameloot et al. and Ritala et al. independently reported the gas phase deposition of ZIF-8 and MOF-5 onto surfaces using the chemical vapour deposition and atomic layer deposition, respectively. The optimisation of vapour phase approaches to the growth of MOF thin films would allow greater control over uniformity and thickness for device manufacturing.^14,15

Resolving the interface structure following epitaxial growth is challenging, especially for hybrid solids that are inherently ‘soft’ materials. Different experimental approaches have been attempted for determining the surface termination of MOFs, including an atomic-force microscopy examination of HKUST-1.¹⁶ Computational approaches to predicting the termination of MOF crystals have also provided useful guidance. Schmid et al. used classical and first principles methods to predict the surface structure and growth mechanism of HKUST-1.¹⁷ The calculated surface formation energies suggested that a ligand termination of the (111) surface is the lowest energy cleavage of the crystal. It should also be noted that complex defects are possible in MOFs, such as missing ligands and metal clusters during growth, further complicating surface structure characterisation.¹⁸

Following an initial screening procedure to identify lattice-matched MOF/substrate combinations, we report a combined classical and first principles investigation of the mechanism of epitaxy of MOF-5 on the (110) surface of rutile TiO₂. Favourable binding positions of the BDC linker forming MOF-5 are calculated with density functional theory (DFT) and used to re-parameterise an existing MOF forcefield to describe the interface. The (011) surface of MOF-5 is then interfaced with the (110) surface of TiO₂ and the resulting chemical interactions and surface reconstruction are reported.

2 Methodology

2.1 Lattice matching procedure

Screening for optimal materials combinations based on minimisation of lattice mismatch was conducted to find viable combinations of MOFs and binary materials including mainly metal oxides. We employed the electronic-lattice-site (ELS) procedure.¹⁹ The process involves cleaving all low-index surfaces of both the MOFs and substrates using the atoms object in the atomic simulation environment (ASE).²⁰ The mismatch of surface cell vectors (u × v) can then be calculated as a percentage difference, taking into account possible orientations and supercell expansions. The definition of surface cell parameters is given by the Zur and McGill scheme of defining a set of primitive vectors of the cell such that they are independent of any rotations or reflections of the lattice, thus allowing a complete identification of compatible surfaces for epitaxial interfacing.^21,22 We define an arbitrary cut-off of 8% mismatch, which takes into account the mechanical softness and flexibility of many hybrid frameworks.

2.2 First principles calculations

Total energy calculations were performed within the Kohn–Sham density functional theory framework using the PBEsol functional for the exchange–correlation potential with a D3 van der Waals correction in the QUICKSTEP module of the CP2K program.^23–25 PBE pseudopotentials of the analytical form of Goedecker, Teter and Hutter (GTH) were used to model the interaction between the valence electrons [Ti(3s²3p⁶4s²3d²), C(2s²2p²), O(2s²2p⁴), F(2s²2p⁵), and H(1s)] and the atomic cores.²⁶ Kohn–Sham orbitals are expanded with local Gaussian functions from the cc-TZ library for C, H, F and O with those from the DZV-GTH-PADE library for Ti.^27,28 The planewave cut-off for the auxiliary density basis set was set to 300 Ry with convergence of the self-consistent field set to 1 × 10⁻⁷ Ha.

2.3 Forcefield calculations

Analytical forcefield calculations used GULP^29,30 and the VMOF forcefield³¹ as parameterised for common MOFs. VMOF represents the interaction between metal and ligands by modified MM3 Buckingham potentials, as given in eqn (1), plus the Coulomb terms, and is an extension to a previous parameterised forcefield, BTW-FF.^32–34

Combination rules, where and , were used for the ε_ij and d⁰_ij values of the MM3 potential, that were parameterised for each metal to reproduce the structural and mechanical properties of a subset of MOFs including: MOF-5, IRMOF-10, MOF-650, UiO-66, UiO-67, MIL-125, NOTT-300, and MOF-74.

Intramolecular bonding parameters of the ligands are taken directly from the CHARMM library and charges derived using the charge equilibrium scheme of Gasteiger.^35–37 Formal charges were used for the metal cations and oxide anions within the metal nodes; we can therefore consider the node and ligands as essentially separate components.

The total internal energy (U) expression can be written as

where, k_r, k_θ and k_Ψ are interatomic force constants, r the distance between a pair of atoms, θ are 3-body and Ψ are 4-body angles, q represents point charges and ε₀ the vacuum permittivity.

Long-range interactions of the ligand were treated as Lennard-Jones functions with combination rules for ε and σ of each individual atom. Lennard-Jones interactions of the ligands were truncated at 12.5 Å and the cut-off of the MM3 Buckingham interactions of the metal node was set to 12.0 Å.

For the TiO₂ surface, a formal charge model was adopted, with Ti and O interacting via a MM3 Buckingham term as was used for the MOF metal node. Combination rules were used for the ε_ij and d⁰_ij values of the MM3 functional, which were parameterised to reproduce structural and mechanical properties of bulk rutile TiO₂.

3 Results and discussion

3.1 Epitaxial matching

A screening procedure was executed to calculate the lattice mismatch between a range of common MOFs and inorganic materials, such as oxides and chalcogenides, with the purpose being to identify systems with small variations of the interfacing surface cell parameters (Fig. 1). The procedure does not consider chemical identity, but focuses on lattice strain, calculated as the percentage difference between interfacing surface parameters. If lattice mismatch of the cleaved surfaces is calculated to exceed 8% between a “soft” MOF and comparatively hard binary material surface, is it unlikely that a uniform coverage of that MOF on the surface would be observed. Furthermore, large lattice mismatch between surface parameters is likely to introduce extended defects at interfacing sites, with mechanical instability and weak chemical bonding.

Fig. 1 Assessment of epitaxial matching between a range of metal–organic frameworks and inorganic materials. The square markers signify lattice mismatch within 8%. The size of the marker is inversely proportional to the mismatch in surface area (u × v expansion) between the interfacing surfaces. The colour of the marker reflects the lattice strain, from blue (low strain) to red (high strain). The most favourable matches therefore appear as large blue squares. The interface with the smallest mismatch is plotted for each combination. The polymorphs chosen are the experimental structures at T = 300 K, with both the anatase (A) and rutile (R) phases considered for TiO₂.

The results of the epitaxial screening are summarised in Fig. 1. Numerical values of lattice mismatch, surface indices and surface expansions are provided as ESI.† Several interesting trends are found. Firstly, common binary materials used as templates for surface growth, including Al₂O₃, TiO₂ and α-SiO₂, are identified to have low lattice mismatch with many MOF topologies. This includes common MOFs such as COF-1, MIL-125, DMOF-1 (orthorhombic and rhombohedral polymorphs) and MOF-649, suggesting these to be functional templating materials. Secondly, we highlight that MOF-649 and DMOF-1 with wine-rack pore topology show the greatest surface compatibility across most binary materials considered. These MOFs are good candidates for porous thin film materials. The screening procedure has also identified ZnO and ZrO₂ to be poor substrate materials for the epitaxial growth of MOFs.

Our results confirm experimental observations made by Hermes et al. who reported the successful growth of MOF-5 on Al₂O₃ surfaces but state that growth was not possible on α-SiO₂ wafers with the same reaction conditions.³⁸ There is a favourable lattice mismatch between the (011) and (110) surfaces of MOF-5 and the (010) and (100) surfaces of Al₂O₃. However, we do not find any favourable lattice planes that would interface effectively between MOF-5 and α-SiO₂ in accord with the experimental observations.

3.2 Epitaxial growth of MOF-5 on TiO₂

We now consider the interface between the (110) surface of rutile TiO₂ and the (011) surface of MOF-5 in more detail as a representative system (see Fig. 2). TiO₂ is a photocatalytic material and popular for its surface reactivity.^39–41 The (110) surface of TiO₂ is known to favourably bind the benzene dicarboxylate (BDC) ligands that form MOF-5.⁴² MOF-5, first synthesised by Yaghi et al., is composed of BDC ligands and Zn²⁺ cations in tetrahedral coordination, with each metal node containing 4 Zn ions and 1 central inorganic oxygen anion.⁴³

Fig. 2 Representation of the unit cell of rutile TiO₂ (left) and MOF-5 (right) with highlighted lattice planes (pink) along (110) and (011), respectively.

3.2.1 Carboxylate attachment on TiO₂: first principles. The (110) surface of TiO₂ was cut with a 6 × 3 surface expansion using GDIS.⁴⁴ The final surface slab contained 540 atoms with cell dimensions of 17.7 and 19.5 Å. The complexity of surface defects and the presence of water adlayers is extensively reported for TiO₂.^45–48 To reduce the cost and complexity of the calculations we consider the bare unreconstructed (110) surface of TiO₂ as a representative interface with MOFs. To create a 3-D periodic model, a vacuum gap of 32 Å was added in the z-direction. Such a large vacuum region is necessary as the length of BDC exceeds 7 Å. Once placed on the oxide surface, the gap must be sufficient that electrostatic interactions of the BDC/oxide surface do not extend across periodic images. Multiple configurations for placing BDC on the (110) surface of TiO₂ were investigated including atomic relaxation (see Fig. 3 and Table 1).

Fig. 3 Initial configurations of the BDC ligand on the (110) surface of TiO₂. Highlighted (purple) is a central region of the TiO₂ surface for depicting the relative positions of Ti polyhedra that the BDC were bound to. For each surface model we highlight (blue) the Ti polyhedra that the BDC is initialised as bound to, with the capping of the carboxylate/carboxylic acid group that was considered at this binding site (lower right of each image). Yellow highlighted oxygen atoms indicate an initial H-bonding interaction with protonated BDC and green shows the initial position of protons when considering a deprotonated ligand at the surface. Arrows across binding sites show the orientation of the BDC ligand.

Table 1 Relative energy (ΔU_R) of the six models considered for binding of BDC on TiO₂ (shown in Fig. 3). The internal energy of binding (ΔU) is calculated with reference to BDC (details in ESI). Note models 5 and 6 include the energy of proton transfer to the surface

Model	ΔUR (eV)	ΔU (eV)
1	1.231	−1.442
2	1.263	−1.410
3	1.202	−1.472
4	Transforms to model 6
5	0.870	−1.804
6	0.000	−2.673

---

Starting from our initial trial structures, we found several configurations that are stationary points on the potential energy surface, as well as the ground-state configuration for BDC on TiO₂. The most favourable site for BDC adsorption (Table 1) is model 6 (Fig. 3) with monodentate binding and deprotonation.

When initialised as monodentate above one Ti surface cation, the configuration is only stable if the proton of the BDC points at a surface oxide anion (as depicted in model 1 in Fig. 3). When initialised as bidentate above two Ti surface cations, several configurations representing stationary points are found. When considering both a protonated and deprotonated BDC ligand in bidentate coordination, the configuration was only stable for the ligand being positioned over face-sharing Ti cations (model 3). Furthermore, when deprotonated, two stable configurations are observed for the ligand where the proton is located on either the 3 or 2 coordinate surface oxide anion (as depicted as model 5 and 6, respectively). The high relative energies of the configurations compared to model 6 (Fig. 3) suggest that for an isolated ligand on the surface of TiO₂, only this ground-state configuration would be accessible. As a final note, when initialised above two edge-sharing Ti cations (model 4), spontaneous deprotonation and rotation of the ligand into the ground-state configuration was observed (as illustrated in Fig. 4).

Fig. 4 Deprotonation of BDC (model 4 in Fig. 3) during optimisation: (a) the initial configuration of the ligand; (b) the loss of the bonding interaction between the protonated carboxylic acid oxygen and Ti as the ligand begins to rotate; (c) the proton on the carboxylic acid oxygen beginning to interact with a bridging oxygen on the surface; (d) the proton transfer from the ligand to the surface.

3.2.2 Carboxylate attachment on TiO₂: forcefield. Exploring the configurational space of the interface between MOF-5 and TiO₂ using first principles techniques is prohibitively expensive from a computational perspective. We therefore developed an analytical forcefield model to describe the interface. The starting point was an existing forcefield derived for metal–organic frameworks, VMOF, which describes the bulk properties of a wide range of metal–organic frameworks.³¹ A forcefield for TiO₂ was fitted to reproduce the lattice parameters and elastic constants of bulk TiO₂ (rutile) the results of which are shown in Table 2.

Table 2 Comparison of structural and mechanical properties of TiO₂ from forcefield (FF) calculations against experimental (Exp.) values. The percentage deviation in the cell parameters is given in parenthesis. The elastic constants are given in GPa, while cell lengths are in Å

Property	Exp.	FF
a	4.587	4.421 (3.61)
b	4.587	4.421 (3.61)
c	2.954	3.068 (3.87)
C11	268.0	235.3
C12	175.0	178.5
C23	147.0	114.6
C33	484.2	403.9
C44	123.8	108.6
C66	190.2	178.5

---

The next step was to model the MOF/oxide interaction, which includes reproducing the first principles binding energies of BDC to the TiO₂ surface. For consistency, experimental lattice parameters of TiO₂ were used for both sets of calculations on the surface, with the atomic positions fully relaxed in each case. Charges in the cross-linking BDC ligand were re-fitted to reproduce the Ti–O(carb) bond lengths and binding energies of the most stable protonated (model 3) and deprotonated (model 6) configurations. The reference states, including a required Coulomb correction, are detailed in the ESI.† When a complete thermodynamic cycle is considered, both the binding energies and proton transfer energies calculated from DFT and FF techniques agree to within 0.1 eV.

3.2.3 MOF-5 on TiO₂: interface structure. The lattice spacing of the MOF-5 epilayer was expanded uniformly to match the more rigid TiO₂ substrate. Note that the calculated bulk moduli are 217.7 GPa for TiO₂ and 8.8 GPa for MOF-5. The uniform expansion in the x and y direction will lead to consistent strain throughout the MOF layer regardless of height of material considered, which may not be consistent with realistic growth, but is unavoidable within periodic boundaries. The atom positions of the vacuum slabs of the protonated oxide and deprotonated MOF were then separately optimised using our forcefield.

The MOF-5 termination consists of two ligands that contact the surface. The ligands are perpendicular, with a difference in rotation of 90° at the surface. Although one ligand (ligand a) was initialised as bidentate between two Ti atoms – as identified as the lowest energy position for the isolated ligand – the other (ligand b) must be initialised as monodentate above one Ti cation. Only the most stable deprotonated carboxylate model was considered.

For an initial model of the interface between the surfaces of MOF-5 and TiO₂, we first consider the minimum layer thickness of MOF-5 to contain one internal pore. During optimisation, ligand b significantly changes geometry. The ligand begins to tilt relative to the surface normal as the carboxylic acid head rotates. The rotation causes electrostatic repulsion between the carboxylic acid oxygen and the protonated bridging oxygen (O(O–H)) of the TiO₂ surface. As a consequence of this repulsion the O(O–H) is displaced, fulfilling the valence of a neighbouring 5-coordinate Ti site, and the carboxylic acid oxygen becomes incorporated into the surface in place of the O(O–H) bridging group. Two configurations were possible depending on the direction that the ligand initially tilts in, as determined by which side of the carboxylic acid oxygen that the proton was initialised on, with one configuration being 0.094 eV more stable than the latter. The two possible models for ligand b as described will be referred to as ligand b model 1 and 2 and are depicted in the ESI.† Ligand a remains in bidentate coordination following optimisation and little change in geometry is observed, further confirming the stability of the ground-state configuration of isolated BDC on TiO₂.

First principles calculations were conducted to verify the predicted reconstruction on three positions of the isolated ligand b on the surface, systems 1–3 correspond to results given in Table 3: (1) ligand b on the surface prior to reconstruction; (2) ligand b model 1 following surface reconstruction; (3) ligand b model 2 following surface reconstruction. Comparison of relative energies (Table 3) prior to and following the reconstruction of the TiO₂ surface, confirm forcefield predictions of the described reconstruction mechanism at the site of ligand b adsorption. Following the reconstruction, the final configuration is the second most energetically favourable position of isolated BDC on the (110) surface of TiO₂, when compared to the identified ground-state configuration (ligand a position).

Table 3 First-principles relative energies prior to and following the predicted reconstruction of the TiO₂ surface at the site of ligand b adsorption. Energies are also given for the three configurations (ΔU) relative to the thermodynamic ground-state configuration (ΔU^ground-state), illustrated as model 6 in Fig. 3

Position	ΔU (eV)	ΔU^ground-state (eV)
1	2.730	3.014
2	0.000	0.284
3	0.037	0.321

---

3.2.4 MOF-5 on TiO₂: interface energetics. Following the identification of the ground-state surface configuration of TiO₂ when interfaced with MOF-5, the strength of interaction between the layers can be assessed. Several models were constructed to simulate the simultaneous growth of the layers of MOF-5 and TiO₂ (Fig. 5). From a plot of layer thickness against total internal energy (see ESI†), the energy of adhesion can be defined. We calculated this to be 1.43 J m⁻² for MOF-5 and TiO₂.

Fig. 5 Structure models considered for different layer thickness of MOF-5 and TiO₂. Ligand a and ligand b are labelled for each model. Note that each structure is periodic in the xy plane.

To put 1.43 J m⁻² into perspective, values of between 2 and 5 J m⁻² are typical for metal/metal oxide interfaces, whilst values less than 1 J m⁻² are expected for weakly bound interfaces dominated by van der Waals interactions. For example, Kohyama et al.⁴⁹ calculate the energy of adhesion between Si-terminated SiC with Ti and Al metal surfaces to be 2.52 and 3.74 J m⁻², respectively. In contrast, the interface between multi-layered graphene and SiO₂ has an energy of adhesion of 0.44 J m⁻².⁵⁰ The value we report for a MOF/oxide interface therefore suggests a reasonably strong chemical binding.

3.2.5 MOF-5 on TiO₂: proton distribution. To further validate the identified ground-state configurations of the interfacing surfaces, the possibility of proton migration across neighbouring oxide sites was considered for the model of interfacing MOF-5 and TiO₂. An extensive analysis of different proton positions was conducted by migrating the protons along neighbouring oxide bridges of the TiO₂ surfaces and comparing relative energies. Owing to the expense of the size of the systems considered, relative energies were calculated with the parameterised forcefield, which can reproduce the binding energy of the proton to the surface as previously predicted by first principles calculations. The movement of the protons from ligand a and b were considered separately for two different structure models with configuration descriptions and relative energies reported in the ESI.†

We calculate for both ligand a and b that the lowest energy positions of the proton (following deprotonation of the BDC ligands), is to remain on a neighbouring oxide site. Fig. 6 highlights the oxide positions surrounding ligand a and b where proton movement from the identified ground-state configurations costs <0.1 eV. Interestingly, at the site of ligand b a further reconstruction of the surface is observed when substituting in a neighbouring site as highlighted (Fig. 6), resulting in the newly protonated oxygen being displaced off-site and tilt of the ligand changing as it rotates into the new position at the surface in its place. The thickness of the MOF layer (as studied by considering the difference in calculated energies for structures 1 and 2) was found not to affect the binding mechanism of MOF-5 on the surface of TiO₂ (see ESI†). The same ground-state configuration is maintained for all layer thicknesses considered in this study.

Fig. 6 The (110) surface of TiO₂ with highlighted oxygen sites where proton migration from ligand a (green) and ligand b (purple) costs <0.1 eV as predicted by forcefield calculations. TiO₄ polyhedra are highlighted if ligand a (green) or ligand b (purple) remain bonding to the Ti sites following optimisation. At the site of ligand b, an oxygen atom is highlighted (yellow), which belongs to carboxylate group of ligand b and has become incorporated into the surface following reconstruction. Ligand orientation (white line) is given for both ligand a and b. The full structure of the MOF has been removed for clarity.

4 Conclusions

A number of viable metal–organic framework/inorganic substrate combinations have been identified that could be used for epitaxial thin-film growth. We have focused on the prototype case of the (011) surface of MOF-5 on the (110) surface of rutile TiO₂, including the initial contact of the ligands, followed by the structure and thermodynamics of complete films. Ideal growth with clean termination of surfaces were considered for this initial analysis of the thermodynamics associated with the growth of MOFs on oxide surfaces. We find that following deprotonation of the BDC ligand, it is thermodynamically favourable for the TiO₂ surface to reconstruct, resulting in the incorporation of a ligand into the surface. The energy of adhesion between MOF-5 and TiO₂ surfaces is calculated to be 1.43 J m⁻², which reflects the significant chemical interaction at the hybrid heterointerface.

5 Data access statement

The VMOF forcefield is available as a library file for GULP from https://github.com/WMD-group/VMOF. The interface structure models and program input/output files are available from https://github.com/WMD-group/MOF-Epitaxy and https://researchdata.ands.org.au.

Acknowledgements

J. K. B. is funded by the EPSRC (grant no. EP/G03768X/1). J. D. G. thanks the Australian Research Council for funding under the Discovery Programme, as well as the Pawsey Supercomputing Centre and NCI for the provision of computing resources. A. W. acknowledges support from the Royal Society University Research Fellowship scheme. K. L. S. is funded under ERC Starting Grant 277757 and K. T. B. is funded under EPSRC grant no. EP/M009580/1. The work benefited from the high performance computing facility at the University of Bath. Access to the ARCHER supercomputer was facilitated through membership of the HPC Materials Chemistry Consortium (EP/L000202).

References

1. O. Shekhah, H. Wang, D. Zacher, R. A. Fischer and C. Wöll, Angew. Chem., Int. Ed., 2009, 48, 5038–5041.
2. L. K. Cadman, J. K. Bristow, N. E. Stubbs, D. Tiana, M. F. Mahon, A. Walsh and A. D. Burrows, Dalton Trans., 2016, 45, 4316–4326.
3. N. Stock and S. Biswas, Chem. Rev., 2011, 112, 933–969.
4. D. Zacher, O. Shekhah, C. Wöll and R. A. Fischer, Chem. Soc. Rev., 2009, 38, 1418–1429.
5. O. Shekhah, Materials, 2010, 3, 1302–1315.
6. D. Zacher, R. Schmid, C. Wöll and R. A. Fischer, Angew. Chem., Int. Ed., 2011, 50, 176–199.
7. S. Hermes, F. Schröder, R. Chelmowski, C. Wöll and R. A. Fischer, J. Am. Chem. Soc., 2005, 127, 13744–13745.
8. B. Liu, M. Ma, D. Zacher, A. Bétard, K. Yusenko, N. Metzler-Nolte, C. Wöll and R. A. Fischer, J. Am. Chem. Soc., 2011, 133, 1734–1737.
9. E. Biemmi, C. Scherb and T. Bein, J. Am. Chem. Soc., 2007, 129, 8054–8055.
10. C. H. Hendon and A. Walsh, Chem. Sci., 2015, 6, 3674–3683.
11. R. Ameloot, L. Stappers, J. Fransaer, L. Alaerts, B. F. Sels and D. E. De Vos, Chem. Mater., 2009, 21, 2580–2582.
12. Y. Yoo and H.-K. Jeong, Cryst. Growth Des., 2010, 10, 1283–1288.
13. M. Tu and R. A. Fischer, J. Mater. Chem. A, 2014, 2, 2018–2022.
14. L. D. Salmi, M. J. Heikkilä, E. Puukilainen, T. Sajavaara, D. Grosso and M. Ritala, Microporous Mesoporous Mater., 2013, 182, 147–154.
15. I. Stassen, M. Styles, G. Grenci, H. Van Gorp, W. Vanderlinden, S. De Feyter, P. Falcaro, D. De Vos, P. Vereecken and R. Ameloot, Nat. Mater., 2015, 15, 304–310.
16. M. Shöâeè, J. R. Agger, M. W. Anderson and M. P. Attfield, CrystEngComm, 2008, 10, 646–648.
17. S. Amirjalayer, M. Tafipolsky and R. Schmid, J. Phys. Chem. Lett., 2014, 5, 3206–3210.
18. J. K. Bristow, K. L. Svane, D. Tiana, J. M. Skelton, J. D. Gale and A. Walsh, J. Phys. Chem. C, 2016, 120, 9276–9281.
19. K. T. Butler, Y. Kumagai, F. Oba and A. Walsh, J. Mater. Chem. C, 2016, 4, 1149–1158.
20. S. R. Bahn and K. W. Jacobsen, Comput. Sci. Eng., 2002, 4, 56–66.
21. A. Zur and T. McGill, J. Appl. Phys., 1984, 55, 378–386.
22. A. Zur, T. McGill and M.-A. Nicolet, J. Appl. Phys., 1985, 57, 600–603.
23. J. VandeVondele, M. Krack, F. Mohamed, M. Parrinello, T. Chassaing and J. Hutter, Comput. Phys. Commun., 2005, 167, 103–128.
24. J. P. Perdew, A. Ruzsinszky, G. I. Csonka, O. A. Vydrov, G. E. Scuseria, L. A. Constantin, X. Zhou and K. Burke, Phys. Rev. Lett., 2008, 100, 136406.
25. S. Grimme, J. Antony, S. Ehrlich and H. Krieg, J. Chem. Phys., 2010, 132, 154104.
26. M. Krack, Theor. Chem. Acc., 2005, 114, 145–152.
27. K. A. Peterson and T. H. Dunning Jr, J. Chem. Phys., 2002, 117, 10548–10560.
28. T. H. Dunning Jr, J. Chem. Phys., 1970, 53, 2823–2833.
29. J. D. Gale, J. Chem. Soc., Faraday Trans., 1997, 93, 629–637.
30. J. D. Gale and A. L. Rohl, Mol. Simul., 2003, 29, 291–341.
31. J. K. Bristow, J. M. Skelton, K. L. Svane, A. Walsh and J. D. Gale, Phys. Chem. Chem. Phys., 2016, 18, 29316–29329.
32. J. K. Bristow, D. Tiana and A. Walsh, J. Chem. Theory Comput., 2014, 10, 4644–4652.
33. N. L. Allinger, Y. H. Yuh and J. H. Lii, J. Am. Chem. Soc., 1989, 111, 8551–8566.
34. N. L. Allinger, F. Li and L. Yan, J. Comput. Chem., 1990, 11, 848–867.
35. B. R. Brooks, C. L. Brooks, A. D. MacKerell, L. Nilsson, R. J. Petrella, B. Roux, Y. Won, G. Archontis, C. Bartels and S. Boresch, J. Comput. Chem., 2009, 30, 1545–1614.
36. J. Gasteiger and M. Marsili, Tetrahedron, 1980, 36, 3219–3228.
37. J. Gasteiger and M. Marsili, Tetrahedron Lett., 1978, 19, 3181–3184.
38. S. Hermes, D. Zacher, A. Baunemann, C. Wöll and R. A. Fischer, Chem. Mater., 2007, 19, 2168–2173.
39. N. Negishi, K. Takeuchi and T. Ibusuki, J. Mater. Sci., 1998, 33, 5789–5794.
40. D. O. Scanlon, C. W. Dunnill, J. Buckeridge, S. A. Shevlin, A. J. Logsdail, S. M. Woodley, C. R. A. Catlow, M. J. Powell, R. G. Palgrave and I. P. Parkin, et al. , Nat. Mater., 2013, 12, 798–801.
41. J. Buckeridge, K. T. Butler, C. R. A. Catlow, A. J. Logsdail, D. O. Scanlon, S. A. Shevlin, S. M. Woodley, A. A. Sokol and A. Walsh, Chem. Mater., 2015, 27, 3844–3851.
42. A. Tekiel, J. S. Prauzner-Bechcicki, S. Godlewski, J. Budzioch and M. Szymonski, J. Phys. Chem. C, 2008, 112, 12606–12609.
43. S. S. Kaye, A. Dailly, O. M. Yaghi and J. R. Long, J. Am. Chem. Soc., 2007, 129, 14176–14177.
44. S. Fleming and A. Rohl, Z. Kristallogr., 2005, 220, 580–584.
45. C. Zhang and P. J. Lindan, J. Chem. Phys., 2003, 118, 4620–4630.
46. J.-M. Pan, B. Maschhoff, U. Diebold and T. Madey, J. Vac. Sci. Technol., A, 1992, 10, 2470–2476.
47. K. Onda, B. Li, J. Zhao, K. D. Jordan, J. Yang and H. Petek, Science, 2005, 308, 1154–1158.
48. B. J. Morgan and G. W. Watson, J. Phys. Chem. C, 2009, 113, 7322–7328.
49. M. Kohyama and J. Hoekstra, Phys. Rev. B: Condens. Matter Mater. Phys., 2000, 61, 2672.
50. S. P. Koenig, N. G. Boddeti, M. L. Dunn and J. S. Bunch, Nat. Nanotechnol., 2011, 6, 543–546.

---

Footnote

† Electronic supplementary information (ESI) available. See DOI: 10.1039/c7ta00356k

---