Periodic energy decomposition analysis for electronic transport studies as a tool for atomic scale device manufacturing

Atomic scale manufacturing is a necessity of the future to develop atomic scale devices with high precision. A different perspective of the quantum realm, which includes the tunnelling effect, leakage current at the atomic-scale, Coulomb blockade and Kondo effect, is inevitable for the fabrication and hence, the mass production of these devices. For these atomic-scale device development, molecular level devices must be fabricated. Proper theoretical studies could be an aid towards the experimental realities. Electronic transport studies are the basis to realise and interpret the problems happening at this minute scale. Keeping these in mind, we present a periodic energy decomposition analysis (pEDA) of two potential candidates for moletronics: phthalocyanines and porphyrins, by placing them over gold substrate cleaved at the (111) plane to study the adsorption and interaction at the interface and then, to study their application as a channel between two electrodes, thereby, providing a link between pEDA and electronic transport studies. pEDA provides information regarding the bond strength and the contribution of electrostatic energy, Pauli’s energy, orbital energy and the orbital interactions. Combining this analysis with electronic transport studies can provide novel directions for atomic/close-to-atomic-scale manufacturing (ACSM). Literature survey shows that this is the first work which establishes a link between pEDA and electronic transport studies and a detailed pEDA study on the above stated molecules. The results show that among the molecules studied, porphyrins are more adsorbable over gold substrate and conducting across a molecular junction than phthalocyanines, even though both molecules show a similarity in adsorption and conduction when a terminal thiol linker is attached. A further observation establishes the importance of attractive terms, which includes interaction, orbital and electrostatic energies, in correlating the pEDA study with the transport properties. By progressing this research, further developments could be possible in atomic-scale manufacturing in the future.

electronic devices [1]. These molecular devices are the building blocks for atomic/close-to-atomic-scale manufacturing. The foundation must be laid so that a mass production of these devices could be developed in the future. The electronic transport studies across a molecular junction can lead to this perspective. Many studies have been performed in revealing the transport properties through a molecule attached to the ends of lead electrodes [2][3][4][5]. A recent technical article shows the importance of contact between the molecule and the electrodes and how it can affect the transport of electrons [6]. Hence, the bonding properties and the interactions at the junction are an inevitable factor to be considered. Even though the transport studies are well established, a quantified analysis that enhances this study is still lacking in the present research fields and it still remains as a research gap. To fill this gap, a new approach called periodic energy decomposition analysis (pEDA), introduced by Raupach and Tonner [7], which is an extension of energy decomposition analysis (EDA) by the consideration of periodic boundary conditions (PBC), could be used to have a detailed study on the bonding interactions between two fragments in a periodic system. Through this study, the total energy by which the molecules get attached to the electrodes could be decomposed into interaction energy, Pauli's repulsion energy, orbital energy and orbital interactions, which could be analysed separately to have a deeper knowledge on the molecule-substrate interactions. By extracting the information regarding interactions occurring between the molecule and the substrate, a connection can be established between the bond strength and the transport properties which forms the motivation for this study.
Phthalocyanine (H 2 Pc), tin (II) phthalocyanine (SnPc), porphyrin (Pr) and tin (II) porphyrin (SnPr) have proved to be excellent candidates for the molecular electronic studies [8,9], hence, they are selected for this research. Phthalocyanines are blue-green coloured aromatic compounds which have applications of molecular switching [10]. This application is due to the rapid charge transfer, less than a picosecond, across the molecule [11]. Apart from that, phthalocyanines and similar porphyrins exhibit stability in most physical and chemical conditions. Phthalocyanines were discovered by Linstead in 1934 [10,12] and the first preparation of porphyrin, called hematoporphyrin, by Johann Ludwig Wilhelm Thudichum was between 1829 and 1901 [13]. The optimized structures of H 2 Pc and Pr are given in figure 1. Similarly, the structures of SnPc and SnPr are given in figure 2. The optimized bond lengths of molecules under study are given in supplementary material (see figures S1-S4, available online at stacks.iop.org/IJEM/2/015401/mmedia), which are comparable with that reported by Zhang et al, albeit some bond lengths differ due to the different exchange-correlation functionals and the software used, details of which can be found in [14].
Orbital overlapping takes place when a molecule is incorporated between the electrodes. This overlapping leads to chemical bonding which in turn determines the coupling strength at the junction. Therefore, importance of chemical bonding should be considered since it forms the basis for molecular chemistry calculations [15]. The effects of chemical bonds between alkyl chains and the metal junctions on tunnelling mechanism have been already established [16]. In quantum calculations, wavefunctions and delocalized electronic structures help to calculate the probability of electrons to present in a particular orbit, whereas the chemical bonding properties are more localized and needed to be established. A summary on different methods to analyse the electronic structures of the compounds can be found in [7].  To establish the link between the above mentioned study and electronic transport analysis, we have presented transport calculations by placing the molecules under study between two gold electrodes by employing non-equilibrium Green's function (NEGF), thereby, obtaining transmission spectrum and current (in 2e/h), where e is the electronic charge and h is Planck's constant Vs bias potential profile. Section 2 gives the theories required for the pEDA and NEGF calculations, which are referred from literatures and the documentations of the software used (ADF 2018). The material and methods including the set-ups used for the study are given in section 3. The results obtained, necessary interpretations and analysis are described in section 4 and finally the concluding remarks are provided in section 5. Throughout the context, the terminology 'phthalocyanines' generally represents H 2 Pc and SnPc. Likewise, the term 'porphyrins' represents SnPr and Pr.

Periodic energy decomposition analysis
The EDA was first studied by Kituara and Morukuma [17] and further developed by Zigler and Rauk [18,19]. Here, intrinsic bond energies of two fragments, A and B, forming AB molecule, are considered and the bond formation is divided into different sub-steps. A summarized theory is given, details of which can be found in [7,20]. A schematic representation of the EDA process and the corresponding energy terms are given in figure 3.
From figure 3, the bond dissociation energy, E bond , is the combination of preparation energy, E prep , and the interaction energy, E int , given as The preparation energy is required since the fragments get excited from its ground state to a reference point when they are geometrically distorted and possibly electronically excited , is the energy of the fragments, A and B, in the ground states respectively and E A B , is the energy of the distorted fragments, respectively. The fragments interact to form the AB molecule and the intrinsic bond energy would be the difference of the combined molecule and the individual fragment energies, i.e. the interaction energy would be partitioned into different terms such as electrostatic energy (E elstat ), Pauli repulsion energy (E Pauli ) and the orbital relaxation energy (E orb ), given as where N is the normalization factor and Â is the antisymmetrization factor. Finally, the orbitals of y¢, which are frozen, can be made to relax and the optimal wavefunction y AB is found with energy E , AB which is the orbital energy, where r AB is the density of y AB and r¢, the density of y¢.
It is important to add a dispersion term to the abovementioned theory which improves the performance of the calculations, which is made possible by implementing the DFT-D3 functional put forth by Grimme et al [21]. Hence, the interaction energy can be further be expressed as the combination of dispersion effect (disp) and the electronic effect (elec) [20], given as The different energy terms are given briefly without going into details [7]. The electrostatic energy term can be given as is the atomic numbers of fragments A and B, respectively, a b R , is the nuclear positions of two fragments, r i j , is the position of electrons of both the fragments and r A B , represents the charge density of fragments A and B, respectively. After the inclusion of orthogonalized wavefunction, y¢ (k), and the expansion of bloch basis, the Pauli repulsion energy term can be given as where D ¢ T represents the occupied and the virtual space of the transformation matrix, DV P Coul , is the change in coulomb interaction due to fragment orthogonalization and D ¢ V XC gives the exchange-correlation energy. Finally, after relaxing the intermediate wavefunction, y¢, and obtaining the final wavefunction, y , AB the orbital relaxation energy can be found out by the equation orb AB All these energy terms can be calculated separately using the above equations.

Non-equilibrium Green's function (NEGF)
NEGF can be used for modelling the transport properties across a molecular junction. For the electronic transport, a voltage should be provided at the drain electrode. For a threeterminal device, voltage at the gate controls the current flow across the junction. When the voltage is provided, there would be a difference in the potentials on both sides and electrons start to flow from the source to drain. This can be solved by using transport equations, which can be used to calculate the electron density and current due to the potential difference [22]. When the electrons start to flow from the electrodes, coherent scattering takes place, and this can be modelled using software. Here, we have modelled this transport properties using SCM DFTB [23]. From this software, it is possible to calculate the transmission function, T(E), which gives the number of electrons with energy, E, transferred across the junction filtering through the scattering region. The current across the junction can be calculated using the Landauer-Buttiker formula is the Fermi-Dirac functions for the left and right electrodes with the chemical potential ( ) m . L R As mentioned above, the transmission function can be found out from the Green's function, G(E), which can be given as Once the Green's function is calculated using the above equations, transmission function can be found out by taking the trace of G(E) and the broadening matrices of two electrodes, G L and G , R which can be given as where the broadening matrices, which are the anti-Hermitian part of the self-energies [24], are given by where I is the imaginary part.

Materials and methods
All calculations were performed by the Amsterdam Density Functional (ADF) molecular modelling suit programme using the BAND [25] and DFTB [23] platforms. All molecules were optimized before performing the calculations. Relativistic corrections were included by means of zeroth-order regular approximation (ZORA) formalism [26]. The exchange-correlation (XC) functional of Perdew [27,28] and Becke [29,30] were used. Triple-zeta polarization (TZP) basis sets were used to describe the valence electrons of carbon, hydrogen, nitrogen, sulfur and tin. For the pEDA calculation, generalised gradient approximation functional with dispersion correction, PBE-D3 [21], was used to add the dispersion term and to carry out the calculations in three dimensional PBC. The frozen core was set to small and normal numerical quality is used in BAND calculation, in order to reduce the calculation time [7]. Regular k-space grids of 1×1 and 5×5 were employed for the calculations of pEDA energy terms. All calculations were assumed to be performed in vacuum since the air resistance is neglected by default in the software. The adsorption of H 2 Pc, SnPc, Pr and SnPr on the Au (111) surface with and without thiol linkers were the set-ups used for the pEDA studies (coordinates of which are given in supplementary tables S3-S10). The convergence criteria for the SCF and the structural optimization calculations are provided in supplementary table S2. A large super cell with 2×3 dimension and three layers of Au for Pr, SnPr and a super cell of 4×4 dimension with three layers of Au for H 2 Pc and SnPc were used to set-up the pEDA calculation, since the latter molecules have a larger dimension. The Au substrate was described as a slab with two dimensional PBC and theoretically optimized lattice parameter of 4.159 Å [31] was used. The Au substrate was optimized using DFTB platform employing DFTB.org parameter directory [32]. For the pEDA calculations without the thiol linker, the terminal hydrogen was attached to the gold with a bond length of 1.980 Å as shown in scheme 1(A). For the molecules with the thiol linker, the bond length between the terminal sulfur atom and the gold was set to 2.611 Å as shown in scheme 1(B). The mentioned bond lengths were chosen in order to avoid the overlap error and to provide enough gap between the atoms. The same results will be obtained even if no bond connection is provided, but maintaining the mentioned distance, since the DFT calculations are independent on the bond characteristics. Two separate fragments were selected and the overlapping of one or more atoms were carefully avoided to rectify the errors of used fragments. The convergence was found to be slow for SnPr, which was rectified by turning on the degenerate occupancy. All calculations were performed thrice to test the reproducibility.

Results and discussion
The highest occupied molecular orbitals (HOMO) and the lowest unoccupied molecular orbitals (LUMO) of all the molecules are obtained from ORCA [33] (supplementary figures S5-S8). The different HOMO and LUMO values along with the HOMO-LUMO gap are given in supplementary table S1 for reference.
The pEDA calculations were performed using the parameters mentioned in section 3. The simulated calculations were quicker for Pr and slower for H 2 Pc than the other molecules, for both k-space grids, since the size is lesser for the former than the latter. The pEDA energy terms were obtained that contributed to the bonding interactions between the molecule without thiol linker and the substrate which is given in table 1. A detailed explanation of each term is beyond the scope of this paper since we have used the pEDA values only for the aid of the following studies and to establish the link between pEDA and electronic transport properties. Readers are advised to refer to recent reviews [20,34] and a study [35] on pEDA analysis to gather an overall idea of the energy terms. In the case of 1×1 k-space grid, from table 1, phthalocyanines have a fewer absolute values of !E elstat and !E orb (except for H 2 Pc), which are due to the low back bonding of electrons from Au surface to the molecule.
The increased !E orb is the only stability available for H 2 Pc. At the same time, for porphyrins, an average increase of ∼17 kJ mol −1 in !E orb can be observed. Special attention should be given to !E int , which is almost comparable for H 2 Pc, SnPc and Pr but is larger for SnPr. Since the former three have a comparable !E int , they might be interpreted to have a similar bonding scenario, but evaluating !E Pauli , which is much higher, about an average of 1.2 times, for porphyrins when compared with phthalocyanines, they appear to lack stability. However, this destabilization is compensated by the increased dispersive effect (about an average increase Phthalocyanine of 14.25 kJ mol −1 ), along with the higher attractive terms (!E elstat and !E orb ) make them more interactive with the surface electrons. It should be observed that the electronic interaction energy, !E int (elec), is actually repulsive but is overruled by the much higher dispersion effects. Hence, a low interaction energy, a weaker orbital and electrostatic energies make phthalocyanines poorly interactive with the gold substrate compared with porphyrins. When the attaching hydrogen was replaced with a single thiol linker, as expected, the adsorption properties got enhanced since the thiol functional group has lone pairs and makes the system more electron rich. Tivanski et al [36] have shown the influence of thiol linkers on the conductivity of molecules over the gold substrate. They have reported the enhanced conductivity with double thiol linkers than single but with different molecules, biphenyl −4,4′-dithiol and biphenyl −4,4′-dicarbodithioic acid. Here, in our study, we have used a single thiol linker to make the calculations simpler. Table 2 shows the corresponding pEDA energy terms for the molecules with the thiol linker. From the table, the contributions of interaction energy from both electronic and dispersive effects as well as, all the pEDA energy terms, repulsive (!E Pauli ) and attractive (!E elstat , !E orb ), are comparable. Hence, the bonding in these molecules are identical and shows almost an equal conduction across a molecular junction. It should be noted that !E Pauli is higher for thiol linked molecules than without and yet shows a high conduction and interaction. The reason for this is the increased dispersion effects and the attractive interactions when a thiol linker is attached. Analysing each molecule with thiol linker, H 2 Pc has ∼36%, SnPc, ∼38%, Pr, ∼32% and SnPr, ∼27% higher electrostatic and H 2 Pc ∼10%, SnPc ∼28%, Pr and SnPr ∼14% orbital contribution than without thiol linker. This dominance increases the back bonding of electrons from the surface to the molecule. Also, the distance from the terminal sulfur to gold contributes to the increased !E Pauli . The influence of distance, yaw and pitch on pEDA energy terms can be referred from [7,20]. Additionally, the increase in size and steric demand of thiol linked molecules contribute to the increased conduction. However, even though the values are comparable, an increased attractive term make porphyrins a slightly more conductive than phthalocyanines.
The pEDA values when performed in a regular k-space grid of 5×5 and a tight k-mesh in correspondence with the lattice space vector, to have more reliable results, are also shown in tables 1 and 2. Here, for molecules without thiol linker (table 1), !E Pauli , is found to decrease for H 2 Pc while an increase in !E Pauli term can be seen for other molecules and for all molecules with thiol linker (table 2). The !E int term got decreased for all molecules in both configurations except for Pr without thiol linker. All the other terms are almost comparable with no change in dispersion energy terms.
Our findings show that the attractive terms that include !E int , !E elstat and !E orb contribute mainly in correlating the transport properties. The !E int (disp) and !E Pauli terms influence the formers in altering their effects. The electronic transport studies performed in the following are in favour of the results obtained with the pEDA calculation. For all these calculations, the relative contributions of dispersion effect to the interaction energy is higher than the electronic effect that can compensate for the other energy terms in decreasing/ increasing the interactions in terms of the dispersion values accordingly. Additionally, an increased size of the adsorbate can also have an impact on the k-space sampling.
To take this study further, the transport properties were studied for all the molecules by placing them in between two 3×3 gold electrodes cleaved at the (111) plane. A single point contact is provided between the terminal atom and the electrode (scheme 2). The calculations employed the usage of DFT and NEGF and were performed using SCM DFTB platform. For the NEGF calculations, QUASINANO2013.1 [37] parameter directory was used. To maintain the consistency, distance between terminal hydrogens and gold atoms were set to 1.980 Å, for molecules without thiol linkers, and the terminal sulfur and gold atoms to 2.611 Å, for thiol linked molecules, (as shown in scheme 2) same as given Phthalocyanine in pEDA calculations (see scheme 1). The required theory for the same is given in section 2.2. With this, a transmission spectrum and current (in 2e/h)-energy characteristic curves can be obtained. The spectra of molecules without thiol linkers at different gate voltages along with the device set up can be referred from supplementary figures S9-S16. From the transmission spectra of each device set-up without thiol linker, H 2 Pc has comparably very higher peaks above the Fermi level, positioned at 0 eV, which represents the LUMO levels (supplementary figure S9). This increase in peaks are due to the higher !E orb obtained for H 2 Pc from the pEDA calculation. However, almost similar levels of peaks are obtained below and above the Fermi level of SnPc (supplementary figure S10). For Pr (supplementary figure S11) and SnPr (supplementary figure S12), peaks are almost equally distributed on both sides of the Fermi level. For different voltages, peaks increase with an increase in the applied gate voltage. Even though the peaks seem to be less distinct for porphyrins, they have higher areas of transmission making them more transmission favourable than phthalocyanines, which is confirmed in I versus E curves ( figure 4(B)). Hence, among the molecules under study, without thiol linkers, as obtained from the pEDA calculations, porphyrins are best adsorbed over the gold substrate which led to higher transmission and conduction of current across Au-molecule-Au molecular junction.
The device set-up and the transmission for different gate voltages when the molecules under study are attached with a thiol linker on both sides can be referred from supplementary figures S17-S24. From the transmission peaks obtained, there is a considerable increase in the conductivity of all molecules removing the energy window. This happens since there are very high overlap of the orbitals, the change in their symmetries and an increased electrostatic influence as explained with the pEDA results. Due to this, evidently, leakage current is much pronounced with the thiol linked molecules. A slight change in the position of the contact could alter the peaks in a different way since this affects the position of the energy levels [38]. It should be noted here that the transmission at  different gate voltages for all molecules are almost similar above and below the Fermi levels. Even though this is the case, a slightly higher peaks can be visible for porphyrin molecules (supplementary figures S19 and S20). A comparison of transmission for all molecules with thiol linker at zero gate voltage in figure 4(C) also proves this situation. Apart from that, the current (in 2e/h) versus energy curves for the molecules at zero gate voltage in figure 4(D) also shows a consistent result.
From figure 4(A), for molecules without thiol linker, the transmission peaks are very well visible for Pr and SnPr. Also, the current conduction is highest for SnPr and the least for Pc showing better conduction for porphyrin and are consistent with the pEDA calculation. This happens because the current is directly proportional to the area under the transmission peaks which is highest for SnPr and least for Pc as shown in figures 4(A) and (B). For molecules with thiol linkers, even though all molecules show similar conductivity and transmission, porphyrins have slightly higher transmission ( figure 4(C)). In addition, the current versus energy plot for porphyrins are slightly enhanced ( figure 4(D)).
For the next generation manufacturing, Manufacturing III, the functional features must be in the atomic scale [1,39]. For this implementation, a stable and steady method for the mass production of atomic scale electronic components are needed. The best and foreseeable method for the mentioned is the application of self-assembling molecules [9]. Phthalocyanines [40] and porphyrins [41] are tested for their selfassembling capabilities for biological and molecular electronic purposes, out of which, porphyrins were self-assembled very fast, within 15 min, and phthalocyanines, in 30 min, on ZnO thin film [42]. With the studies performed here, further molecules could be analysed and calculated theoretically. Once a solid method is achieved for durable transistor actions with these molecules, future studies can be concentrated on connecting the molecular leads without altering the contact points with the outside circuitry. Also, a consistent and solid current conduction can be made possible with a strong contact, thereby, linking the microscopic with the macroscopic world.

Conclusions
In conclusion, we have presented a periodic EDA on four different molecules, derivatives of phthalocyanines and porphyrins, to prove their adsorption capabilities over a gold substrate. Through this study, we were able to provide an idea regarding separate energy terms contributing to the bond formation of the mentioned molecules at a molecular junction. Besides, the electronic transport properties across each of the molecules were studied and analysed. Porphyrins were found to be much more interactive with the gold surface according to the pEDA results which are confirmed by the consistent favourable results obtained with NEGF calculations, thereby, we have established a connection between the pEDA and electronic transport studies across a molecular junction. All results corresponding to both studies were consistent throughout. Hence, we have successfully quantified a method to correlate the bonding analysis with transport properties across a molecular junction. Apart from that, comments on further studies of these molecules such as the self-assembly and the stability of contact, a mass production of atomic-scale devices with these molecules are provided, which would be of help in forming the foundation of Manufacturing III in the future. Furthermore, questions such as 'is this connection between pEDA and electronic transport studies consistent for all the molecules?', 'will this study help to identify potential molecules for the fabrication of molecular devices?', 'will this help to authentically study the interactions at the metalmolecule interface?', 'is experimental verification possible?', and 'would surface planes other than (111) gold substrate and a higher k-space grid enhance the study?' need to be addressed in the future.