Field emission: calculations supporting a new methodology of comparing theory with experiment

This paper provides a demonstration-of-concept of a new methodology for comparing field electron emission (FE) theory and experiment. It uses the parameter κ in the mathematical equation Im = CVmκ exp[–B/Vm] (where B and C are weakly varying or constants) that is taken to describe how measured current Im depends on measured voltage Vm for electronically ideal FE systems (i.e. systems that (i) have constant configuration during voltage application and (ii) have Im(Vm) given by the emission physics alone). Experimental parameter values (κm) are used to compare two alternative FE theories, for which allowable (but different) κ ranges have been established. At present, contributions to the ‘total theoretical κ’ made by voltage dependence of notional emission area are not well known: simulations reported here provide data about four commonly investigated emitter shapes. The methodology is then applied to compare 1928/1929 Fowler–Nordheim (FN) FE theory and 1956 Murphy–Good (MG) FE theory. It is theoretically certain that the 1956 theory is ‘better physics’ than the 1928/1929 theory. As in previous attempts to reach known correct theoretical conclusions by experimentally based argument, the new methodology tends to favour MG FE theory, but is formally indecisive at this stage. Further progress needs better methods of establishing error limits and of measuring κm.

varying or constants) that is taken to describe how measured current I m depends on measured voltage V m for electronically ideal FE systems (i.e. systems that (i) have constant configuration during voltage application and (ii) have I m (V m ) given by the emission physics alone). Experimental parameter values (κ m ) are used to compare two alternative FE theories, for which allowable (but different) κ ranges have been established. At present, contributions to the 'total theoretical κ' made by voltage dependence of notional emission area are not well known: simulations reported here provide data about four commonly investigated emitter shapes. The methodology is then applied to compare 1928/1929 Fowler-Nordheim (FN) FE theory and 1956 Murphy-Good (MG) FE theory. It is theoretically certain that the 1956 theory is 'better physics' than the 1928/1929 theory. As in previous attempts to reach known correct theoretical conclusions by experimentally based argument, the new methodology tends to favour MG FE theory, but is formally indecisive at this stage. Further progress needs better methods of establishing error limits and of measuring κ m .
As part of this wider project, it is wished to (i) improve existing methodologies for analysing FE currentvoltage data and (ii) improve existing procedures for using analysis results to derive useful scientific conclusions.
Sections 2-4 set out some detailed background to the present work; the remaining sections then present technical details and conclusions.
2. General background 2.1. Issues related to overall field emission system behaviour The term FE system is defined to include all aspects of the experimental system that can affect the relationship between measured current I m and measured voltage V m , including: emitter 'configuration' (composition, geometry and surface condition); the mechanical, geometrical and electrical arrangements in the vacuum system; all aspects of the electronic circuitry and of all electronic measurement instruments; the emission physics; and ALL relevant physical processes that might be taking place (for example, the generation of field-emitted vacuum space charge, Maxwell stressinduced reversible changes in emitter geometry and adsorbate atom dynamics).
This I m (V m ) relationship can be discussed as follows. It is assumed that there exists a definite (albeit unknown) relationship I m (F C ) between I m and the magnitude F C of the electrostatic field at some defined characteristic location 'C' on the emitter surface, and that the relationship between F C and V m can be written where ζ C is a parameter called the (characteristic) local (measured) voltage conversion length (LVCL). The LVCL ζ C (discussed further below) is a system characterization parameter. The description 'local' indicates that it applies to some specific location on the emitter surface and involves the local field magnitude at that location. An LVCL has the dimensions of length, but does not in practice correspond to any geometrical length in the FE system. An FE system is termed electronically ideal if (i) there are no changes in system configuration as a result of application of voltage and (ii) there are no voltage-or current-dependent changes in the value of any LVCL as a result of 'system complications'. The most commonly discussed system complications are series resistance, field-emitted vacuum space charge and system geometry change due to Maxwell stresses, but there are many others. For a list of the 12 types of system complication currently known, see https://doi.org/10.13140/RG.2.2.10294.57927, and for a discussion of some of the resulting difficulties, see [2].
Paradigm examples of electronically ideal systems are the 1920s and 1930s systems where the emitter and all aspects of the emitter-mounting arrangements are metal, and the system is operated under conditions of good high vacuum and at fields where there is no significant field-emitted vacuum space charge. Many modern systems are electronically ideal, over all or part of their voltage working range, but many others are not.
Normally, some form of validity check needs to be applied to establish whether experimental data are likely to have been derived from an electronically ideal system and hence can be validly interpreted using standard data-analysis techniques such as Fower-Nordheim (FN) plots. There currently exist three forms of validity check, namely: (i) apparent linearity of an FN or similar data-analysis plot; (ii) the Orthodoxy Test [3] and (iii) the so-called 'magic emitter test' (see ResearchGate presentation just cited). With electronically non-ideal systems, the interpretation of measured FE current-voltage characteristics can be extremely complex and often impossible to carry out exactly in the present state of research knowledge.
For electronically ideal FE systems, the interpretation problem becomes part of FE emission physics, with the value of the characteristic LVCL determined in principle by the zero-current electrostatics of the system geometry (but in practice nearly always by experiment).
We make the point that this paper is NOT directly about predicting the outcome of FE experiments of the type usually carried out. Rather, it is aimed at the question of how to design experiments (and experimental apparatus) so that useful and reliable comparisons can be made between theory and experiment. Further, it describes the start of a research journey, not the conclusion of one.
Obviously, normal scientific thinking about journeys of this kind is that one needs to avoid using experimental systems and conditions that are more complicated than necessary. In particular, one should in principle avoid experiments on emitters, such as so-called large-area field electron emitters (LAFEs), that have very many individual emission sites, with a statistical distribution of characterization parameters.
Hence, also, in a journey of this kind in FE contexts, it is desirable to only use experimental data from electronically ideal systems. Although the LVCL definition above applies to all systems, it will in this paper only be applied to electronically ideal systems (for which it can be assumed that ζ C is constant).
We also note that the parameter ζ C (which is best measured in nm) is the reciprocal of the parameter β apparently first introduced by Dyke and Trolan in 1953 [4] and best measured in m -1 .
We prefer using LVCLs rather than values of the Dyke-Trolan β, partly in order to avoid potential confusion resulting from the widespread use (in FE literature) of the symbol 'β' to denote dimensionless field enhancement factors, partly because we think a parameter measured in nm is easier to work with than one measured in m -1 . In particular, if equation (2.1) is inverted into the alternative form then ζ C has a simple physical interpretation as the factor that links the onset measured voltage V m,on to the onset value F C,on of the characteristic local field magnitude. For example, for a 4.50 eV emitter, F C,on is around 2 V nm −1 , so if ζ C had the value 1000 nm, then the onset measured voltage would be predicted as around 2000 V.

Issues related to emission physics and data analysis
For the last 90 years, since the publication in 1928 [5] of the ground-breaking (but seriously flawed [6,7]) paper of FN, most FE theory has been developed in the context of so-called smooth planar metal-like emitter (SPME) methodology. This methodology disregards the existence of atomic structure (including the possible role of surface atom orbitals in FE theory), uses Sommerfeld-type free-electron theory and treats emitters as having smooth, structureless planar surfaces of large lateral extent. However, it is possible to treat needle-shaped emitters that are 'not too sharply pointed' by using the planar emission approximation [1] and carrying out an integration of local emission current density (LECD) over the emitter surface. As is well known, the main theoretical LECD equations in SPME methodology have been: (i) the original 1928 FN FE equation (as corrected in 1929 [8]), which used an exactly triangular (ET) transmission barrier; and (ii) the 1956 Murphy-Good (MG) FE equation [9], which is based on the Schottky-Nordheim (SN) (or 'planar image rounded') transmission barrier.
The most common method of analysing I m (V m ) data has been the well-known FN plot, introduced by Stern et al. in 1929 [8] because these authors thought that an FN plot would be exactly straight. The 1956 MG FE theory is better physics than the 1928/1929 FN FE theory (see [7] for a discussion) and predicts that an FN plot will be slightly curved. This curvature makes it difficult to reach accurate quantitative conclusions from the vertical-axis intercept of the straight line fitted to the experimental data plot (see [10] and Appendix A in [11]). However, a different plot form, the so-called MG plot (see [11] for details) should (in SPME methodology) be 'very nearly straight'.
However, even for electronically ideal FE systems, it is not expected that an experimental MG plot should be exactly straight. This is because other factors, neglected in MG FE theory (and more generally in SPME methodology), are known to contribute additional voltage dependence and thus to cause departure from strict MG plot linearity. The main factors of this type are (i) the long-known fact [12] that voltage dependence in notional emission area (see below) will be introduced when an experimental or theoretical integration of LECD over the needle surface is carried out and (ii) the need to take surface atom atomic orbitals into account in a complete FE theory.
There are significant fundamental physical difficulties involved in accurately incorporating atomic structure into FE theory and this has persuaded the authors to take a strategically different approach to making comparisons between FE theory and experiment. Rather than formulating some specific advanced FE theory and then attempting to test this against experiment, there has been a proposal [13,14] that attempts should be made to fit experimental data for FE systems proven to be electronically ideal to the mathematical form where B, C and κ (also written 'k' in some past work) are parameters that in the first instance have been treated as constants. Two methods of carrying out this fitting have been discussed in [14]. Our present expectation is that C and κ will be weakly varying functions of field, and that B will be a constant or a very weakly varying function of field. Thus, depending on the data-analysis method used, the extracted values of C and κ may be average values for the range of voltages (and hence the range of characteristic field values) used.
We emphasize that equation (2.3) is a mathematical form used for fitting, not a theory of FE. In particular, it should not be assumed that exp[-B/V m ] is an expression for the Gamow factor (or 'barrier strength') that is used in FE theory. We also emphasize again that equation (2.3) is a mathematical form that is being applied in this paper to experimental data from FE systems assumed to be electronically ideal (i.e. systems for which the LVCL is constant), not to FE I m (V m ) data in general--for which the LVCL may be a function of measured voltage. Equation (2.3) has previously been called the 'empirical FE equation', but this name now seems potentially misleading. The main previous users of this mathematical form were Abbott & Henderson in 1939 [12] (who used only integral values of κ), Forbes in 2008 [13], and Forbes et al. in 2019 [14]. Both [13] and [14] allowed κ to be non-integral. Thus, this paper employs the first letters of the family names of the main past users and refers to equation (2.3), with non-integral values allowed, as the AHFP mathematical form (for measured voltage).
The parameter 'κ' has been called the 'pre-exponential voltage exponent' [14]. However, for an electronically ideal FE system, both the characteristic field F C and the related characteristic scaled-field f C (defined below) are directly proportional to V m . Thus, when AHFP-type equations are formulated in terms of these alternative variables, factors of F C κ or f C κ will appear in the equation. Thus, particularly when discussing electronically ideal FE systems, it is now considered preferable to formally call κ the AHFP equation pre-exponential parameter. However, for simplicity here, we will sometimes refer to it as the 'total κ-value'.
In principle, our long-term aim in the context of emission theory is as follows. After confirming that an FE system under investigation is electronically ideal, and then extracting experimentally derived values of B, C and κ, our aim is to use these values to deduce more detailed information about the nature of FE theory. It is expected that the values of C and κ will be more useful than the value of B, but that trying to interpret extracted values of C and κ simultaneously would be a difficult endeavour. Current efforts have therefore concentrated on issues related to the extraction and interpretation of values of κ alone.
It has been shown by Forbes et al. [14], by simulations, that the best average value of κ found (for an electronically ideal system) from a plot of ln{I m /V m Κ} versus 1/V m (a so-called power-κ plot) is a sensitive function of both (i) the details of emission theory and (ii) the shape of the field emitter. Hence, if I m (V m ) measurements of high quality were available for electronically ideal FE systems, and the effects of emitter shape and surface condition could be disentangled from these, then a measured value of κ should provide a method of comparing different FE theories with experiment. Here we report some steps towards at least assessing emitter-shape effects.
An experimental alternative would be to make reliable direct measurements of LECD, but this seems to be significantly more difficult experimentally.
As part of the overall process, we need to know: (i) what physical effects contribute to the experimental value (κ m ) of κ; (ii) what is the size of each contribution (or what are the limits within which each contribution size is likely to lie); and (iii) what information is already available in the literature. Existing knowledge about these things is reviewed in §3, where it will become apparent that a gap in our knowledge is good understanding of how the voltage dependence of the notional emission area depends on the emitter shape and on an assumed work-function value. An aim of this paper is to report the result of simulations relating to this knowledge gap. This is done in § §3 and 4. However, we need first to provide, in §3.1, some additional details of existing theory.
Finally, here, we need to point out that equation (2.3) is not the only mathematical form to have been suggested for fitting to experimental results. On the basis of a theoretical discussion based on so-called fractional-dimensional space, Zubair et al. [15] have suggested the mathematical form where V and I are a voltage and a current that appear in their theory, and a is a parameter that they call the fractional space dimensionality parameter.
royalsocietypublishing.org/journal/rsos R. Soc. Open Sci. 9: 220748 The present authors have significant hesitations about the relevance of the theoretical methodology used to derive this formula. It would be outside the scope of this paper to discuss details here, but we plan to raise them privately with Zubair et al. at a later time. However, whatever the merits of its derivation, it is clear that equation (2.4) can be regarded as an alternative mathematical form to which experimental data might be fitted.
Zubair et al. argue that this is an appropriate mathematical form for experimental data taken from rough surfaces. Equation (2.4) would also, in principle, be an alternative form for fitting data from electronically ideal 'flat-smooth-surface' FE systems of the kind that we are interested in. In this context, we make two points. First, it is possible to envisage that we shall have to explore (in sequence) a number of mathematical forms of gradually increasing complexity. Form (2.4) is more complex than form (2.3), but it is possible to argue that a better choice of 'second preferred form' might be where a and b are arbitrary mathematical parameters. Second, the best approach seems to be to thoroughly explore the usefulness of equation (2.2) before moving on to anything more complicated.

Summary of extended Murphy-Good field electron emission current-voltage theory
For our simulations, we use the extended Murphy-Good (EMG) formulation (see [11]) of FE theory. The local kernel current density for the SN barrier (J kL SN ) is given by where the symbols have their usual meanings [16]. This equation can be put into scaled format by defining three parameters: where c S [≡ (e 3 /4πε 0 ) 1/2 )] is the Schottky constant and f L is the local scaled-field for a barrier of zero-field height ϕ. Using these yields the scaled format expression: In EMG theory, the LECD J L EMG is written where the parameter λ L SN is a prediction uncertainty factor of unknown value and functional form. It formally takes account of all relevant physics not included in 1956 MG FE theory, including atomic-level wave function effects and band structure effects, and also incorporates the usual MG pre-exponential correction factor t F -2 and the MG temperature correction factor.
If the kernel current density J kL SN is integrated over the surface of a needle-shaped emitter, using the planar emission approximation, then the resulting calculated notional emission current I n SN can be written in the form where J kC SN is the value of the kernel current density at some surface location 'C' that characterizes the emitter. (In modelling, 'C' is usually taken as the emitter apex). The notional emission area A nC SN is defined by equation (3.7); it depends on the chosen emitter shape and work-function distribution, on the assumed nature of the barrier (here, the SN barrier), and on the choice of location 'C'.
royalsocietypublishing.org/journal/rsos R. Soc. Open Sci. 9: 220748 We may suppose that, for the chosen emitter model, there would be a true (but unknown) model emission current I tm , and that we can write: where λ J SN is a prediction uncertainty factor for current density, of the same general kind as λ L SN . The factor λ J SN will depend on many things, but we specifically show the barrier label. This is because the interpretation of experimental results (in particular, FN and related plots) requires an assumption about the nature of the transmission barrier.
In reality, when attempting to accurately model the 'predicted' emission current I p from a real emitter, it is inevitable (at least at present) that the model emitter will not be an accurate model for the real emitter. Hence, a second source of prediction uncertainty is introduced, and we can write where λ EM SN is a prediction uncertainty factor related to emitter-model inadequacy.
Using equation (3.7) yields Leakage current, if present, can be modelled by an additional correction factor included within λ EM SN .
We now assert that (for the electronically ideal systems under discussion) the 'formally predicted' emission current I p can be identified with the measured current I m , and hence (using equation (3.5)) we can write where f C is the scaled-field at location 'C'. In equation (3.13), I m is well defined and J kC SN is well defined in principle (if location 'C' is known, and the values of ϕ and f C are stated). Hence, the formal area A fC SN is well defined in principle. It is this formal emission area that is extracted from experiments, typically by using FN plots.
For an electronically ideal FE system, the parameter f C can also be written f C = V m /V mR, where V mR is the reference measured voltage needed to pull the top of the SN barrier down to the Fermi level, at location 'C'.
Our assumption is that, in a first approximation (at least over defined ranges of measured voltage V m ), A SN nC can be written in the form where C A and κ A are parameters associated with the shape and work-function characteristics of the emitter. As already indicated, we expect these parameters either to be constants, or (more likely) to be weakly varying functions of apex field and hence of measured voltage. The discussion above is given specifically for the SN barrier and EMG theory, but generalized versions of the theory can be given for other barrier forms.

Individual contributions
For clarity, in this section, the individual contributions to the total κ-value are denoted by the symbol Δκ. The factors currently known to affect the experimental κ-value are as follows.
Ã A. Emitter band structure. A three-dimensional free-electron model of emitter band structure, as used in deriving the 1928/1929 FN and 1956 MG FE equations, makes a contribution Δκ = 2. Models involving lower dimensionality will make, and situations involving quantum confinement may make, a contribution less than 2. It is not clearly known what the effects of more complicated real band structures, such as the 'many-valley' band structure of silicon [17], are expected to be.
royalsocietypublishing.org/journal/rsos R. Soc. Open Sci. 9: 220748 Ã B. Barrier form (i.e. shape). As has been shown in [13,14,18,19], there exist formulations of FE theory in which an SN barrier makes a contribution to κ. (For an emitter with local work function equal to 4.50 eV, Δκ ≈ -0.773.) Barriers similar in form to the SN barrier may also behave in this way. Ã C. Atomic-level structure. In Oppenheimer's approach [20,21], FE is interpreted as the field ionization (FI) of the atomic-level orbitals of surface metal atoms. The existing standard theory of the FI of a hydrogen atom [21] is equivalent to setting Δκ = -1. It might reasonably be assumed that an FE theory that in effect involved FI out of surface metal atoms might make a contribution somewhere in the range -1 ≤ Δκ < 0. Some relevant work (e.g. [20][21][22][23][24]) currently exists, but greater depth is probably needed. Ã D. Voltage dependence in the notional emission area. As discussed above, a formula for notional emission current I n will involve a notional emission area A n , and this will usually have field and voltage dependence. Details will depend on the chosen emitter shape and the work-function distribution.
There is some relevant work, e.g. [16,18,19,25,26], but no good general knowledge about the range of Δκ-values involved. The present work aims to help fill this gap. Ã E. Other physical factors included in the prediction uncertainty factor (λ J ) related to emission theory. As indicated above, correction factors related to other physical effects (for example, a temperature correction factor) are formally included in λ J . These will or may contribute to κ. Our present understanding is that--certainly for bulk metals with conventional three-dimensional band structures--contributions due to these 'other factors' are likely to be relatively small and hence can be disregarded for the time being. It is currently not clear whether this is also likely for materials with non-conventional band structures. Ã F. The prediction uncertainty factor (λ EM ) related to emitter-model inadequacy. When predicting emission currents, if an inadequate emitter model is used, then it is likely that an incorrect total κ-value will be predicted and that a correction will be necessary. At present, there is no useful classified knowledge about λ EM or this correction, and this paper does not discuss this issue in detail. Table 1 provides a summary of existing knowledge relating to the AHFP pre-exponential exponent κ. Experimental data has been included only if the related FE experimental system is known or thought to be electronically ideal.
A further account of how notional area depends on apex field has recently been presented by Ramachandran & Biswas [26]. This is discussed in §6.

The concept of notional cap-area efficiency
As already indicated, it is of considerable theoretical interest to know how the notional emission area A nC varies as a function of characteristic field F C , and how this functionality depends on the emitter shape and work-function characteristics. In practice, we use cylindrically symmetric models and take 'C' as located at the emitter model apex 'a'. Also, it is more useful to investigate the behaviour of the socalled notional cap-area efficiency g n (also called the 'area factor') given by where r a is the apex radius of curvature of the model. Further, it is better to investigate g n as a function of the apex scaled-field f a . Results will depend on what assumption is made as to the form of the transmission barrier. All discussion here assumes an SN barrier and that the planar emission approximation is being used. (Except that to compile the contributions table below, additional calculations have been carried out, using an ET barrier.) For the special case of a hemisphere-on-plane (HSP) case, Jensen ([17], see eq. (94), or [18], see eq. (30.19)) has obtained an analytical approximation, which he writes in the form g(F) ≈ 1/[b + 4-ν]. In our notation, this becomes where the bracketed acronym indicates the shape model under discussion. However, what we actually need is the form of the relationship ln{g n } versus ln{ f a }, so that we can establish the power law dependence of g n on f a and hence (for electronically ideal systems) on V m . Straightforward algebraic manipulation leads to the result  Lepetit's advanced theoretical treatment does not easily fit into the classification system used above. f The argument is also made by Biswas that the behaviour of LAFEs is much more complicated than that of single-tip emitters and involves a contribution to the main exponent. g The value 2.23 applies to an emitting surface with local work function 4.50 eV. Table 2. For the apex scaled-field range 0.15 ≤ f a ≤ 0.45 and the four emitter shapes investigated (assuming an SN barrier), the limits of the variation of the contribution κ A ( f a ) across the range and its average value (κ A ) av . Combining the individual four ranges provides an estimate of the total likely range of variation.
shape Jensen has also given a formula ( [18], see eq. (113)) for the case of a hemi-ellipsoid-on-plane (HEP). He writes his formula as g(F ) ≈ 1/[b + 1-ν], which in our notation becomes The formula is stated to hold when the aspect ratio (the ratio of the semi-major to the semi-minor axis, i.e. height divided by base-radius) is large; however, when the aspect ratio approaches unity, then factors that were neglected in the derivation come into play. This means that equation (3.21) as written is not expected to reduce to equation (3.17) when the aspect ratio becomes unity.
Since these are approximate analytical results, we have thought it helpful to check the accuracy of equation (3.17) by means of finite-element electrostatic simulations, as discussed below. We then went on to carry out numerical analyses for emitter shapes of greater practical interest, for which analytical treatment looks either impossible or impracticable.

Simulation details
In this and the following sections, for historical reasons, the contribution (to the total value of κ) made by field/voltage dependence in the notional emission area is denoted by κ A, rather than by Δκ.

Procedures
To solve the Laplace equation numerically, a commercial finite-element software package (COMSOL v5.3) was used. In addition to the HSP model, we studied three well-known emitter shape models, namely the hemisphere-on-cylindrical-post (HCP), the spherically rounded cone (SRC) and the HEP. Recently [36], it has been shown that the apex field enhancement factors for the hemi-ellipsoidal, paraboloidal and hyperboloidal (with half-angle 5°) posts are similar, so we do not consider paraboloidal and hyperboloidal posts here. Details of our simulation procedure are similar to those found in [14,36].
The emitters (other than the HSP model) are given the following geometric parameters: total height h = 4 µm and apex radius of curvature r a = 50 nm. This implies that the apex sharpness ratio σ a = h/r a = 80. The SRC emitter has an additional parameter θ, which defines the half-angle of the vertex (here, θ = 5°). The work function of the emitters was taken as constant across the surface and set to 4.50 eV.
A cylindrical simulation cell was used. The cell radius and height, the meshing geometry, and the electrostatic boundary conditions were chosen in accordance with the 'minimum domain dimensions', meshing criteria and electrostatic arguments outlined in [37]. As an illustration, figure 1a shows the resulting field distribution for a particular case of the HSP model.
Integration of current density over the emitter surface, to evaluate the emission current, was carried out directly in the software package using the planar emission approximation and the kernel current density for the SN barrier.

Consistency check with Jensen's analytical treatment
For an initial comparison between our present results for the HSP model and those of Jensen, we show (in figure 1a) our results in the form used by Jensen in Fig. 30.3 of [19], namely a direct plot of g n versus the apex field magnitude F a . (But note that Jensen uses a quantity measured in eV nm −1 to represent the behaviour of electrostatic fields: his quantity 'F' equals our 'eF'.) This is done for the three workfunction values Jensen used. Our plot and Jensen's plot are closely similar.
For ϕ = 4.50 eV, figure 1b shows a more detailed comparison. The percentage relative difference δ between the analytical and simulation results is The inset to figure 1c shows how δ varies with f a : δ decreases almost linearly as the apex field magnitude F a (and hence the scaled-field f a ) increase. Over the 'Pass' range of the Orthodoxy Test (0.15 ≤ f a ≤ 0.45), the percentage relative difference lies approximately in the range (-3% ≤ δ ≤ -1%). royalsocietypublishing.org/journal/rsos R. Soc. Open Sci. 9: 220748 Obviously, there is a question of whether this small difference is due to 'approximation error' in the analytical calculations or to 'simulation error'. However, we know from the work of de Assis & Dall'Agnol [37], on minimum simulation domain dimensions, that the simulations ought to be highly precise. Thus, our expectation is that most of the difference δ arises from the mathematical approximations made in the analytical treatment.
For ϕ = 4.50 eV, figure 1c confirms the 1/g n is a good linear function of 1/f a . A further comparison can be made between the plots of κ A ( f a ) derived from equation (3.16) and from our finite-element calculations. This is shown in figure 2. The degree of agreement is very satisfactory. More generally, figure 2 clearly predicts that κ A will be a noticeable function of f a (and hence of apex field and of measured voltage). For the circumstances of these calculations, in the range 0.15 ≤ f a ≤ 0.45, κ A takes values (for the SN barrier) in the approximate range 0.9 ≥ κ A ≥ 0.75.

Simulation results and discussion
The HSP model is the easiest to use when making comparisons between analytical models and simulations, and provides a useful degree of verification for the simulations. However, the HSP model is obviously not a realistic field emitter model. Other field emitter shapes, closer to experimental reality, are of greater interest. As indicated earlier, we have investigated three other emitter models: the HCP model, the HEP model and the spherically rounded (truncated) cone (SRC) model. For completeness, the data for the HSP model discussed in §4 are also included in diagrams below.
The current-voltage characteristics were calculated using the general procedure described in §4.1. As before, an emitter height of 4 µm and apex radius of 50 nm (corresponding to an apex sharpness ratio of 80) were used, and calculations were carried out for apex scaled-fields in the range 0.15 ≤ f a ≤ 0.45. Figure 3 shows the results of the current density integration over the emitter surface, but with the current plotted as a function of the macroscopic field F M that is applied between the top and bottom surfaces of the simulation box. The corresponding apex field enhancement factors (γ a ) are also indicated.
As illustrated by the graphs, the lowest threshold values of macroscopic field (and hence of applied voltage) correspond to sharply pointed elongated shapes, as is well known. The shapes of individual carbon nanotubes correspond most closely to the HCP model, which is the model most effective at field enhancement. This is part of the reason why carbon nanotubes have been extensively investigated as field emitters in recent years. Figure 4 shows logarithmically how the notional cap-area efficiency g n varies as a function of apex scaled-field f a , for the four emitter shapes investigated. Clearly, these results fall into two groups. It royalsocietypublishing.org/journal/rsos R. Soc. Open Sci. 9: 220748 might be expected intuitively that emitters with similar apex shapes would exhibit similar behaviour. However, the highest g n values were in fact obtained for the HCP model. This result was not obviously expected. However, closer examination shows that the HCP emitter not only has the largest apex FEF, but also has higher field values along a line along the surface, in a plane containing the model axis. This leads to greater integral current values and therefore to larger values of notional emission area (and hence to larger values of g n ).
The quantity κ A represents the contribution (to the total κ-value) that is due to emitter-shape effects. To investigate numerically how κ A varies with f a (and hence for electronically ideal emitters, with measured voltage), plots of dln{g n }/dln{ f a } versus f a were made for the studied tip shapes. These plots are shown in figure 5. Figure 5 suggests that, for all tip shapes, there will be a contribution κ A (to the total κ-value) that is a function of the apex scaled-field f a and hence of the independent variable in an AHFP-type equation. For the range 0.15 ≤ f a ≤ 0.45 under discussion, table 2 records the variation in κ A and also its average value (κ A ) av over the range (obtained by linear regression to the data plots in figure 5).

A methodology for using κ-values in theory-experiment comparisons
With the theory and simulations described, we now present a new methodology for comparing FE theory and experiment. The methodology seeks to use an experimental κ-value to decide (if possible) which of two specific FE theories is 'better science'. As an example, we shall seek to establish experimentally that the SN barrier used by Murphy and Good in 1956 is 'better science' than the ET barrier used by Fowler and Nordheim in 1928. Theoretically, there is no doubt about this (see [7]), but establishing a decisive proof from experiment has proved elusive (e.g. [38]).
We emphasize that the remarks here aim to establish a methodology and provide a 'proof of concept', but do not aim to produce a definitive result at this stage. We anticipate that getting a decisive result may be a long and messy process, but we hope that establishing this methodology may stimulate the development of techniques able to provide good experimental values of κ, with reliable error limits.    the assumption of constant local work function), but at present we have no good knowledge of effects of this kind.
In each case, on the last line of table 3, we arrive at an estimate of the range within which κ is predicted to lie if the equation in question is an adequate representation of experimental reality. These two ranges can then be used to generate figure 6, which shows how an experimental estimate κ m of the AHFP equation pre-exponential parameter κ is to be interpreted, for emitters of the kind under discussion. It is stressed that the numbers in this table are 'first estimates', based on our current theoretical assumptions, and may change as theoretical knowledge improves.
It is of interest to apply this decision table to the small number of relevant experimental κ m values listed in table 1.
The Lauritsen value of 0 (on row 2) was a highly important early empirical result but has limited accuracy: it is also of note that the emission was coming from a small area ( probably a pointed protrusion or low-work-function patch) on a cylindrical wire [28].
An objective of Stern et al.'s 1929 paper [8] was to determine whether the Lauritsen empirical value κ = 0 or the FN theoretical value κ = 2 better represented experimental reality. They concluded that the most appropriate existing dataset (taken over the largest voltage range) was that in Fig. 3 of the Millikan & Eyring's 1926 paper [29]. This dataset was also taken with a cylindrical wire emitter. Stern et al. concluded that taking κ = 2 gave the better fit to the experimental data. This conclusion was also supported (albeit less convincingly) by careful re-examination of some results of Gossling et al. [30] taken from a tungsten needle-like emitter (see Fig. 1  In 1939, the 'precise measurements' by Abbott & Henderson reported in Part II of their paper led them to conclude that the best integral value was κ m = 4 (see row 6). According to the decision table in figure 6, this result is not compatible with either FN FE theory or MG FE theory. This finding has never been satisfactorily explained and has usually been ignored. However, Abbott & Henderson were the first to point out (more than 80 years ago) that the value of κ m must depend on the 'exact microscopic shape of the point (emitter)'--as now well demonstrated by our present calculations.
Basic problems with all this early work are the limited accuracy of the experiments and the lack of realization that κ m could legitimately be non-integral (e.g. [13]).
Turning to discuss modern work on κ m , the only table entry (on row 11) is to the Ioffe Institute work on LAFEs built from carbon nanotubes. The use of LAFE data is not ideal, but these are the only experimental datasets that we currently have access to that have been controlled and recorded electronically (which is necessary for precise extraction of κ m -values). If we make the usual assumption that traditional metal emission theory is 'good enough' for CNT data analysis, then applying figure 6 thinking to the row 11 result that κ m = 1.65 yields the conclusion that these experiments are marginally indecisive.  royalsocietypublishing.org/journal/rsos R. Soc. Open Sci. 9: 220748 In summary, the existing experimental evidence, as listed in table 1, does not enable us to reach a decisive experimentally based decision that 1956 MG FE theory is a better description of experimental reality than 1928/1929 FN FE theory, notwithstanding that it is theoretically certain that the former is 'better physics' than the latter. Thus, there is a strong need for appropriate and accurate new experiments, preferably on metal emitters.

Summary and discussion
The main focus of this research has been the prediction of values for the pre-exponential parameter κ in the AHFP mathematical equation, and on preliminary comparison of predicted values with such experimental results as currently exist. So-called EMG FE theory has been outlined, six main theoretical contributions to κ (for current-voltage characteristics taken from electronically ideal FE systems) have been identified and existing reliable knowledge about these contributions (as far as we are aware of it) has been collated into table 1. A gap in this knowledge, namely the detailed effect of emitter shape on the related contribution κ A, has been filled (at least partially) by simulations on four different emitter tip shapes, for the specific common work-function value 4.50 eV. These shapes include the HSP model, for which there is an analytical treatment by Jensen. For this model, our simulations and Jensen's analytics are in very good agreement.
A new methodology for using an experimental value κ m of the pre-exponential exponent κ in the AHFP mathematical equation to decide which of two specific theories of electronically ideal FE current-voltage characteristics is 'better science' has been outlined and applied to the choice between the 1928/1929 FN and 1956 MG FE theories. Based on current theoretical knowledge, the relevant decision table (figure 6) has been constructed and applied both to historical values of κ m and to the one modern value that we consider worth discussing.
As indicated above, an interesting result relates to the Stern et al. analysis of one dataset reported by Millikan & Eyring and two datasets reported by Gossing (all taken in 1926). Stern et al. concluded that taking κ m = 2 gave lower residual error than taking κ m = 0, and hence that the results were compatible with 1928/1929 FN FE theory. As already indicated, the whole apparatus of FN plots, still in use 90 years later, was built on this finding. But the decision table in the present paper shows that the 1929 conclusion that κ m = 2 is also compatible with 1956 MG FE theory and hence indicates from experimental analysis that FN plots are not the only suitable option for analysing FE current-voltage data. (An alternative is the so-called MG plot [11], introduced on the basis of theoretical arguments.) The single modern estimate of κ m in the table tends to favour 1956 MG FE theory rather than 1928/ 1929 FN FE theory, but is marginally indecisive. Thus, FE science still remains in the state in which it has been since the late 1950s: it is theoretically certain that MG FE theory is 'better physics' than FN FE theory, but demonstrating this experimentally by a decisive argument remains elusive. The most important need looks to be for better experimental methodologies for accurate measurement of κ m for metal emitters and for determining reliable error limits.
We make again the general point that this paper represents 'first steps' toward a new methodology. As theoretical knowledge improves, and experimental results get better defined and more accurate, the decision table will change numerically (and perhaps new contributions to κ will be added), and it is expected that its interpretation will change qualitatively. When we have reached the point that this methodology can clearly show, in alignment with theoretical understanding, that 1956 MG FE theory is 'better science' than 1928/1929 FN FE theory, then it will be timely to construct one or two new decision tables that compare 1956 MG FE theory with one or two more sophisticated FE theories-for example, the 2014 Kyritsakis-Xanthakis theory [39] of emission from an earthed sphere.
We note that our simulations have so far been carried out only for one value of work function and one value of apex sharpness ratio. This has been sufficient to establish a proof of concept for our new methodology, but we also plan to carry out further simulations for ranges of values of these parameters.
Finally, two recent papers by other authors need comment. For a model involving a hemi-ellipsoidon-cylindrical-post (HECP model), Ramachandran & Biswas [26] have investigated how the dependence of the notional emission area on the apex field varies with the ratio 'hemi-ellipsoid height divided by apex-radius'. Detailed investigation of the relationship between their results and ours is outside the scope of this paper, but our intention is to incorporate an appropriately customized version of their results into updated versions of tables 2-4, in due course.
In the second recent paper, Ayari et al. [40] discussed whether methodology based on measuring κ m could be used to extract a value for the local work function ϕ from FE current-voltage measurements. We royalsocietypublishing.org/journal/rsos R. Soc. Open Sci. 9: 220748 agree with their conclusion that this would be a very challenging task, but we think that the discussion in the present paper shows that the complete theoretical picture is more complex than their arguments have assumed (although the conclusion would be the same).
Finally, we state again our belief that the decision table methodology introduced here could become a useful tool in developing better FE science. Data accessibility. Computations described in this paper have been carried out using two standard commercial Software packages. COMSOL has been used to carry out electrostatic modelling. LABVIEW has been used for calculation of emission current densities and emission currents, and related analyses.
When COMSOL is used, a large internal data file is generated by COMSOL and passed to LABVIEW. Appropriate modelling information (as described in the main text) and appropriate sequences of input instructions have to be entered into the software packages, but there are no large data input files and there has been no need for low-level code to be written. It is assumed that researchers able to drive the software packages used (or packages with equivalent functionality) are able to generate the package input instructions needed to implement the procedures described in the main text or in the electronic supplementary material, and to generate associated diagrams.
Additional details concerning the electrostatic modelling, beyond those provided in the main text, are essentially similar to those provided as electronic supplementary material to [11], which can be found via the link: https://doi. org/10.6084/m9.figshare.c.5324958.
The main difference is that in the present paper four different emitter shapes were analysed, rather than the single shape treated in [11].