Linear-scaling density functional simulations of the effect of crystallographic structure on the electronic and optical properties of fullerene solvates

Abstract

In this work, the crystal properties, HOMO and LUMO energies, band gaps, density of states, as well as the optical absorption spectra of fullerene C₆₀ and its derivative phenyl-C₆₁-butyric-acid-methyl-ester (PCBM) co-crystallised with various solvents such as benzene, biphenyl, cyclohexane, and chlorobenzene were investigated computationally using linear-scaling density functional theory with plane waves as implemented in the ONETEP program. Such solvates are useful materials as electron acceptors for organic photovoltaic (OPV) devices. We found that the fullerene parts contained in the solvates are unstable without solvents, and the interactions between fullerene and solvent molecules in C₆₀ and PCBM solvates make a significant contribution to the cohesive energies of solvates, indicating that solvent molecules are essential to keep C₆₀ and PCBM solvates stable. Both the band gap (E_g) and the HOMO and LUMO states of C₆₀ and PCBM solvates are mainly determined by the fullerene parts contained in solvates. Chlorobenzene- and ortho-dichlorobenzene-solvated PCBM are the most promising electron-accepting materials among these solvates for increasing the driving force for charge separation in OPVs due to their relatively high LUMO energies. The UV-Vis absorption spectra of solvent-free C₆₀ and PCBM crystals in the present work are similar to those of C₆₀ and PCBM thin films shown in the literature. Changes in the absorption spectra of C₆₀ solvates relative to the solvent-free C₆₀ crystal are more significant than those of PCBM solvates due to the weaker effect of solvents on the π-stacking interactions between fullerene molecules in the latter solvates. The main absorptions for all C₆₀ and PCBM crystals are located in the ultraviolet (UV) region.

---

1 Introduction

Conjugated polymer/fullerene organic solar cells (OSCs) have high potential as low-cost but efficient renewable energy sources and have attracted increasing attention in recent years, as they are relatively lighter, cheaper and easier to produce from renewable sources in contrast with most inorganic solar cells.^1–6 Today, power conversion efficiencies (PCEs) of over 10% for solar cells based on blend films of conjugated polymers with fullerenes have been reached.^7,8 In these high performance devices, fullerenes are used as electron-accepting materials. Furthermore, it has been found that fullerenes crystallize in some solid blends prepared from organic solvents such as benzene, chlorobenzene (CB), ortho-dichlorobenzene (oDCB), and chloroform at high concentration or following annealing treatment.^9–23

The crystallite formation of fullerenes can have a significant effect on the properties of blend films because the properties of crystalline fullerenes may be different from those of amorphous fullerene molecules. For instance, the band gaps (E_g) of both C₆₀ and its derivative phenyl-C₆₁-butyric-acid-methyl-ester (PCBM) crystals are lower when compared to those of the isolated C₆₀ and PCBM molecules.²⁴ This E_g reduction has a significant effect on their optical properties. Moreover, the large crystallite formation of fullerenes in some systems may aid charge mobility relative to the amorphous phase due to the molecular ordering and the lack of grain boundaries.^25,26 Some previous studies^{14,16,17,20,21,23} have also shown that devices based on poly-3-hexylthiophene (P3HT) and PCBM blends display a better performance (e.g., a large increase in photocurrent and carrier lifetime) following thermal annealing which alters the fraction of amorphous and microcrystalline PCBM domains, with the formation of PCBM nano-crystalline domains that are crucial for high PCEs.^14,16,27

To date, a large number of studies both in experimental implementations^{9–23,27–43} and in theoretical aspects^24,44–49 have been carried out on the preparation, characterization of morphology and structure and physical properties of solvent-free or solvated C₆₀ and PCBM crystals. All these studies can not only help to improve the performance of these crystals in OSCs but can also help to provide a deeper insight into their fundamental properties. It is found that both solvent-free C₆₀ and PCBM crystals crystallize with four molecules, in space group Pa3̄ of the cubic system^34,35 and P2₁/c of the monoclinic system^37,38 respectively [see Fig. 1(a) and (e)], while the crystal structure of solvated C₆₀ and PCBM co-crystals strongly depends on the solvent because it can diffuse into the crystal lattice and play an important role in modifying the interactions between the molecules; thus fullerenes can assume different arrangements when using different solvents. For instance, Meidine et al.³³ found that the crystal structure of benzene-solvated C₆₀ is triclinic with a P1̄ space group symmetry where the ordered molecules of C₆₀ are in an approximately hexagonal close-packed arrangement and separated by benzene molecules. Pénicaud et al.³⁰ yielded some single crystals of the formula [(C₆₀)][(C₆H₅)₂] with a monoclinic structure by slow co-crystallization of saturated solutions of C₆₀ and biphenyl in toluene at ambient temperature. Rispens et al.²⁷ obtained two PCBM packing structures from single crystals grown from CB and oDCB by using extremely slow solvent evaporation conditions and maintaining a solvent saturated atmosphere. When PCBM crystallizes from CB, the co-crystal structure is triclinic with four PCBM and two solvent molecules in each unit cell, while a monoclinic structure with four PCBM and four solvent molecules in the unit cell is formed when using oDCB (the asymmetric units of these crystal structures are shown in Fig. 1). From a structure–property relationship point of view, different crystal structures may result in different electronic and optical properties. Both of these are very important for the optimization of solar cell devices, in particular since the maximum open-circuit voltage (V_oc), one of the most important parameters for the solar cell efficiency, largely depends on the gap between the highest occupied molecular orbital (HOMO) energy of the donor material and the lowest unoccupied molecular orbital (LUMO) of the acceptor material.^50,51 In addition, the photocurrent is also determined to some extent by the absorption coefficient of the photoactive layer.⁵² In contrast to the case of solvent-free C₆₀ and PCBM crystals, whose electronic and optical properties have been widely investigated,^31,40–49 there is little information available on the electronic and optical properties of solvated fullerene crystals so far, to the best of our knowledge, and therefore, in order to tailor the efficiency of polymer/fullerene based solar cells in a suitable way for higher PCEs, it is necessary to study the electronic and optical properties of fullerenes co-crystallised with different solvents. To this end, a suitable approach is to use simulations with first-principles density functional theory (DFT) that can provide accurate ground state properties such as energy eigenstates and also allows the calculation of optical absorption spectra. Conventional first-principles DFT codes such as CASTEP, VASP, PWSCF and ABINIT require a computational effort that scales cubically with the model size. As the solvated co-crystals we are interested in here are large in scale, involving as many as 216 to 2208 atoms, calculations using these codes can be too computationally expensive. Hence, we have employed a linear-scaling (LS) DFT method as implemented in the ONETEP (Order-N Electronic Total Energy Package) program⁵³ which is particularly suitable for studying large systems, since the computational cost only increases linearly with respect to the number of atoms.

Fig. 1 Unit cells for pristine C₆₀ and PCBM crystals and their solvates. Fullerenes are shown in a CPK representation with solvent molecules drawn as bonds.

2 Computational details

All structural properties of fullerene crystals such as the space group and lattice parameters in this paper (summarised in Table 1) are taken from the Cambridge Structural Database (CSD), except for the 1-methylnaphthalene-solvated PCBM co-crystal whose crystallographic structure is reported herein for the first time (see ESI†). One should note that these structures may not be the only possible ones since fullerene crystals usually can have a number of polymorphs close in energy, and there is no guarantee that these polymorphs whose X-ray structures are solved in bulk powder samples are prevalent or even present in films. Although this is not the focus of this work, in order to predict the possible or the most stable structures of these fullerene crystals, computational methods for organic crystal structure prediction based on global minimisation of the lattice energy⁵⁴ can be employed, since they have been successfully applied to the prediction of single-⁵⁵ and multi-component⁵⁶ organic molecular crystals and their reliability has been tremendously improved in recent years.⁵⁷ However, XRD is the most reliable method for the experimental determination of crystal structures and the R₁ and wR₂ refinement indices for the structures we considered are within accepted limits.

Table 1 Refcode (CSD entry identifier), composition, space group, lattice parameters and density (ρ) of the C₆₀ and PCBM crystals as well as the temperature (T) at which these parameters are determined

Refcode	Composition (per unit cell)	Space group	a (Å)	b (Å)	c (Å)	α (°)	β (°)	γ (°)	ρ (g cm^−3)	T (K)
SOCTOT03	4C60 (fullerene)	Pa3̄ (cubic)	14.06	14.06	14.06	90.00	90.00	90.00	1.722	153
JUGCET	2C60, 8C6H6 (benzene)	P1̄ (triclinic)	9.96	15.07	17.50	65.26	88.36	74.94	1.495	173
MULFOO	4C60, 4C12H10 (biphenyl)	I2/a (monoclinic)	10.40	16.91	20.46	90.00	100.87	90.00	1.644	155
YOLSOH	8C60, 96C6H12 (cyclohexane)	F4̄3c (cubic)	28.17	28.17	28.17	90.00	90.00	90.00	1.029	RT
PESJII01	4PCBM (fullerene)	P21/c (monoclinic)	13.47	15.14	19.01	90.00	106.90	90.00	1.621	100
—	4PCBM, 4C11H10 (1-methylnaphthalene)	P21/c (monoclinic)	18.83	13.27	18.43	90.00	105.96	90.00	1.580	100
EKOZUZ	4PCBM, 2C6H5Cl (chlorobenzene)	P1̄ (triclinic)	13.83	15.29	19.25	80.26	78.56	80.41	1.649	90
EKOZOT	4PCBM, 4C6H4Cl2 (ortho-dichlorobenzene)	P21/c (monoclinic)	13.76	16.63	19.08	90.00	105.29	90.00	1.669	90

---

During the geometry optimization we kept the dimensions of the unit cell as well as the atomic positions fixed, with the only exception of the hydrogen atoms, since it is not usually possible to accurately determine in experiments their position.

As is well known, the virtual orbitals in the exact Kohn–Sham (KS) model of DFT are genuinely bound one-electron states at exactly the same (local) potential as the occupied orbitals, and represent excited electrons in the non-interacting KS reference system.⁵⁸ Their energies, under certain conditions, can be interpreted as ionization energies from the corresponding states since the occupied KS orbital energies have already been interpreted as approximate ionization energies.⁵⁹ As a result, we can easily interpret the excitation spectrum using the one-electron KS model. In contrast, according to the canonical Koopmans' theorem, the energies of virtual orbitals in the Hartree–Fock (HF) model are interpreted as approximate electron affinities to an extra electron based on the frozen-orbital HF approximation.⁶⁰ Moreover, in many cases the HF virtual orbitals are not bound states and are shifted to much higher energies than the KS virtual orbitals, with their shapes being very diffuse.^58,60 Therefore, the HF virtual orbitals are often unphysical. The ONETEP program, which employs the one-particle density matrix reformulation of KS DFT, can describe the occupied (valence) KS states very well as compared to the conventional DFT implementations. However, it is unable to represent the unoccupied (conduction) equally well because the local orbitals used in ONETEP are specifically optimised to describe the valence states. In ONETEP calculations, these localized orbitals are represented by a set of non-orthogonal generalized Wannier functions (NGWFs)⁶¹ which are atom centered and strictly localized within a set radius. The NGWFs are represented by a basis set of periodic cardinal sinc (psinc) functions,⁶² which is equivalent to a plane-wave basis set, and is optimized self-consistently during the calculations. Due to the self-consistent optimisation of the NGWFs and the psinc basis set, the distinguishing feature of ONETEP is that it is able to perform linear-scaling DFT calculations⁶³ with near-complete basis set accuracy, as is possible in conventional cubic-scaling DFT approaches.

Typically, the unoccupied states that lie 1–2 eV above the Fermi level obtained from the NGWFs in a standard ONETEP calculation are higher in energy than those obtained from a conventional cubic-scaling DFT code such as CASTEP,⁶⁴ and some unoccupied states are even lost.⁶⁵ Therefore, a recently developed methodology in which a second set of localized orbitals (conduction NGWFs) is optimized to accurately describe the unoccupied states has been employed. Further details of this new methodology have been given elsewhere.^65–67 Once the optimization of the set of conduction NGWFs has finished, the valence and conduction NGWF basis sets can be combined into a new joint basis set which is able to well represent both the occupied and unoccupied KS states of the system. The optical matrix elements of the electronic transitions between valence and conduction states can then be calculated using this new joint basis using Fermi's golden rule, which has been implemented in ONETEP. The electronic and optical properties obtained using the new joint basis set are expected to agree well with those obtained using the conventional cubic-scaling DFT formulation.

Using the dipole approximation in which the exponential is expanded in a Taylor series and only first order terms are included, the imaginary part ε₂(ω) of the frequency-dependent dielectric function is defined as,

where Ω is the cell volume, the indices v and c refer to valence and conduction bands respectively,|ψⁿ_k〉 and Eⁿ_k are the nth eigenstate and corresponding energy at a given k-point, and q̂ is the polarization direction of the photon with an energy of ℏω. Using ONETEP the imaginary part of the dielectric function of many systems such as the metal-free phthalocyanine and C₆₀-conjugated polymer hybrids has been successfully calculated.^65,66 One can then also obtain the real part ε₁(ω) of the dielectric function from ε₂(ω) using the Kramers–Kronig relation as follows,where P denotes the principal value, and it is generally set to 1 in calculations.

The optical absorption coefficient α(ω) indicates how far light of a particular wavelength can penetrate into a material before it is absorbed. It can be calculated directly from the dielectric function by using the following formula:

where c is the speed of light.

All calculations, including geometry optimizations, ground states and conduction states optimizations, have been performed using a dispersion–corrected PBE exchange–correlation (PBE-D2) functional which is based on the correction according to Grimme,⁶⁸ with Γ point sampling only and periodic boundary conditions. The reason why Grimme's DFT-D2 was used is that D2 dispersion correction, which only includes two-body energies and neglects many-body terms and any effects originating from an atom's environment, has been widely used due to its very low computational cost.^69,70 In addition, since many results in the available literature have also been obtained by using PBE-D2,⁴⁷ the use of the same dispersion-corrected exchange–correlation functional would make our calculated results more comparable to them. Moreover, the D3 correction is not currently implemented in ONETEP.

The projector augmented wave (PAW) method was used in geometry optimizations due to its many advantages such as the higher computational efficiency, all-electron precision, and a lower psinc kinetic energy cutoff to converge compared with norm-conserving pseudopotentials (NCPPs). However, the NCPPs were employed in ground-state and conduction calculations because the PAW method is not currently available for computing ε₂(ω) during conduction calculations in ONETEP. The total energies were converged with respect to the psinc cut-off energy and NGWF radii. The psinc grid spacing was set to be equivalent to a kinetic energy cutoff of 800 eV for geometry optimizations and 1000 eV for ground-state and conduction calculations and no truncation was applied to the density kernel. Both geometry optimizations and ground-state calculations were performed using 1 NGWF per H atom and 4 NGWFs each per C, O and Cl atom with a fixed radius of 8.0 Bohr. Conduction calculations were performed with 5 NGWFs per H atom and 13 NGWFs each per C, O and Cl atom, and a fixed radius of 9.0 Bohr. The number of conduction states explicitly optimized was increased linearly with the number of atoms in the system. A Gaussian smearing width of 0.1 eV was used for the calculations of density of states (DOS) and for the dielectric function. The imaginary part of the dielectric function was calculated using the momentum operator formalism since the position operator is ill-defined under periodic boundary conditions.

3 Results and discussion

3.1. Crystal properties

C₆₀ substitution and solvates. Fullerenes are highly symmetric molecules: their spherical shape aids in obtaining an ordered packing structure, leading to the formation of crystals in the solid state. Chemical substitution of a fullerene molecule breaks the symmetry often producing “awkward” molecular structures which present difficulties in organized and symmetrical solid-state packing. Although many examples exist where crystals of C₆₀ have been produced including many solvated crystals, few examples of crystals have been reported for substituted fullerenes, specifically those commonly used in molecular electronics such as the ubiquitous PCBM.

Growth of single crystals of PCBM was attempted using slow crystallization conditions to maximize our chances of success. We employed slow vapour diffusion of isopropanol into a concentrated solution of 1-methylnaphthelene (1-MN). Various concentrations of PCBM in 1-methylnapthalene were explored in order to balance the rate of crystal formation to achieve diffractable single crystals. Attempts to grow single crystals from tetralin and ortho-xylene did not yield success.

Non-planar molecules are known to promote solvated structures.^71–73 This is presumably due to the significant void space commonly observed in crystals formed from non-planar molecules. The void volume and void fraction were calculated for the studied fullerene crystals using Crystal Explorer 3.1 (isovalue: 0.002, quality: high), as shown in Table 2. In each case, lattice inclusion of solvates reduced the void fraction of the crystals in comparison to the non-solvated crystals, the only exception being cyclohexane. However, a closer inspection of the cyclohexane solvate crystal revealed that the increased void space in this case is due to the inefficient molecular packing of the cyclohexane layer formed between fullerene layers within the crystal and not due to an increased void space within the fullerene layers (See Fig. S1, ESI†).

Table 2 Calculated void volumes and void fractions for the studied fullerene crystals

Crystal	Solvent	Void volume (Å^3)	Volume of unit cell (Å^3)	Void volume fraction (%)	Solvent-solute close interaction	Close interaction distance (Å)
SOCTOT03	None	466.50	2780.03	16.78	None	—
JUGCET	Benzene	278.64	2294.84	12.14	π⋯C60	3.25
MULFOO	Biphenyl	427.13	3533.97	12.09	π⋯C60	3.48
YOLSOH	Cyclohexane	5716.50	22347.13	25.58	None	—
PESJII01	None	451.69	3708.70	12.18	None	—
PCBM-1-MN	1-Methylnaphthalene	478.69	4426.36	10.82	π⋯C60	3.20
EKOZUZ	CB	420.51	3894.90	10.80	None	—
EKOZOT	oDCB	402.15	4210.81	9.55	π⋯C60	3.29

---

Void filling may not be the only driving force for solvate inclusion within the crystal lattice. It appears that solvent–solute interactions are prevalent in the crystal structure of both C₆₀ and PCBM in combination with aromatic solvents. The host–guest inclusions of aromatic molecules display close contact between the π centroid of solvents and the fullerene cage (π⋯C₆₀). Fullerene, being an electron acceptor and π-acidic, demonstrates a closer molecular contact to π-basic solvents (1-MN, benzene) than π-acidic solvents (oDCB). A deeper understanding of the strength of these host–guest interactions can be explained by studying their cohesive energies.

Cohesive energies of solvent-free and solvated C₆₀ and PCBM crystals. As shown in Table 3, the cohesive energy E^coh(total) of a fullerene crystal was calculated as E^coh(total) = −(1/m)[E_cryst − mE_mol(solvent) − nE_mol(fullerene)], where E_cryst is the total energy per unit cell of a fullerene crystal under periodic boundary conditions, E_mol(solvent) and E_mol(fullerene) are the energies of an isolated solvent and fullerene molecules, respectively, and m and n are the number of solvent and fullerene molecules contained in an unit cell. In order to study the contributions of the fullerene–fullerene, solvent–solvent and fullerene–solvent interactions to the cohesive energies of fullerene crystals, we have also calculated the cohesive energies of the fullerene and solvent parts in fullerene crystals as they represent the fullerene–fullerene and solvent–solvent interactions, respectively. Here, the E_cryst of the fullerene (solvent) parts were calculated by removing the solvent (fullerene) molecules contained in the unit cells of crystals and by keeping the same lattice constants and atomic positions as their parent unit cells. The contributions of the fullerene–solvent interactions to the cohesive energies of fullerene crystals were calculated as E^coh(fullerene–solvent) = E^coh(total) − E^coh(fullerene–fullerene) − E^coh(solvent–solvent).

Table 3 Cell volumes and cohesive energies E^coh(total), E^coh(fullerene–fullerene) and E^coh(solvent–solvent) for C₆₀ and PCBM crystals, pure fullerene and solvent parts contained in fullerene crystals compared with the available experimental and theoretical results in parentheses; values in brackets are percentage ratios of E^coh(fullerene–fullerene), E^coh(solvent–solvent) and E^coh(fullerene–solvent) to E^coh(total)

Fullerene	Solvent	Cell volume (Å^3/fullerene molecule)	E ^coh (eV/fullerene molecule)
			Total	Fullerene–fullerene	Solvent–solvent	Fullerene–solvent
a Ref. 74, experiment. b Ref. 46, LDA. c Ref. 47, PBE-D2.
C60	None	695.01	1.64	1.64 (1.74,^a 1.60^b)	0.00	0.00
	Benzene	1147.42	2.65	0.48	0.43	1.74
				[18.09%]	[16.30%]	[65.61%]
	Biphenyl	883.49	2.22	0.75	0.08	1.39
				[33.72%]	[3.67%]	[62.61%]
	Cyclohexane	2793.39	4.94	0.02	2.52	2.40
				[0.36%]	[51.09%]	[48.55%]
PCBM	None	927.18	1.90	1.90	0.00	0.00
	1-Methylnaphthalene	1106.59	2.91	1.49	0.16	1.26
				[51.22%]	[5.37%]	[43.40%]
	Chlorobenzene	973.72	2.34	1.83 (1.91^c)	0.00	0.51
				[78.37%]	[0.00%]	[21.61%]
	ortho-Dichlorobenzene	1052.70	2.67	1.62 (1.76^c)	0.12	0.93
				[60.58%]	[4.67%]	[34.75%]

---

From Table 3 it is possible to observe that the cohesive energies E^coh (fullerene–fullerene) of the solvent-free C₆₀ crystal and PCBM parts in CB- and oDCB-solvated PCBM co-crystals compare quite well with the experimental and theoretical values in the literature.^46,47,74 The E^coh(fullerene–fullerene) of the PCBM part in the CB-solvated PCBM co-crystal is higher than the corresponding one in the oDCB-solvated PCBM co-crystal, which is also in good agreement with previous calculations using the PBE-D2 functional.⁴⁷ We find that the E^coh(fullerene–fullerene) of pure fullerene parts increases with the decrease in the cell volume. As the solvent-free C₆₀ and PCBM crystals do not contain any solvent molecules, they are closest-packed and their cell volumes are the smallest among pure C₆₀ and PCBM parts, respectively. Their cohesive energies are the highest because the cell volume (per fullerene molecule) can represent to some extent the average distance between fullerene molecules under periodic boundary conditions for these pure fullerene parts, and therefore the average distance will decrease as the cell volume per fullerene molecule decreases, while the interactions between fullerene molecules become stronger. The cohesive energy, which is the energy required to separate a solid into its elementary “building blocks”, can become higher within a certain distance range, and this can also explain the computational results reported by Gajdos et al.:⁴⁷ the cohesive energy of different PCBM crystal structures increases in the order body-centered-cubic < hexagonal < simple cubic < monoclinic < triclinic, with the cell volume of these structures decreasing in the same order.

From the cohesive energies E^coh (fullerene–fullerene) in Table 3 we can also conclude that the pure fullerene parts originated from fullerene solvates are more likely to decompose due to their lower cohesive energies, in contrast with the solvent-free fullerene crystals that are able to exist stably under realistic conditions; thus these pure fullerene parts would be metastable or unstable, while the fullerene solvates are stable as observed in experiments. This indicates that the solvents contained in the fullerene solvates are essential for the stability of these solvates. Due to the existence of solvent molecules, the π-stacking interactions between fullerene and solvent molecules can be formed in the fullerene solvates. From the cohesive energies E^coh (total), E^coh (fullerene–fullerene) and E^coh (solvent–solvent) in Table 3 it can be seen that the cohesive energies E^coh (total) of these C₆₀ and PCBM crystals are mainly determined by the fullerene–fullerene and fullerene–solvent interactions therein, except for the cyclohexane-solvated C₆₀ whose E^coh (total) is mainly determined by the interactions between solvent–solvent molecules and fullerene–solvent molecules, as the number of solvent molecules is much higher than that in other fullerene solvates. To sum up, the contributions of the interactions between fullerene–solvent molecules to the cohesive energies E^coh (total) of C₆₀ and PCBM solvates are significant in particular for the benzene- and biphenyl-solvated C₆₀, and this implies that solvent molecules are indispensable to fullerene solvates.

3.2 Electronic properties

HOMO, LUMO energies and energy gaps of solvent-free and solvated fullerene crystals. From Fig. 2 and the densities of C₆₀ and PCBM crystals shown in Table 1, we found that both the HOMO and LUMO energies for C₆₀ and PCBM crystals shift upward as the crystal density increases, which is in agreement with the experimental results that suggest that higher HOMO and LUMO levels of the PCBM film could be linked with the increased film density.⁴¹ As can be seen in Fig. 2, both the HOMO and LUMO energies of C₆₀ solvates are lower than those of the solvent-free C₆₀ crystal, and a similar behaviour can also be seen in the 1-methylnaphthalene-solvated PCBM co-crystal; however the HOMO and LUMO energies of CB- and oDCB-solvated PCBM co-crystals are higher compared to the solvent-free PCBM crystal energies. As a result, CB- and oDCB-solvated PCBM co-crystals should be the most suitable electron acceptor materials to increase the driving force of charge separation in solar cells due to their higher LUMO energies. In addition, the magnitude of the changes in the HOMO and LUMO energies of C₆₀ solvates relative to the solvent-free C₆₀ crystal is clearly larger than that of PCBM solvates compared to the solvent-free PCBM crystal. This is probably because the π-stacking interactions between C₆₀ cages in C₆₀ solvates are largely reduced relative to the strong π-stacking interactions in the pristine C₆₀ crystal, as the existence of solvents therein separates the C₆₀ cages and enlarges the distances between the cages. As shown in Table 3, this kind of interaction only takes up a small portion (less than 50%) of all the interactions in C₆₀ solvates. In contrast, the solvents contained in PCBM solvates do not seem to significantly reduce the π-stacking interactions between PCBM molecules, with these being the strongest among the interactions in PCBM solvates (see Table 3).

Fig. 2 HOMO and LUMO energies of (a) solvent-free, (b) benzene-solvated, (c) biphenyl-solvated, (d) cyclohexane-solvated C₆₀ crystals as well as (e) solvent-free, (f) 1-methylnaphthalene-solvated, (g) chlorobenzene-solvated, (h) ortho-dichlorobenzene-solvated PCBM crystals.

As shown in Fig. 3, the calculated band gap (E_g) of the pristine C₆₀ crystal is 1.36 eV, and despite being underestimated compared to the experimental value of 2.3 eV,⁷⁵ it is close to the published theoretical values of 1.09 eV obtained using PBE⁴⁴ and 1.5 eV obtained using LDA functionals⁴⁶. In addition, the E_g of 1.28 eV for monoclinic solvent-free PCBM crystal found in this work is also comparable with those of 1.21 and 1.32 eV for simple cubic and body-centered-cubic PCBM crystals, both calculated by using the LDA functional in ref. 24. Hybrid functionals such as B3LYP,^76,77 which includes a fraction of exact Hartree–Fock exchange, would provide an even more accurate prediction of band gaps, with values closer to the corresponding experimental ones. However, hybrid functionals are not currently implemented for periodic systems in ONETEP.⁷⁸ Both the pristine C₆₀ and PCBM crystals have a lower E_g compared to the isolated C₆₀ and PCBM molecules, and this means that E_g decreases for C₆₀ and PCBM from the gas phase to solid phase. The E_g for solvated C₆₀ and PCBM co-crystals are different from those of the solvent-free C₆₀ and PCBM crystals.

Fig. 3 Energy gaps of isolated C₆₀ and PCBM molecules, solvent-free C₆₀ and PCBM crystals as well as solvated fullerene co-crystals and their corresponding fullerene and solvent crystals which only retain the fullerene and the solvent parts of the co-crystal via removing the solvent and fullerene parts, respectively.

To shed light on the effect of fullerenes and solvents contained in the fullerene solvates on the E_g of solvates, we then calculated the E_g of the solvent and fullerene parts of the solvates, as shown in Fig. 3. Both the pure solvent and fullerene parts have the same lattice constants as their parent solvated co-crystals, as they are obtained by deleting the fullerene and solvent parts of the co-crystals, respectively. From Fig. 3 it is possible to observe that the E_g of pure fullerene parts derived from solvates are also lower than the E_g of the corresponding isolated fullerene molecules, except for the one derived from benzene-solvated C₆₀ whose E_g is equal to that of the isolated C₆₀ molecule. By comparing the E_g of fullerene solvates with those of the corresponding pure solvent and fullerene parts, we found that the E_g of fullerene solvates are slightly smaller than those of the corresponding pure fullerene parts, this probably being due to the effect of the solvents contained in the solvates.

The reductions of E_g of the biphenyl-solvated C₆₀ and 1-methylnaphthalene-solvated PCBM co-crystals are the largest among the C₆₀ and PCBM solvates, respectively. However, the difference in the band gap between the fullerene solvates and their corresponding fullerene parts is so small that the E_g of fullerene solvates is nearly the same as that of the corresponding fullerene parts. In contrast, the difference in E_g between solvates and their corresponding pure solvent parts is significantly larger: this implies that the E_g of solvated fullerene co-crystals is mainly determined by the fullerene part contained therein. Furthermore, the difference in E_g among solvated fullerene co-crystals should be attributed to the different E_g of fullerene parts. The reason why the fullerene parts in solvates have different E_g is that their packing structures are different. In general, the packing structure of fullerene molecules in fullerene solvates is affected by the solvent: the solvent molecules with which the fullerene molecules co-crystallise play an important role in modifying the interaction between fullerene cages via forming π-stacking interactions with fullerene molecules. Thus, the packing structure of fullerene molecules in fullerene solvates is a result of the interaction between solvent and fullerene molecules. Different solvents may result in different fullerene packing structures, which can consequently lead to variations in E_g of fullerene solvates.

Electronic density of states of fullerene crystals. Fig. 4 shows the comparison of the occupied parts of the total density of states (TDOS) for the pristine PCBM crystal and isolated PCBM molecule with the experimental and theoretical data in ref. 40. As the TDOS of the PCBM film measured by photoelectron spectroscopy (PES) in ref. 40 gave the binding energies as positive numbers, we have changed the energy of the TDOS to be positive to match the results of PES, and also shifted the energy to align the lowest energy peak of TDOS. In Fig. 4, one should note that 0 eV on the energy scale is the Fermi level, not the vacuum level since photoemission spectra are usually displayed in binding energies relative to the Fermi edge of the sample. As can be seen, the TDOSs calculated in this work agree well with the experimental photoelectron spectrum and theoretical results which were calculated by using the DMol³ package using the BLYP functional, thus proving that our results are reliable.

Fig. 4 Comparison of the occupied parts of the total density of states (TDOS) for the solvent-free PCBM crystal and isolated PCBM molecule between our results and the available experimental photoelectron spectrum and theoretical data in the literature. The TDOSs reported here were generated with a bigger Gaussian smearing width of 0.2 eV to make them more comparable with the results in the literature.

In Fig. 5, we only show the DOS of fullerene crystals near the Fermi level because it is the most important part in the DOS and governs the properties of many materials. From the total and local DOS it is clear that the LUMO states of all fullerene solvates are entirely derived from the fullerenes C₆₀ and PCBM, but the HOMO states of these solvates are not. The contributions of biphenyl and 1-methylnaphthalene to the HOMO states of their fullerene solvates are obvious, thus also indicating that the effects of biphenyl and 1-methylnaphthalene on the HOMO states of their fullerene solvates are distinctly stronger than the others. Similar effects of these two solvents on the E_g of their fullerene solvates have also been previously found, however the contributions of these two solvents are minor: the HOMO states of their fullerene solvates are still mainly determined by the fullerenes. Therefore, we can conclude that fullerenes C₆₀ and PCBM contribute greatly to the HOMO and LUMO states of their solvates, which is in agreement with the results we found above: the E_g of fullerene solvates are mainly determined by the fullerene parts contained in solvates.

Fig. 5 Total and local density of states of (a) pristine, (b) benzene-solvated, (c) biphenyl-solvated, (d) cyclohexane-solvated C₆₀ crystals as well as (e) pristine, (f) 1-methylnaphthalene-solvated, (g) chlorobenzene-solvated, (h) ortho-dichlorobenzene-solvated PCBM crystals near the Fermi level.

3.3 Absorption spectra of solvent-free and solvated fullerene crystals

The ultraviolet-visible (UV-Vis) absorption spectra of solvent-free and solvated C₆₀ and PCBM crystals are shown in Fig. 6. From Fig. 6(a) it is possible to see that the absorption peaks of the solvent-free C₆₀ crystal occurring at 449.7 nm, 394.0 nm and 334.5 nm shift about 104 nm or more towards longer wavelengths as compared to the absorption peaks of the C₆₀ thin film (346 nm, 268 nm and 220 nm) measured by Zhou et al.;³¹ this red-shift is most probably caused by the GGA-PBE functional which usually underestimates the energy gap. As the intensity of the absorption spectrum of the C₆₀ thin film is dimensionless, we are not able to make a quantitative comparison with that of the absorption spectrum of the pristine C₆₀ crystal. In Fig. 6(e), a similar red-shift of the absorption peak (at 417.3 nm) of the solvent-free PCBM crystal can also be found due to the use of the GGA-PBE functional when compared to the absorption spectrum of the PCBM thin film, as recorded by Cook et al.⁴³ However, the general trend of the absorption spectra is not affected by the systematic red-shift error introduced by the PBE functional, and therefore for pristine C₆₀ and PCBM crystals, the spectra in the present work are similar to those for C₆₀ and PCBM thin films shown in ref. 31 and 43. As shown in Fig. 6(e), a weak tail extending as far as 900 nm can be found in the absorption spectrum of both the pristine PCBM crystal and thin film. This absorption tail has also been observed in the absorption spectra of the pristine C₆₀ crystal and thin film at long wavelengths [see Fig. 6(a)], which confirms the results obtained by Cook et al.: the origin of this absorption tail for the PCBM thin film is from its parent C₆₀ cages rather than the light scattering or other experimental artifacts.⁴³ In Fig. 6(a) the weak absorption in the region from 500 to 900 nm for the pristine C₆₀ crystal and thin film is due to the high number of optical transitions that occur in this region for the C₆₀ molecule: for instance, nearly all transitions to the first excited state such as the HOMO–LUMO transition are symmetrically forbidden because of the high symmetry of the C₆₀ molecule. Only transitions to higher energy bands such as the second and third lowermost conduction bands that lie above the LUMO band, which corresponds to the absorptions of shorter waves, can be allowed.

Fig. 6 Ultraviolet-visible (UV-Vis) absorption spectra of pristine and solvated C₆₀ and PCBM crystals compared with experimental data^31,43 as well as the contributions of solvent molecules contained in fullerene solvates to the total absorption spectra.

As shown in Fig. 6, the main absorptions for all fullerene crystals studied here locate in the UV region. Comparing the absorption coefficients of C₆₀ solvates with those of the solvent-free C₆₀ crystals, we found the absorption peaks of the former are lower than the corresponding ones of the latter, with the exception of the enhanced absorption peak of benzene-solvated C₆₀ in the middle UV (MUV) region (200 to 250 nm) and the one of biphenyl-solvated C₆₀ near 200 nm. The enhanced MUV absorptions for benzene- and biphenyl-solvated C₆₀ co-crystals should be attributed to the contributions of the solvent molecules contained in these two solvates, as shown in Fig. 6(b) and (c). From Fig. 6(d) it is possible to see that there is no contribution from the cyclohexane molecules to the absorption spectrum of the corresponding C₆₀ solvate. The absorption coefficients of benzene-, biphenyl- and cyclohexane-solvated C₆₀ co-crystals in the visible region (450 to 570 nm) are much smaller than those of the solvent-free C₆₀ crystal, thus indicating that these three solvated C₆₀ co-crystals have a rather weak ability to absorb light in this region. All of the four C₆₀ crystals show no absorption when the wavelength is above 600 nm. In Fig. 6(f)–(h), similar reductions can also be found for the absorption peaks of PCBM solvates in the region from 300 to 900 nm when compared to the solvent-free PCBM crystal. However, their reductions are small relative to the ones of C₆₀ solvates, since the decrease of π-stacking interactions between PCBM molecules caused by solvents in PCBM solvates is smaller than that in C₆₀ solvates. In particular, the CB-solvated PCBM co-crystal has a very similar absorption spectrum to the solvent-free PCBM crystal, which implies that CB is the optimal solvent to keep the absorption property of the PCBM crystal at its strongest. The absorption coefficients of the three PCBM solvates in the visible region (420 to 760 nm) are almost the same. The difference in computed absorption coefficients between CB- and oDCB-solvated PCBM co-crystals is so small that it would be difficult to observe experimentally. From Fig. 6(f)–(h), it can also be seen that the absorptions of PCBM solvates around 250 nm increase because of the significant contributions of the solvents in this region.

4 Conclusions

In conclusion, linear-scaling density functional theory (LS-DFT) calculations, which allowed us to simulate large systems with thousands of atoms, have been carried out in order to clarify the crystallographic, electronic and optical properties of the electron-acceptor materials for organic photovoltaics (OPVs) C₆₀ and its derivative PCBM crystals, as well as their solvates co-crystallised with different solvents, also including the PCBM-cyclohexane solvate composed of 2208 atoms in its unit cell. Our results can be summarized as follows:

(i) The relatively low cohesive energies of pure C₆₀ and PCBM parts which originate from the C₆₀ and PCBM solvates compared with those of the solvent-free C₆₀ and PCBM crystals with closest-packing structures indicate that the fullerene parts are unstable and prone to decompose without solvents. The significant contributions of the interactions between fullerene and solvent molecules to the cohesive energies of C₆₀ and PCBM solvates also imply that solvent molecules are essential for C60 and PCBM solvates.

(ii) With an increase in the crystal density, both the HOMO and LUMO energy levels for C₆₀ and PCBM crystals move upward, in agreement with experimental observations. We have examined a variety of solvents and found that chlorobenzene and ortho-dichlorobenzene are the best from our set of solvents for the preparation of OPVs due to the relatively high LUMO energies of their PCBM solvates, which are helpful to increase the driving force for charge separation in OPV devices.

(iii) Both the energy gap (E_g) and the HOMO and LUMO states of C₆₀ and PCBM solvates are mainly determined by the fullerene parts contained in solvates. The solvents play an important role in determining the packing structures of fullerenes contained in C₆₀ and PCBM solvates, which lead to the different E_g of C₆₀ and PCBM solvates. Thus, the effect of solvents on the E_g of fullerene solvates is not due to the π-stacking interactions between solvent and fullerene molecules but merely due to the different fullerene packing structures in solvates.

(iv) In contrast to the significant changes in the intensity of the UV-vis absorption spectra of C₆₀ solvates as compared with the one of the solvent-free C₆₀ crystal, the intensity of the absorption spectra of PCBM solvates changes slightly, relative to that of the solvent-free PCBM crystal. This is because the reduction of π-stacking interactions between fullerene molecules caused by solvents in PCBM solvates is not as strong as in the case of C₆₀ solvates. The main absorption peaks for these fullerene crystals are located in the UV region. The contributions of solvent molecules to the absorption spectra of their respective fullerene solvates mainly occur in the middle UV region (200 to 300 nm) with the exception of cyclohexane whose contribution cannot be found in the entire absorption spectrum of the C₆₀ solvate.

This work was motivated by the need to examine the properties of solvated C₆₀ and PCBM co-crystals, as it is possible that they can be used instead of solvent-free fullerene crystals. We expect our obtained results to help in the design and optimisations of new OPV devices.

Acknowledgements

H. X. thanks Lanzhou University of Technology (LUT) for a scholarship from the LUT Training Program of the National and Gansu Provincial Excellent Doctoral Dissertations. G. B. thanks the EPSRC Doctoral Training Centre (grant number: EP/L015382/1) and Merck for a PhD studentship. This work was financially supported by the National Natural Science Foundation of China (11164014 and 11364025). The authors are grateful for helpful discussions with Dr Laura E. Ratcliff and Jolyon Aarons. The calculations in this work were performed on the Iridis4 Supercomputer of the University of Southampton.

References

1. G. Yu, J. Gao, J. C. Hummelen, F. Wudl and A. J. Heeger, Science, 1995, 270, 1789.
2. M. Reyes-Reyes, K. Kim and D. L. Carroll, Appl. Phys. Lett., 2005, 87, 83506.
3. C. J. Brabec, A. Cravino, D. Meissner, N. S. Sariciftci, T. Fromherz, M. T. Rispens, L. Sanchez and J. C. Hummelen, Adv. Mater., 2010, 22, 3839.
4. C. B. Nielsen, R. S. Ashraf, N. D. Treat, B. C. Schroeder, J. E. Donaghey, A. J. P. White, N. Stingelin and I. McCulloch, Adv. Mater., 2015, 27, 948.
5. H. Hu, K. Jiang, G. Yang, J. Liu, Z. Li, H. Lin, Y. Liu, J. Zhao, J. Zhang, F. Huang, Y. Qu, W. Ma and H. Yan, J. Am. Chem. Soc., 2015, 137, 14149.
6. S. Berny, N. Blouin, A. Distler, H.-J. Egelhaaf, M. Krompiec, A. Lohr, O. R. Lozman, G. E. Morse, L. Nanson, A. Pron, T. Sauermann, N. Seidler, S. Tierney, P. Tiwana, M. Wagner and H. Wilson, Adv. Sci., 2015, 3, 1500342.
7. Y. Liu, J. Zhao, Z. Li, C. Mu, W. Ma, H. Hu, K. Jiang, H. Lin, H. Ade and H. Yan, Nat. Commun., 2014, 5, 5293.
8. M. C. Scharber, Adv. Mater., 2016, 28, 1994.
9. G. Paternó, F. Cacialli and V. García-Sakai, Chem. Phys., 2013, 427, 142.
10. A. L. Ayzner, S. C. Doan, B. T. de Villers and B. J. Schwartz, J. Phys. Chem. Lett., 2012, 3, 2281.
11. F. Piersimoni, S. Chambon, K. Vandewal, R. Mens, T. Boonen, A. Gadisa, M. Izquierdo, S. Filippone, B. Ruttens, J. D'Haen, N. Martin, L. Lutsen, D. Vanderzande, P. Adriaensens and J. V. Manca, J. Phys. Chem. C, 2011, 115, 10873.
12. B. A. Collins, J. R. Tumbleston and H. Ade, J. Phys. Chem. Lett., 2011, 2, 3135.
13. R. M. Beal, A. Stavrinadis, J. H. Warner, J. M. Smith, H. E. Assender and A. A. R. Watt, Macromolecules, 2010, 43, 2343.
14. Y. Kim, J. Nelson, T. Zhang, S. Cook, J. R. Durrant, H. Kim, J. Park, M. Shin, S. Nam, M. Heeney, I. McCulloch, C.-S. Ha and D. D. C. Bradley, ACS Nano, 2009, 3, 2557.
15. A. L. Ayzner, C. J. Tassone, S. H. Tolbert and B. J. Schwartz, J. Phys. Chem. C, 2009, 113, 20050.
16. T. A. Bull, L. S. C. Pingree, S. A. Jenekhe, D. S. Ginger and C. K. Luscombe, ACS Nano, 2009, 3, 627.
17. L. Li, G. Lu, S. Li, H. Tang and X. Yang, J. Phys. Chem. B, 2008, 112, 15651.
18. P. A. C. Quist, T. J. Savenije, J. M. Schins, J. E. Kroeze, P. A. Rijkers and L. D. A. Siebbeles, Phys. Rev. B: Condens. Matter Mater. Phys., 2007, 75, 195317.
19. A. Swinnen, I. Haeldermans, M. vande Ven, J. D'Haen, G. Vanhoyland, S. Aresu, M. D'Olieslaeger and J. Manca, Adv. Funct. Mater., 2006, 16, 760.
20. X. Yang, J. Loos, S. C. Veenstra, W. J. H. Verhees, M. M. Wienk, J. M. Kroon, M. A. J. Michels and R. A. J. Janssen, Nano Lett., 2005, 5, 579.
21. T. J. Savenije, J. E. Kroeze, X. Yang and J. Loos, Adv. Funct. Mater., 2005, 15, 1260.
22. X. Yang, J. K. J. van Duren, M. T. Rispens, J. C. Hummelen, R. A. J. Janssen, M. A. J. Michels and J. Loos, Adv. Mater., 2004, 16, 802.
23. J. K. J. van Duren, X. Yang, J. Loos, C. W. T. Bulle-Lieuwma, A. B. Sieval, J. C. Hummelen and R. A. J. Janssen, Adv. Funct. Mater., 2004, 14, 425.
24. J. M. Nápoles-Duarte, M. Reyes-Reyes, J. L. Ricardo-Chavez, R. Garibay-Alonso and R. Lopez-Sandoval, Phys. Rev. B: Condens. Matter Mater. Phys., 2008, 78, 35425.
25. C. Menelaou, S. Tierney, N. Blouin, W. Mitchell, P. Tiwana, I. McKerracher, C. Jagadish, M. Carrasco and L. M. Herz, J. Phys. Chem. C, 2014, 118, 17351.
26. L. Zheng and Y. Han, J. Phys. Chem. B, 2012, 116, 1598.
27. M. T. Rispens, A. Meetsma, R. Rittberger, C. J. Brabec, N. S. Sariciftci and J. C. Hummelen, Chem. Commun., 2003, 16.
28. C. D. Reddy, Z. G. Yu and Y.-W. Zhang, Sci. Rep., 2015, 5, 12221.
29. J. J. Richards, A. H. Rice, R. D. Nelson, F. S. Kim, S. A. Jenekhe, C. K. Luscombe and D. C. Pozzo, Adv. Funct. Mater., 2013, 23, 514.
30. A. Pénicaud, O. Y. Carreón, A. Perrier, D. J. Watkin and C. Coulon, J. Mater. Chem., 2002, 12, 913.
31. W. Y. Zhou, S. S. Xie, S. F. Qian, T. Zhou, R. A. Zhao, G. Wang, L. X. Qian and W. Z. Li, Sci. China, Ser. A: Math., Phys., Astron., 1996, 39, 864.
32. M. Jansen and G. Waidmann, Z. Anorg. Allg. Chem., 1995, 621, 14.
33. M. F. Meidine, P. B. Hitchcock, H. W. Kroto, R. Taylor and D. R. M. Walton, J. Chem. Soc., Chem. Commun., 1992, 1534.
34. H.-B. Bürgi, E. Blanc, D. Schwarzenbach, S. Liu, Y.-J. Lu, M. M. Kappes and J. A. Ibers, Angew. Chem., Int. Ed. Engl., 1992, 31, 640.
35. S. Liu, Y.-J. Lu, M. M. Kappes and J. A. Ibers, Science, 1991, 254, 408.
36. G. M. Lazzerini, G. M. Paternò, G. Tregnago, N. Treat, N. Stingelin, A. Yacoot and F. Cacialli, Appl. Phys. Lett., 2016, 108, 53303.
37. G. Paternò, A. J. Warren, J. Spencer, G. Evans, V. G. Sakai, J. Blumberger and F. Cacialli, J. Mater. Chem. C, 2013, 1, 5619.
38. M. Casalegno, S. Zanardi, F. Frigerio, R. Po, C. Carbonera, G. Marra, T. Nicolini, G. Raos and S. V. Meille, Chem. Commun., 2013, 49, 4525.
39. Y. Yang, C. Liu, S. Gao, Y. Li, X. Wang, Y. Wang, T. Minari, Y. Xu, P. Wang, Y. Zhao, K. Tsukagoshi and Y. Shi, Mater. Chem. Phys., 2014, 145, 327.
40. C.-Q. Sheng, P. Wang, Y. Shen, Y.-J. Li, W.-H. Zhang, F.-Q. Xu, J.-F. Zhu, G.-Q. Lai and H.-N. Li, Chin. Phys. B, 2012, 21, 17102.
41. Y. Zhong, S. Izawa, K. Hashimoto, K. Tajima, T. Koganezawa and H. Yoshida, J. Phys. Chem. C, 2014, 119, 23.
42. S. Falke, P. Eravuchira, A. Materny and C. Lienau, J. Raman Spectrosc., 2011, 42, 1897.
43. S. Cook, H. Ohkita, Y. Kim, J. J. Benson-Smith, D. D. C. Bradley and J. R. Durrant, Chem. Phys. Lett., 2007, 445, 276.
44. S. Jalali-Asadabadi, E. Ghasemikhah, T. Ouahrani, B. Nourozi, M. Bayat-Bayatani, S. Javanbakht, H. A. R. Aliabad, I. Ahmad, J. Nematollahi and M. Yazdani-Kachoei, J. Electron. Mater., 2016, 45, 339.
45. N. Troullier and J. L. Martins, Phys. Rev. B: Condens. Matter Mater. Phys., 1992, 46, 1754.
46. S. Saito and A. Oshiyama, Phys. Rev. Lett., 1991, 66, 2637.
47. F. Gajdos, H. Oberhofer, M. Dupuis and J. Blumberger, J. Phys. Chem. Lett., 2013, 4, 1012.
48. F. Frigerio, M. Casalegno, C. Carbonera, T. Nicolini, S. V. Meille and G. Raos, J. Mater. Chem., 2012, 22, 5434.
49. D. L. Cheung and A. Troisi, J. Phys. Chem. C, 2010, 114, 20479.
50. C. J. Brabec, S. Gowrisanker, J. J. M. Halls, D. Laird, S. Jia and S. P. Williams, Adv. Funct. Mater., 2001, 11, 374.
51. M. C. Scharber, D. Mühlbacher, M. Koppe, P. Denk, C. Waldauf, A. J. Heeger and C. J. Brabec, Adv. Mater., 2006, 18, 789.
52. F. Zhang, K. G. Jespersen, C. Björström, M. Svensson, M. R. Andersson, V. Sundström, K. Magnusson, E. Moons, A. Yartsev and O. Inganäs, Adv. Funct. Mater., 2006, 16, 667.
53. C.-K. Skylaris, P. D. Haynes, A. A. Mostofi and M. C. Payne, J. Chem. Phys., 2005, 122, 84119.
54. H. P. G. Thompson and G. M. Day, Chem. Sci., 2014, 5, 3173.
55. M. A. Neumann, F. J. J. Leusen and J. Kendrick, Angew. Chem., Int. Ed., 2008, 47, 2427.
56. P. G. Karamertzanis, A. V. Kazantsev, N. Issa, G. W. A. Welch, C. S. Adjiman, C. C. Pantelides and S. L. Price, J. Chem. Theory Comput., 2009, 5, 1432.
57. M. Vasileiadis, A. V. Kazantsev, P. G. Karamertzanis, C. S. Adjiman and C. C. Pantelides, Acta Crystallogr., Sect. B: Struct. Sci., 2012, 68, 677.
58. E. J. Baerends, O. V. Gritsenko and R. van Meer, Phys. Chem. Chem. Phys., 2013, 15, 16408.
59. O. V. Gritsenko and E. J. Baerends, Can. J. Chem., 2009, 87, 1383.
60. D. B. Cook, Int. J. Quantum Chem., 1996, 60, 793.
61. C.-K. Skylaris, A. A. Mostofi, P. D. Haynes, O. Diéguez and M. C. Payne, Phys. Rev. B: Condens. Matter Mater. Phys., 2002, 66, 35119.
62. A. A. Mostofi, P. D. Haynes, C.-K. Skylaris and M. C. Payne, J. Chem. Phys., 2003, 119, 8842.
63. J.-L. Fattebert and J. Bernholc, Phys. Rev. B: Condens. Matter Mater. Phys., 2000, 62, 1713.
64. S. J. Clark, M. D. Segall, C. J. Pickard, P. J. Hasnip, M. I. J. Probert, K. Refson and M. C. Payne, Z. Kristallogr. – Cryst. Mater., 2005, 220, 567.
65. L. E. Ratcliff, N. D. M. Hine and P. D. Haynes, Phys. Rev. B: Condens. Matter Mater. Phys., 2011, 84, 165131.
66. L. E. Ratcliff and P. D. Haynes, Phys. Chem. Chem. Phys., 2013, 15, 13024.
67. L. E. Ratcliff, Optical absorption spectra calculated using linear-scaling density-functional theory, Springer, 2013.
68. S. Grimme, J. Comput. Chem., 2006, 27, 1787.
69. J. F. Dobson, A. White and A. Rubio, Phys. Rev. Lett., 2006, 96, 073201.
70. J. B. A. Davis, F. Baletto and R. L. Johnston, J. Phys. Chem. A, 2015, 119, 9703.
71. S. F. Alshahateet, R. Bishop, D. C. Craig and M. L. Scudder, CrystEngComm, 2003, 5, 417.
72. G. E. Morse, I. Gong, Y. Kawar, A. J. Lough and T. P. Bender, Cryst. Growth Des., 2014, 14, 2138.
73. I. S. Neretin and Y. L. Slovokhotova, Russ. Chem. Rev., 2004, 73, 455.
74. C. Pan, M. P. Sampson, Y. Chai, R. H. Hauge and J. L. Margrave, J. Phys. Chem., 1991, 95, 2944.
75. R. W. Lof, M. A. van Veenendaal, B. Koopmans, H. T. Jonkman and G. A. Sawatzky, Phys. Rev. Lett., 1992, 68, 3924.
76. A. D. Becke, J. Chem. Phys., 1993, 98, 5648.
77. C. Lee, W. Yang and R. G. Parr, Phys. Rev. B: Condens. Matter Mater. Phys., 1988, 37, 785.
78. J. Dziedzic, Q. Hill and C.-K. Skylaris, J. Chem. Phys., 2013, 139, 214103.

---

Footnote

† Electronic supplementary information (ESI) available. CCDC 1519619. For ESI and crystallographic data in CIF or other electronic format see DOI: 10.1039/c6cp08165g

---