1. Introduction

The use of organic materials as the active elements in a diverse range of plastic electronic products has become commercially attractive due to the possibility of synthetically optimising the constituent materials to the target application. Novel applications have been demonstrated in niche markets such as displays, imaging systems and chemical sensors. Many plastic electronic products demand that the organic materials are deposited as a thin-film over a large area using either solvent based casting or vacuum deposition methods. The viability of candidate organic materials for device applications is normally assessed using such thin-films and a benchmark parameter that is often experimentally measured is the carrier mobility using transient photoconductivity or current injection methods[1]. Such experiments have generated a valuable database of mobility parameters for a diverse range of organic materials. An interesting subset of this data is focussed upon organic materials whose transport occurs via hopping processes which are controlled by energetic disorder due to background dipole electric fields[2-4]. Such hopping systems frequently involve redox reactions between pole moment. Examples of such systems include pure dipolar glasses (DGs) where the molecules are usually condensed into a thin-film using thermal evaporation, and molecularly doped polymers (MDPs) where the small molecules are used to dope a supportive polymer using spin-coating from a common solvent. It is normally desirable to deposit such DGs and MDPs in an amorphous form to ensure that the electrical and mechanical properties of the film are globally maintained. The amorphous structure ensures that the orientation of internal dipole directions is random and that the associated Gaussian distribution of energetic sites through which carriers hop is correlated[5,6]. Monte Carlo simulations of hopping transport through such correlated Gaussian landscapes has resulted in the Gaussian disorder model (GDM)[7]. The GDM has been successfully used to characterise the dependence of the mobility to electric field and temperature for many DG and MDP materials[7,8].

The assumption that the dipole distribution remains entirely random throughout DG and MDP systems may be invalidated, however, under circumstances where local ordering of neighbouring dipole moments is favoured. Local ordering is likely to arise as a result of dipole-dipole interactions that result in molecular agglomeration and crystallisation[9-12]. Under these circumstances the ability of the GDM model to parameterise the associated hopping transport may be compromised. A particularly important test for the GDM model under these circumstances is to examine to what extent it is capable of detecting the underlying presence of local dipole ordering through the parameters that are returned. To investigate the GDM sensitivity to local dipole order a series of Monte Carlo hopping simulations have therefore been conducted using initially amorphous landscapes into which a number of locally ordered regions were introduced. Simulations with local ordering have been performed to model both DGs and MDPs. This work extends independent previous Monte Carlo studies which have either considered short-range order in MDPs using non-correlated random orientation vectors[12], or have examined the effect of introducing ordered regions into DG energy landscapes by using pre-determined Gaussian distributions to perform the allocation of site energies[13-15]. By contrast the generation of energy landscapes in the present study has been achieved using a fundamental electrostatic calculation of the dipole potential, where local variations of the potential have been introduced using dipoles that are locally ordered within cubic regions that have a specified volume. Using this approach, the sensitivity of hopping transport to the total concentration of ordered dipoles, and the associated ordered volume size that is employed, has been independently evaluated.

2. Simulation Details

Monte Carlo simulation of carrier hopping through DG and MDP structures was performed using cubic lattices where the fundamental cube size comprised NxNxN sites with a lattice constant of size a. For DGs dipoles having an identical dipole moment (p) were randomly oriented at each of the fundamental cube lattice sites. The resulting electronic energy for each site was then determined using an Ewald summation of the electric potentials due to all other surrounding dipoles[16]. To implement periodic boundary conditions and ensure that site energies near the faces of the fundamental cube would remain correlated the Ewald calculation included the electric potential contribution from the surrounding 27 identical copies of the fundamental cube. The energy scale of landscapes generated using this procedure was proportional to p/4πεa² where ε is the DG permittivity.

A subset of smaller nxnxn ordered cubes (OCs) were introduced into these DG landscapes to evaluate the influence of local dipole ordering. These OC regions were randomly seeded throughout the fundamental cube and for each OC location the dipoles were anti-aligned with their nearest neighbours to effect cancellation of dipole potentials and so introduce lower local site energies relative to the surrounding DG landscape. The anti-alignment orientation assigned to each OC was, however, randomly selected. The individual volume of the OCs introduced was restricted to (n/N)³ ≤ 0.125 and the number (m) of such cubes was controlled to give an overall concentration (C) of ordered to total fundamental sites. Due to possible overlap of the OC site locations due to the random seeding procedure it was generally found that C ≤ m(n/N)³ with equality only holding for dilute OC lattices. For example, to achieve C = 0.5 it was found that m ~ 80 for n/N = 0.2 whereas m ~ 5 for n/N = 0.5.

Monte Carlo hopping through the generated energy landscapes was conducted using a standard Miller-Abrahams selection criteria[17] for nearest and next-nearest neighbour sites which ignores the dependence of the hopping rate upon the mutual orientation of the transport sites. To ensure that energetic equilibration of carriers was established during hopping a macroscopic transport distance of 2000 lattice planes along the applied electric field direction was used to determine individual carrier transit times. Statistical averaging of transit times using a minimum population of 20 carriers was used to calculate the hopping mobility (µ) associated with the applied electric field strength (F). The electric field strengths used during simulation generally spanned a broad range with F < ~ 1.5 F_max where F_max represents the upper field limit where the energy landscape is “flattened” and carrier velocities saturate[18]. To fully parameterise the mobility results with the GDM, simulations were furthermore conducted across a broad range of temperatures for a selected energy landscape. Appropriate temperature ranges were dictated by the energetic characteristics of the landscape and the achievement of practical simulation times.

3. Results and Discussion

3.1. Simulation of Dipolar Glasses

Previous work has demonstrated that the distribution of site energies in pure DG lattices is well described by a single Gaussian density-of-states (DOS)[5,19]. The standard deviation (σ_DOS) of the DOS is consequently a convenient parameter to characterise the underlying energetic disorder of the associated lattice and may be determined through application of the GDM to mobility data[3,4,7]. In the present work the continued ability of a Gaussian distribution to describe the DOS of DG lattices into which OC regions were introduced was first examined. To generate a sufficient number of site energies whose statistical distribution could be realistically examined through binning the minimum fundamental cube size used for this study corresponded to N = 10. The use of larger N sizes generated increasingly smoother distributions for a fixed bin interval (typically 30 bins across the entire energy spectrum) but the additional computational overhead encountered limited the maximum fundamental cube size to N = 15. Fitting of the lattice DOS to a Gaussian distribution was performed using a least-squares procedure to generate σ_DOS parameters and this process was successively repeated over five sets of random number seeds to allow σ_DOS averaging. The results of these DOS studies are presented in Figure 1 which shows the dependence of σ_DOS upon the overall concentration of inserted OCs for various n/N ratios.

In Figure 1 the data have been normalised relative to the reference σ_DOS for the pure (m = 0) DG lattice. As anticipated the amount of energetic disorder is observed to progressively fall as the total concentration of OC sites becomes larger. The rate of fall of σ_DOS depends upon the n/N ratio where it is noted that the insertion of large OCs is more effective than small OCs in reducing the overall energetic disorder at high C values. The relative inefficiency of smaller n/N OCs to reduce σ_DOS is understandable, however, since with increasing m the probability of overlap of OC locations is enhanced and dipole anti-alignment at the overlap locations is lost. This phenomenon is further illustrated in Figure 1 for the n/N = 0.5 (placed) and single cube data sets. In the n/N = 0.5 (placed) situation the seed locations for the OCs have been specified to ensure that no overlap occurs and in this case C = 1 is achieved when m = 8 exactly. Compared to the randomly seeded n/N = 0.5 data the n/N = 0.5 (placed) results show an enhanced σ_DOS fall rate with C. The rate of fall for the placed results is comparable with the single cube data which has been obtained for DG lattices into which a single n/N cube has been inserted. Completely anti-aligned ordering of all dipoles in the fundamental cube may be achieved under these circumstances for n/N = 1 where it is confirmed in Figure 1 that σ_DOS → 0.

Figure 1. Dependence of the normalised DOS energetic disorder of DGs upon the total concentration of inserted OCs for various n/N cube ratios. The dashed lines represent the best fit of the data to Equation 2

To assess the effect of OC insertion upon hopping transport Monte Carlo simulations were performed using selected energy landscapes from Figure 1. Examples of simulated mobility data are given in Figure 2 for a DG where n/N = 0.5, m = 3 and C = 0.36.

The electric field dependence of the mobility plotted in Figure 2(a) continues to exhibit the characteristic µ ∝ √F dependence that is predicted by the GDM for F < F_max[6,7,10] and extrapolation of the data in this field region to F = 0 allows zero-field mobility values (µ₀) to be determined at a given temperature (T). Analysis of µ₀(T) within the GDM framework enables the underlying transport energy disorder (σ_GDM) to be established since[3,7]:

(1)

In Equation 1 k_B = Boltzmann’s constant and µ₀₀ is the lattice mobility in the absence of disorder when T → ∞_. From Figure 2(b) it would appear that the application of Equation 1 to the µ₀(T) data remains appropriate across the temperature range that was typically employed for the simulations (σ_GDM / 5k_B < T < 2σ_GDM / k_B) allowing σ_GDM to be extracted. For pure DGs it is expected that σ_GDM returned by the GDM analysis should be equal to σ_DOS found using the DOS analysis.

Figure 2. Examples of simulated Monte Carlo mobility data for a DG where n/N = 0.5, m = 3 and C = 0.36. The energy scale used corresponded to p/4πεa2 = 43 meV (a) electric field dependence of the mobility at temperatures indicated (b) temperature dependence of µ0 where the solid line represents the best fit of the data to Equation 1

The general behaviour illustrated by the Figure 2 plots was similarly found when Monte Carlo simulations were subsequently conducted using other landscapes across the entire spectrum of n/N and C parameters. As expected the µ₀₀ parameter was found in all cases to be constant for the chosen fixed lattice size a. A summary of the energy disorder found using the GDM analysis of the Monte Carlo results is given in Figure 3(a) which plots σ_GDM against σ_DOS.

For convenience the data in Figure 3 have been normalised to the reference σ_DOS and σ_GDM values for the pure DG. It is immediately evident from this plot that despite the apparent DOS control suggested from Figure 1 the accompanying energy disorder experienced by carriers as they hop through the landscapes is ultimately less dramatic. There appears to be little if any direct correlation between σ_GDM and σ_DOS although closer inspection of the Figure 3(a) data appears to suggest that there is a systematic shift of the data relative to the σ_GDM = σ_DOS reference line as n/N increases. This may indicate that as the effectiveness of local ordering is improved through the introduction of larger and more numerous OCs the attempt to model the entire lattice DOS using a single Gaussian distribution that is characterised by σ_DOS is invalidated.

Figure 3. Comparison of normalised GDM and DOS energy disorder magnitudes for DGs with inserted OC regions as indicated. (a) σGDM against σDOS; dashed line represents σGDM = σDOS as a reference guide (b) σGDM against σDOS(NOC) ; dashed line represents σGDM = σDOS(NOC) as a reference correlation guide

The DOS fitting procedure was consequently re-examined and this revealed that for some of the energy landscapes generated using high n/N and high C parameters the values for σ_DOS that were being returned were influenced by a strong energetically narrow central peak relative to a much weaker distribution of energetically broader wing states. Similar observations have been noted in independent work[5] for pure DG landscapes where dipole orientations were restricted. It is tempting to associate these peak states and wing states with the OC sites and the non-OC (NOC) sites respectively. To investigate this possibility the Gaussian fitting procedure was therefore successively applied to site energies that were strictly associated with first the OC and then the NOC locations. No systematic difference between the mean energies of the OC and NOC sites was found to exist and the DOS results from these separate analyses, again normalised to the pure DG reference, are presented in Figure 4 for various n/N ratios. Whereas there is comparatively little difference in the σ_DOS values for OC and NOC sites when n/N =2 (Figure 4(a)) a significant difference is revealed for higher n/N values (Figures 4(b,c)). Across a wide range of C σ_DOS for the OC sites (σ_DOS(OC)) appears to be insensitive to the unordered dipole fields from surrounding NOC sites which would be expected given the short-range nature of the dipole field. By contrast σ_DOS for the NOC sites (σ_DOS(NOC)) starts to converge towards the lower σ_DOS(OC) magnitude as the disordered dipole background is progressively diluted and the lattice contains only locally aligned dipoles. It is noted that in a simplified situation where dipole field interference between the OC and NOC is eliminated, both σ_DOS(OC) and σ_DOS(NOC) will assume constant values that are independent of C. Under these circumstances the standard deviation for the energy of a randomly selected site in the fundamental lattice will be given as;

(2)

Figure 4. Normalised DOS energy disorder parameters of DGs upon the total concentration of inserted OCs. The open symbols correspond to σDOS(NOC) data, and the closed symbols to σDOS(OC) data. The dashed lines correspond to the fitted σDOS results in Figure 1 (a) n/N = 0.2 (b) n/N = 0.3 (c) n/N = 0.5

Equation 2 predicts that the reduction of σ_DOS between the (higher) OC and (lower) NOC magnitudes should be linear with C. The data given in Figure 1 exhibit such linearity and the σ_DOS(NOC) values returned through the application of Equation 2 are found to be consistent with those found through the independent NOC site analyses performed in Figure 4.

Whilst the partitioning of the lattice sites is conceptually simple to perform from a DOS perspective it is unclear how the separate σ_DOS parameters for the OC and NOC regions influence the ensuing σ_GDM magnitude that is returned by the Monte Carlo simulations. However, as previously noted from Figure 3(a) the response of σ_GDM is not as sensitive as σ_DOS to the introduction of local order where it is typically found that to achieve a 50% reduction in σ_GDM relative to the DG reference magnitude it is necessary for C to exceed about 0.9 for n/N > 0.3. Similar behaviour has been reported in independent Monte Carlo studies[13-15] which was attributed to weakly percolating transport where the majority of carriers follow low-mobility paths through NOC sites. The preference to hop via the NOC sites is clearly understood with reference to Figure 4 where it is noted that after the injected carriers have had an opportunity to thermalise they will tend to occupy the lowest energy states associated with the σ_DOS(NOC) tail. The energy mismatch with the higher σ_DOS(OC) states would then be expected to constrain hopping within the NOC associated sites and this situation will be imposed until σ_DOS(OC) and σ_DOS(NOC) converge which for a high n/N demands a high C (> 0.9). The disorder parameter σ_GDM returned by the GDM analysis of mobility data is thus anticipated to correspond to σ_DOS(NOC) and to investigate this the σ_GDM data have been replotted in Figure 3(b) against the σ_DOS(NOC) data from Figure 4. Compared to the original plot against σ_DOS in Figure 3(a) the correlation of these disorder parameters is found to be significantly improved. The important implication of these results is that GDM analysis of experimental mobility data will be completely insensitive to the presence of local order within DGs across a broad range of lower C conditions where σ_GDM ≈ σ_DOS(NOC). Direct information concerning the amount of disorder that is associated with the OC regions themselves may only be accurately determined under extremely high C conditions when σ_GDM → σ_DOS(OC).

3.2. Experimental Dipolar Glass Data

An interesting DG material that is prone to crystallisation when deposited as a thin-film is N,N′-Bis (3-methylphenyl) -N,N′-diphenylbenzidine (TPD). A considerable literature exists concerning the determination of the mobility of hole carriers in TPD using direct time-of-flight methods[20-25], and indirect techniques such as current injection[25,26], photo-conductivity[27,28] and field-effect transistor configurations[25,29,30]. A survey of the TPD hole mobility values that have been obtained as F → 0 using these approaches reveals that the highest values (~ 10^-3 cm²V^-1s^-1) are consistently returned from the direct time-of-flight method. Other indirect methods may sometimes return mobility values that surprisingly may range up to almost two orders of magnitude lower[27,30]. Very little information regarding the associated structural condition of these TPD films is available but it seems unlikely that the highest mobility values reflect significant dipole ordering as the direct time-of-flight measurements are generally performed using thicker films which are less prone to crystallisation. A more realistic explanation is that the indirect methods demand the formation of ohmic contacts and these may be difficult to achieve for thinner films where crystallisation is more easily initiated through surface humidity exposure[25,31,32]. Confirmation of the relative insensitivity of the fundamental TPD hole mobility to crystallisation in thin-film structures is suggested from novel surface time-of-flight results[24] which give low-field mobility values that compare well with the conventional time-of-flight method performed through the film thickness. The stability of the TPD mobility to crystallisation would indicate that a significant fraction of dipoles sites remain unordered according to the present Monte Carlo characteristics and that µ₀ ~ 10^-3 cm²V^-1s^-1 is determined by the underlying σ_DOS(NOC) magnitude of these disordered TPD dipoles.

3.3. Simulation of Molecularly Doped Polymers

The controlling influence of the NOC sites upon the transport mechanism in DGs may be further investigated by removing the entire NOC dipole population. The removal procedure leaves a hopping lattice that is populated only by the OC sites. Such a pure OC lattice could be used to represent a MDP system (MDP(OC)) in which all of the dopant dipoles have agglomerated and anti-aligned within an inert polymer matrix that does not facilitate transport. Alternative MDP systems would include a lattice in which the dipoles have again fully agglomerated into cubes but remain disordered (MDP(DC)), and a reference MDP (MDP(R)) where no agglomeration occurs and the dipoles are randomly dispersed throughout the polymer matrix. For each of the MDP systems described the effects of agglomeration and dipole alignment may be explored as a function of the fraction (C) of total sites in the fundamental cube that are occupied by dopant dipoles. Earlier Monte Carlo work using MDP(R) lattices has been conducted to investigate the dependence of F_max upon C[18].

The dependence of σ_DOS upon C for each of the MDP scenarios is summarised in Figure 5 where data for the MDP(OC) and MDP(DC) lattices have been generated using several n/N cube ratios. All of the σ_DOS results given in Figure 5 employed identical magnitudes for p and a to ensure a consistent energy scale, and the Gaussian fitting procedure was only applied to the distribution of lattice energies at occupied dipole sites as these are the relevant sites through which hopping may proceed. A single Gaussian provided an excellent fit in all MDP cases and the data have been normalised to the σ_DOS magnitude for the MDP(R) system when C = 1.

For the reference MDP(R) system σ_DOS shows the expected increase with C (σ_DOS ∝ C^2/3[5, 18]) as the nearest-neighbour dipole distance becomes progressively shorter. A similar dependence is observed for the MDP(DC) data at higher concentrations where the effect of agglomeration upon nearest-neighbour distance becomes insignificant compared to the MDP(R) case. However, for lower concentrations (C < 0.5) agglomeration maintains a much shorter nearest-neighbour distance compared to the reference situation which results in higher agglomerated site energies. The independence of the DOS to the agglomerated cube size in the MDP(DC) system reflects the short-range nature of the random dipole potential. By contrast the combined effects of agglomeration and dipole alignment is seen to have a dramatic effect in Figure 5 upon the DOS of MDP(OC) systems. Dipole anti-alignment is highly efficient at reducing the effective dipole potential within individual cube sites so that the overall energy disorder is consequently smaller for MDP(OC) systems relative to both the MDP(DC) and MDP(R) systems across the entire C range. The effectiveness of disorder reduction is generally improved for larger n/N cube sizes as these offer less opportunity to upset the anti-alignment of dipoles in regions of OC overlap.

Figure 5. Dependence of the normalised DOS energetic disorder of MDPs upon the total concentration of dopant dipoles for various n/N cube ratios.. The open symbols correspond to MDP(DC) systems, the closed symbols to MDP(OC) systems, and the cross symbols to MDP(R) systems.The dashed lines represent a guide for the eye in all cases

Figure 6. Comparison of normalised GDM and DOS energy disorder magnitudes for MDPs using various n/N cube ratios. The open symbols correspond to MDP(DC) systems, the closed symbols to MDP(OC) systems, and the cross symbols to MDP(R) systems. The dashed line represents σGDM = σDOS as a reference correlation guide

Monte Carlo simulations of hopping transport using the above MDP landscapes were again analysed within the GDM framework. The results of these simulations are summarised in Figure 6 where normalised values of σ_GDM are plotted against the corresponding normalised σ_DOS values from Figure 5.

It is evident that the correlation between these parameters for MDPs is much stronger compared to the equivalent plot (Figure 3(a)) for the DGs. The improved correlation in the MDP systems is attributed to the absence of a parallel transport path which constrains carriers to hop via the DOS associated with occupied dipole sites and thus ensures that σ_GDM = σ_DOS. In principle experimental mobility data should therefore be capable of detecting the existence of agglomeration in MDP systems by comparing measured σ_GDM values with the MDP(R) reference magnitude. It is noted from Figure 5, however, that such a comparison is likely to be most effective for MDP(DC) systems when C < 0.5, and for low n/N MDP(OC) systems when C > 0.5, which may not correspond to actual doping ranges in which agglomeration is initiated for a selected combination of dopant molecule and polymer. The results that are presented in Figure 5 and Figure 6 highlight the extreme effect of agglomeration in MDPs whereby all of the dipoles in the MDP(R) system are forced to fully agglomerate into either the MDP(DC) or MDP(OC) structures. Other possibilities may of course be envisaged that involve only partial agglomeration of the MDP(R) dipole population, or MDP systems in which there exists a mix of DC and OC agglomerations. For all of these alternative MDP structures the possible mix of disordered and ordered dipole regions is expected to introduce parallel transport paths as found in the preceding study of DGs. The associated correlation between σ_GDM and σ_DOS will be weaker and the detection of agglomeration using mobility measurements will consequently be less sensitive.

3.4. Experimental Molecularly Doped Polymer Data

A further experimental complication that may exist in some MDP materials originates from the choice of host polymer which may introduce an additional background fluctuation (σ_P) to that arising from the dopant dipoles (σ_DOS). The existence of a finite σ_P will depend upon whether the polymer itself possesses an inherent dipole moment or may be induced to become dipolar through internal van der Waals interactions. Under the normal assumption that σ_P and σ_DOS are uncorrelated the effective MDP energy disorder detected by GDM is then σ_GDM = √( σ_P² + σ_DOS²)[3,5]. Provided σ_P is not dependent upon the internal structure of the dopant dipoles, the relative locations of the DOS characteristics given in Figure 5 (where σ_P = 0) remain invariant, however, and the detection of agglomeration is still possible. Polystyrene based MDPs offer a promising opportunity to study agglomeration effects as they provide a relatively low σ_P background and appear to induce clustering of small dopant molecules such as TPD[11].

5. Conclusions

Hopping in DGs continues to proceed mainly via a transport path that is associated with disordered dipole sites until a significant fraction of these sites are replaced by ordered dipole regions. A significant concentration of ordered dipole regions may therefore remain undetectable in DGs via GDM analysis of mobility data. The effect of dipole ordering is much more apparent in MDPs where the dopants are forced to agglomerate. Agglomeration effectively removes the transport path that is associated with disordered dipole sites in DGs and constrains carriers to hop exclusively via the agglomerated sites. Mobility measurements are consequently anticipated to provide a more sensitive detection method for the presence of ordered dipole regions in fully agglomerated MDPs compared to DGs.