High-throughput Identification and Characterization of Two-dimensional Materials using Density functional theory

We introduce a simple criterion to identify two-dimensional (2D) materials based on the comparison between experimental lattice constants and lattice constants mainly obtained from Materials-Project (MP) density functional theory (DFT) calculation repository. Specifically, if the relative difference between the two lattice constants for a specific material is greater than or equal to 5%, we predict them to be good candidates for 2D materials. We have predicted at least 1356 such 2D materials. For all the systems satisfying our criterion, we manually create single layer systems and calculate their energetics, structural, electronic, and elastic properties for both the bulk and the single layer cases. Currently the database consists of 1012 bulk and 430 single layer materials, of which 371 systems are common to bulk and single layer. The rest of calculations are underway. To validate our criterion, we calculated the exfoliation energy of the suggested layered materials, and we found that in 88.9% of the cases the currently accepted criterion for exfoliation was satisfied. Also, using molybdenum telluride as a test case, we performed X-ray diffraction and Raman scattering experiments to benchmark our calculations and understand their applicability and limitations. The data is publicly available at the website http://www.ctcms.nist.gov/~knc6/JVASP.html.

using the PBE functional, which has been found to be non-suitable for characterizing 2D materials and require more sophisticated vdW functionals than just using PBE functional. Ataca et al. 30 also performed stability analysis to identify stable compounds among monolayer materials using an LDA functional, and they found 52 such materials. For all such compounds they also calculated elastic coefficients using LDA but only in the 2H and 1T hexagonal monolayer structures. Rasmussen et al. 22 studied the electronic structures of 51 semiconducting monolayer transition-metal dichalcogenides and -oxides in the 2H and 1T hexagonal phases using DFT and taking into account many-body effects using GW approximation to alleviate DFT band-gap problem. Although this work gives a huge thrust towards 2D materials repository, it consists of prototype structures (2H and 1T) only, and their experimental fabrication is a debatable issue. In addition, it didn't provide a robust method to identify 2D materials among all other bulk materials. Moreover, most of their computations are carried out with PBE for geometric optimization within DFT. Although, GW can be reinforced to obtain better band-gap predictions, important effects such as strain dependence of the bandgap and other properties are not reliable for 2D materials when computed with PBE. Very recently, Materialsweb repository has found hundreds of 2D materials using topological scaling algorithm 31 . Similar algorithm has been implemented by Mounet et al. 32 to identify novel layered materials. DFT databases such as materials project (MP) 33 , Open Quantum Materials Database (OQMD) 34 , AFLOW 35 have performed thousands of material data calculations with the PBE functional. These databases are quite homogenous among themselves and encompasses many of the ICSD materials data. While they provide a large amount of electronic structure, surface, interface and mechanical properties data, the choice of the PBE functional, fixed plane-wave energy cutoff and fixed K-point selection may not be well suitable for characterizing 2D materials 36 . Many 2D materials such as graphene require plane wave cut-off more than 600 eV 37 , but present DFT databases takes 520 eV cut-off for all materials uniformly. Additionally, these DFT databases do not contain any information for layered (such as single and bi-layer) material properties.
As evident from the above discussion, while there are many DFT-based databases available nowadays, there is still a significant need for a systematic evaluation of 2D material properties in monolayer, bilayer, and multi-layer forms, and for a comparison to their 3D bulk counterpart. Such properties should be calculated for all experimentally observed structures, not just conventional 2H and 1T, using DFT-functionals more suitable to the 2D case, such as optB88 functional 38,39 . Also, sufficient convergence in plane-wave energy cut-off and K-point mesh is needed for reliable results. Lastly, a computationally-light criterion to identify whether a 3D material can be made into a 2D is still missing.
In this work, we developed such a criterion and then computed a vast array of properties for all the materials that our criterion identified as good 2D candidates. Specifically, our screening methodology is based on combining the finding from the Ataca et al. 30 work and existing DFT databases. We notice above that PBE generally overestimates the lattice constant, especially along the direction of vdW interactions. We decided to use such a peculiarity as a basis for screening of materials. All major DFT databases use PBE and they also all contain the corresponding experimental lattice constants (from ICSD databases), for comparison purposes. So, we simply extract the PBE and the experimental lattice constants (a, b and c), and compute their relative error. The relative error in a or b or c-axis equal or greater than 5% is chosen as a preliminary screening criterion. After the preliminary identification, we repeat the DFT procedure for all these materials with better vdW functionals, such as optB88 to compute energetics and all other properties mentioned above.
For all possible 2D material candidates, our database provides materials properties for single layer as well as for their bulk counterparts. While the properties of many possible 2D materials have already been posted into our database, many of our calculations are still running, and their results will be automatically updated to the website once completed. All the possible 2D material candidates that we found are described in the results section below. The properties in our database include energetics, geometrical (computational X-ray diffraction) and electronic properties (band structure, density of states and work-function) and elastic constants for both the bulk and the single layer counterparts of materials. Most of the bulk materials have Inorganic Crystal Structure Database (ICSD) IDs 40 , so experimentalists can easily compare collected XRD patterns and look up how the actual materials are fabricated. Validation for all the 2D materials data with previously published work is beyond the scope of this paper, however we will discuss few examples of common 2D materials as benchmarking of the database. The data is distributed on National Institute of Standards and Technology (NIST) webpage for public use at http://www.ctcms.nist.gov/~knc6/JVASP.html. At present our database consists of 1012 bulk materials and 430 single-layer materials, and the database is still populating.

Results
Screening criterion and classification of vdW materials. As mentioned above, semi-local density functionals such as PBE are known to give incorrect results of lattice constants for vdW-bonded materials. All the publicly available DFT databases such as MP, AFLOW and OQMD uses PBE for property calculation. Hence, materials with significant error in lattice constants can be expected to be vdW bonded. We chose δ greater than or equal to 5% as screening criterion (for details see the discussion section) to identify possible 2D materials, where δ is the relative difference between ICSD (experimental data) and MP-PBE (DFT database computed using the PBE functional) lattice constants (a, b and c lattice constants chosen as in the ICSD and MP databases) for bulk materials in non-cubic crystal systems: PBE ICSD ICSD Using this screening criterion as an initial step, we identified 1356 possible 2D materials out of all the systems in the Materials Project database. Please note that we didn't perform any of the DFT calculations needed to compute δ, we simply obtained the necessary quantities from the Materials Project database through REST API 41 .
Scientific RepoRts | 7: 5179 | DOI:10.1038/s41598-017-05402-0 An overall picture of the classification of these possible 2D materials is given in Fig. 1. To begin with, we found significant relative error in all lattice constants as well as their combinations (Fig. 1a). Specifically, the Venn diagram shows that the most of the materials have a relative error of 5% or more in lattice constant a (868 systems), followed by lattice constant b (679 systems) and, lastly, lattice constant c (651 systems). Relative errors equal or greater than 5% in all three lattice constants (a, b, and c) were found for 260 materials, while similar relative errors in only a, only b, and only c were found for 461, 139 and 174 materials, respectively. Lastly, relative errors equal or greater than 5% in both a as well as b, a as well as c and b as well as c were found for 105, 42 and 175 materials, respectively. We clearly observe that the Vander Waal bonding can have effect not just in z but also in x and y directions of the crystal.
With respect to the nature of these possible 2D materials, we found that most of the layered materials are chalcogenide in nature (Fig. 1b). However, halide, pnictides and their combinations are also possibilities. In Fig. 1c we show the structural prototype distribution of the materials in our database. While the major contribution is from AB 2 structures, most of which are already known and have been fabricated 17 , other dominating stoichiometry are AB, ABC, ABC 2 and AB 3 . These are less known potential 2D materials, possibly deserving further investigation. More details on the materials belonging to different structure-prototypes are given in the supplementary information (Table S1), where each prototype is listed with the corresponding materials. Like the prototype classification, in Fig. 1d we show the space group distribution of the predicted layered structures in our database. Space groups can change between bulk materials and corresponding single layer ones due to the broken   Fig. 1d shows that 2D materials can have counterintuitive crystal structure space groups rather than just 2H and 1T structures. Recently, quite a few 2D materials have been experimentally fabricated which have crystal structures beyond 2H and 1T, such semi-metal WTe 2 (space group Pmn2 1 ) 42,43 . In single layer materials space groups have direct implication on their properties. For instance, crystals with an inversion center cannot display piezoelectric behavior, and hence our space group results could be used to screen for 2D piezoelectric materials. Such a screening criterion has been already used by Materials Project to identify potential bulk piezoelectric materials 44 . Additionally, crystal structure space group information can be used in Raman experiments and second-harmonic generation (SHG) experiments, in which the spectral peaks are sensitive to space groups 45 . Interestingly, materials with lack of inversion symmetry can lead to valley-selective excitation of charge carriers 46,47 . For this purpose, a detailed list of the materials in various space-group is given in supplementary information (Table S2). We also provide the bulk and layered space group information in our database so that it can be used in characterizing and identifying the materials experimentally.
In Fig. 1e we show the crystal system distribution of the bulk layered materials. Most of the layered materials belong to orthorhombic crystals system. Other major crystal systems are monoclinic and tetragonal systems. Up to this work bulk layered materials were mostly considered as having honey-comb structures in hexagonal system. However, our results clearly show that layered materials often belong to other crystal systems as well. A complete list of predicted layered materials in different crystal-systems are given in supplementary information (Table S3). In Fig. 1f, we show that most of our predicted layered materials are binary and ternary compounds. Our results are similar to the work by Ashton et al. 31 . As binary and ternary materials are generally easer to fabricate experimentally than others, the ICSD database is biased towards such compounds. Specifically, it contains 2,033 elementals, 34,785 binaries, 68,730 ternaries, 68,083 among quaternaries and quintenary compounds, as of year 2016-17. As our database imports structures from the Materials Project, and the Materials Project was built on the ICSD database, the same bias towards binary and ternary compounds can be observed.
As the choice of a 5% relative-error for our screening criterion is somewhat arbitrary, we recomputed all the quantities discussed in Fig. 1 using cutoff values of 3%, 7% and 10%. The complete results of these calculations are shown in Fig. S1 and in Tables S4-S5 in the supplementary information. As expected, the total number of screened compounds decreases as δ increases (3116, 1356, 819 and 375, for δ equal to 3%, 5%, 7% and 10%, respectively). However, in all cases the largest number of materials satisfying the screening criterion has large relative error in lattice constant a (Table S4). Similarly, the number of compounds with large relative error in lattice constant b is in all cases very close to the number of compounds with large relative error in lattice constant c (Table S4). We find similar behavior for chemical composition and stoichiometry. For all δ, the most prevalent chemical composition is chalcogenides, followed by halides and pnictides (Table S5). AB 2 is always the most common crystal prototype, and AB, ABC are the second and third most likely crystal prototypes in all cases except for δ = 3% (Fig. S1). It is worth pointing out that the probability of the selected material to be in the AB 2 form grows monotonically as δ increases. The distribution of space groups shows less systematic behavior among all the investigated quantities. Orthorhombic crystal system is always the most commonly found structure, and monoclinic is the second most likely in all cases except for δ = 10%, where it appears as the third most likely one.
Conventionally, the accepted error in lattice constants computed using DFT with the PBE functional is about 2% 36,48 with respect to experimental values. As we found a lot of compounds with relative error in lattice constants greater than 5% (Fig. 1a), we noticed that PBE functional may not be the ideal functional to investigate 2D materials, and tested if the optB88 exchange-correlation functional 38, 39 could be a better choice. We therefore relaxed the 1356 materials that our screening criterion had identified using the optB88 functional. Our results comparing lattice constants from ICSD to those obtained using PBE (Materials-Project data) and our optB88 functional-based data are shown in Fig. 2.
ICSD-data are given along the x-axis, and computational data are given along the y-axis (Fig. 2a-c and j-l). The DFT data is plotted using dots, while the solid lines indicate the results for a perfect agreement between experimental and computational DFT data. We find some degree of deviation from ideality in all cases, but such deviation is clearly more significant for PBE data than optB88. To better quantify such a deviation, we calculated the mean absolute error (MAE) in the data using equation 2: where x i is the optB88-based or the PBE-based DFT data and y i is the corresponding ICSD data. Comparing MAE values, we find that optB88-MAE are always lower than PBE-MAE, up to 4 times lower in the case of lattice constant a ( Fig. 2(a) and (g)), validating the assumption that the optB88 functional describes Vander Waal bonding better than PBE. The error in individual data points is also visualized in Fig. 2(d-f) and (j-l) as relative percentage error. We find that the relative error can be up as large as 53%, as it is the case for FeS (P4/nmm) computed using PBE (materials-project-id mp-20311). Examples of other materials with high relative error in lattice constants are SnS (Cmcm, mp-1379), RbMnAs (P4/nmm, mp-20242), KMnP (P4/nmm, mp-20422), NiBi (P6 3 / mmc, mp-22861), and SnCl2 (P2 1 /c, mp-29179). Although, optB88 can reduce the error in lattice-constants for the predicted layered materials, there still is plenty of room for improvement, as neither of these two functionals predicts lattice constants with 100% accuracy for all the materials. After identification of possible 2D materials Scientific RepoRts | 7: 5179 | DOI:10.1038/s41598-017-05402-0 with the lattice-constant criterion, we validate our findings by computing their exfoliation energy. The exfoliation energy for 2D materials is computed as equation 3: Here, E 1L and E bulk are the energies of the single layer and 3D bulk materials and N 1L and N bulk are the number of the atoms in the single layer and bulk systems respectively. Although PBE results are available for the 3D bulk structures in the Materials Project database, we carried out new bulk calculations with the optB88 functional to determine the exfoliation energy, for consistency between bulk and single layer energetics treatment. Both bulk and single layer data and metadata are available publicly through an easy-to-use web interface, which is discussed in the next section. The conventional computational criterion [49][50][51][52] to predict feasibility of exfoliation suggests that exfoliation energy for a feasible 2D materials should be less than 200 meV/atom. As for equation (2), both bulk and monolayer energies must be computed, and at present, we have carried out 1012 out of 1356 bulk calculations, and 430 single layer calculations. All the remaining calculations are currently in progress and the database is continuously expanding. It is to be noted that the computational cost of single layer material calculations can be much higher than bulk due to inclusion of vacuum padding in the simulation box. We have 371 materials in our database which have both single and bulk calculations completed, and out of these 330 (88.9%) have satisfied the exfoliation energy criterion (exfoliation energy below 200 meV) as shown in Fig. 1.

Figure 2.
Comparison of lattice constants from Materials-project PBE calculations (a-f) and our optB88 functional calculations (g-l) with experimental ICSD data. Lattice constant a is examined in (a) PBE data and (g) optB88 data, while lattice constants b and c are given for PBE data in (b) and (c), respectively, and in (h) and (i) for optB88 data, respectively. In all cases the corresponding mean absolute error (using Eq. 2) is also reported, quantifying the improvement in lattice constants using optB88 over PBE. The data is plotted using dots, while the solid line is given just to indicate what the results should be for a perfect agreement between experimental and computational data. Similarly, relative percent errors with respect to the ICSD data are given in (d-f) for PBE a, b, and c lattice constants, respectively, and in (j-l) for optB88 equivalent data.
Even more interestingly, most of the predicted layered structures were found to have exfoliation energies in the 60-100 meV range. This result clearly suggests that our simple criterion to distinguish vdW bonded structures based on the relative difference δ in lattice constants is a reasonable criterion, at least for the initial screening. As an added advantage, this criterion is very easy to apply and anybody can use it. To facilitate experimentalists requiring options for easy-to-fabricate 2D materials, we have grouped materials with respect to their energy of exfoliation in our database as shown in Fig. 3 and Table 1. Specifically, we have a group for exfoliation energies below 0-20 meV/atom, one for energies between 20 and 40 meV/atom, and so on. Along with the chemical composition, we provide the space group information for bulk and single-layer material in parenthesis, which are hyperlinked to our database webpage. Clicking on these space group navigate to database webpage directly. While there has not been experimental work on all of these 2D materials yet, many of them have been very recently predicted using DFT by Materialsweb database 31 also. Materials such as TiIN, Sc 2 NCl 2 , OsCl 2 O, MnCl 2 , MgBr 2 are found to be common between ours and Materialsweb database. A detailed comparison of these two databases needs to be made for a comprehensive outlook for 2D materials. Our criterion, however, suggests a few more (such as BCl 3 , Te 2 Br, AlI 3 , TiI 3 and so on) materials that have not been identified yet. Moreover, new emerging 2D materials like BiTeI have already been proven candidate of 2D materials with very interesting physical phenomenon such as giant Rashba effect 53,54 , hence other novel 2D materials could also be subjected to experimental fabrication. Although, it is not recommended to invest lot of effort in experimental fabrication of layered materials with E f greater than 200 meV/atom, it is to be noted that large formation energy generally indicates that the 2D material could be formed only if a suitable substrate is found as for silicene and GeO 49,55 in which the exfoliation energy can take the value up to 600 meV/atom. Therefore, a layered material with high exfoliation energy cannot be completely discarded as candidate 2D materials.
Applicability and validation of the database. Our database provides the users with immediate access to several properties that are key in describing a 2D material and determining its technological applicability. All such properties have been computed consistently, to enable meaningful comparison between different materials. The properties that we focused on in this database are energetics, structural, elastic and electronic properties for both monolayer and bulk materials. While some of the calculated properties can deviate from experimental data, a consistent database like ours will enable at least a meaningful qualitative comparison of materials. We provide the applicability and limitation of our database as discussed below.
X-ray diffraction pattern of bulk and monolayer materials. XRD patterns act as key signatures for materials. Using DFT optimized crystal structure, XRD patterns for bulk and monolayer are obtained in a consistent way. As an example, XRD patterns for bulk and single layer 2H-MoTe 2 are shown in Fig. 4. To validate our DFT-based XRD, we also carried out experiments (Fig. 4a) to benchmark the bulk 2H-MoTe 2 structure. As shown in Fig. 4b, there is great agreement between the DFT and the experimental XRD peak positions, while the intensities do not always match in the case of bulk 2H-MoTe 2 , which can be attributed to the experimental sample preparation techniques. For single layer MoTe 2 , the XRD pattern is shown in Fig. 4c. Peak positions are quite different between single layer and bulk, which could be explained in terms of the structure modifications occurring as we switch from bulk to single layer MoTe 2 . Strong diffraction lines for single layer materials are a characteristic feature of 2D materials, as they are also observed in the case of MoS 2 56 . Although experimental XRD patterns are commonly available for bulk materials, obtaining experimental single layer XRD is much more challenging. This adds to the value of our database, as we are providing such a quantity and it could be used by experimentalist as a reference to compare and analyze their findings. Once again, it is important to point out that vdW functionals such as optB88, or conventional LDA functional produce XRD patterns for bulk materials that well match experimental ones, while semi-local functionals like PBE give erroneous XRD patterns, as they are incorrect in  Table 1.
Scientific RepoRts | 7: 5179 | DOI:10.1038/s41598-017-05402-0   determining lattice constants to begin with. We examined at least 20 different material cases using three different exchange correlation functionals, to guarantee we used the most appropriate one for our DFT calculations. Results clearly show that the peak positions in computational XRD using optB88 are very close to experimental ones, while their intensities may not match. A comparison of XRD experimental results to XRD DFT spectra obtained using 3 different exchange correlation functionals is given in the supplementary information (Fig. S2), together with 5 more examples of good match between optB88 spectra and experimental ones (Fig. S3).
Band-structure of bulk and monolayer materials. We compute band structures, density of states and plane-average electrostatic potential to characterize the electronic properties of bulk and single layer materials. An example for 2H bulk and single layer MoTe 2 is shown in Fig. 5. As expected, we notice a change in the nature of the bandgap, from indirect to direct, and a change in the bandgap value going from bulk to single layer. In the figure the red dots indicate the conduction band minima and the green dots show the valence band maxima. In addition to results computed using the optB88 functional, our database also contains HSE06 functional-based band structures for few materials. HSE06 functionals are important to obtain correct description of bandgaps of materials because all the local and semilocal DFT such as LDA and GGA are prone to underestimate the bandgaps. However, HSE06 are computationally very expensive, hence only few HSE06 for bulk and single layer materials are calculated and the database for HSE06 calculations is still populating. Nevertheless, in absence of better results, our optB88 based bandgaps could still be used to predict qualitative trends. To validate our high throughput electronic property data for single-layer and bulk materials, we compared our results to previous experimental and DFT data for quite a few 2D materials, as shown in Table 2. We notice that, although the trends in the bandgaps look consistent, their exact values are less satisfactory. We also observe that optB88 generally underestimates, and HSE06 somewhat overestimates the bandgaps. A detailed list of bandgaps using optB88 functional is given in the supplementary information, to enable an appropriate selection of materials (Table S6). DFT generally underestimates band-gaps of materials, hence if a material has a zero bandgap, it may or not be metallic compared to experimental data. Thankfully if DFT predicts non-zero bandgap then the material must be semiconducting or insulator in nature. We find 61% optB88-based bandgaps of single layer materials are greater than 0 or in other words they are non-metallic. This suggests that most of the layered materials in our database could be semiconducting or insulating in nature. HSE06 based bandgaps for most of the materials are underway, which can give better prediction of bandgaps. Interestingly, the metallic phase materials can also be of technological importance due to emergence of topological states, which are currently an area of intense research 57 . A detailed analysis of whether a material has a topological state or not is beyond the scope of present paper. Lastly, we also computed the work-function of the single-layer materials using the z-direction averaged potentials, as shown in Table 2. These values can be used to characterize the n/p nature of the 2D heterostructures. For a heterostructure made with two single layer materials, if the work-function for the first material is higher than the second one, then the first one should act as an n-type and second one as a p-type semiconductor. The plane average potential and hence work-functions are calculated for all the 2D materials. As evident from Table 2  nature instead of volume nature used in calculation. It is to be noted that during the calculation the vacuum size is taken as z-dimension of the simulation box so to make the elastic constant z-dimension free, we rescale the elastic constants with the z-simulation box size giving the elastic constant in the unit of N/m 58,59 . Moreover, the bulk and the monolayer elastic properties are not comparable but all the monolayer data could be compared for their elastic properties among themselves. A table for comparison of elastic constant for some reported elastic constant in the literature in Table 2. The DFT data for elastic constants are found to be within 10% tolerance of the previously reported values. It correctly reflects our data are comparable to previously reported data. Similar to the  structural data, elastic constants for 2D materials are sensitive to the selection of functional. A comparison of the elastic constants using different functionals are given in supplementary information (Table S7). We detect that elastic constants could deviate upto 27% with PBE functional are used such as for bulk MoS 2 (176 GPa 60 compared to 238 GPa from experiment 61 ), however with optB88 (214 GPa) the error is generally within 10%.
Another key property that can help experimentalists characterize their data is the phonon description. These phonon modes act as the signatures of materials and can be easily compared with experimental data. One of the key advantage characterizing the 2D or any other materials with phonon is that the experimental spectra is less affected by impurities in the materials than X-ray photoelectron spectroscopy (XPS) and other experimental techniques. However, we can only compute Gamma-point phonons using the conventional cell for each material for computational cost reasons. Highly negative phonon modes generally imply structural instability or phase-transition in the material. It is to be emphasized that elastic constants and phonon modes generally require a careful investigation of finite size effects on properties, but that is beyond the scope of this work right now. However, our phonon results can be used to screen layered materials with highly negative phonon modes. As an example, we compare the phonon modes and their representation for 1T'-MoTe 2 in Table 3. The experimental results are obtained with polarization-resolved Raman experiment 45 . We notice that our 1 × 1 × 1 conventional cell phonon results are comparable to experiments. This could be attributed due the fact that conventional cell for MoTe 2 was large enough (more than 1.2 nm) enabling proper description of atomic vibrations. However, this may not be case for all the materials, because their conventional cells could be relatively small. The phonon density of states of the system is given in the supplementary information (Fig. S4). To ensure the effect of cell size, we also compute phonons for 2 × 2 × 2 conventional cell. We notice the 2 × 2 × 2 and 1 × 1 × 1 data are comparable and close to experiments. Therefore, consistency in different sizes ensures the validity of extension of the calculation into phonon description. In addition to providing the phonon modes, we also provide the phonon representation of the phonons using symmetry of the eigenvectors of phonon mode for most of the materials with the help of Phonopy package 62 . However, we manually implement all the crystal point groups in the Phonopy package that were missing before along with their Infrared and Raman activity as available at Bilbao server 63 . The phonon representation modes assignment is challenging task in experiments. As mentioned above, the experimental mode assignment is carried with angular Raman spectroscopy, which are much expensive. In this way we can assign the phonon modes with proper representation for unknown modes in literature such as for 1T'-MoTe 2 64 . Comparison of the experimental and our DFT phonon modes show that our phonon modes could deviate up to 8% from experiment. However, our phonon representations have great agreement with experiments. Additionally, we compared the elastic constants obtained from 1 × 1 × 1 and 2 × 2 × 2 supercell for BN, BiClO and AuI as additional test cases for validation in the supplementary section (Table S8). We didn't notice significant change in elastic constants due to size.
Web-interface. The database is accessible through an easy to use web interface. A snapshot of the website is given in Fig. 6.
Here, for example clicking on Mo and Te and clicking "Search" button shows the options of calculations. In the first column of the results is the JARVIS-ID which leads to as detailed web-page of the calculation, second column is the formula of the structure, third column in the functional used in the calculation, fourth column is the calculation type for example 3D bulk or 2D single layer, fifth column is bandgap and its nature such as direct or indirect using LDA/PBE or optB88 functional, sixth column is the same information but for HSE06 functional, wherever available, seventh column is the Voigt bulk modulus, eighth column is the shear modulus, ninth column is the space group of the structure used in calculation, tenth column is the energy per atom of the system which can be post processed to calculate heat of formation. The energy per atom value difference between bulk and its corresponding single-layer material with consistent functional can be used to calculate the exfoliation energy. Clicking on any of the JARVIS-ID navigates to a new webpage containing detailed information of the calculation. The input files used in the calculation are provided as zipped files to ensure reproducibility of the data. Visualization of the structure is provided with the help of JSMOL 65 . The atomic structure could be rotated, and other valuable information such bond-distance, angles could be obtained interactively by simple click operations as facilitated by JSMOL. Next, convergence of the energy-cutoff and K-points are shown to ensure proper DFT optimization of the structure. Then XRD and RDF analysis of the structure is given. After that electronic structure information, such as density of states and band structure are given with total, orbital and element dos buttons to facilitate the visualization of orbital and elemental character of electronic states. After this, z-plane averaged potential is provided which can be used to calculate work-function for layered materials. Interestingly, the work-function of the 2D single layer material is close to the bulk z-plane averaged potential. This can be due to the interaction in the vdW bonding-direction is small as expected. Moreover, the elastic constant and allied phonon data are also provided for the cases where calculation has been completed. We also provide online visualization of the phonon modes using JSMOL. We provide ICSD-IDs 29 for bulk materials at respective webpages. These ICSD-IDs can be used to find out how the materials were experimentally fabricated. All the above mentioned easy to use web-tools are provided for free on NIST's robust CTCMS web-server http://www.ctcms.nist.gov/~knc6/JVASP.html.

Discussions
According to functional comparison investigation of Blaha et al. 36,48 , PBE generally shows 2.5% maximum error tolerance in lattice constant, thus here we assume any error that is greater than this value may be due to the extra degree of freedom such as vdW bonding. We set quite a conservative screening criterion of 5%, but it can be shown that even values as low as 3.5% can be used, such as in case of GeS, for example. If we consider transition metal chalcogenides and oxides as already known 2D-materials, hundreds of new 2D materials were found here, and the database is still growing. We expect to complete all the exfoliation energy calculations within 5 months, while computing the energetics, electronic and elastic properties will take longer, possibly up to 1 year. Our database clearly shows the 2D materials could not only be transition metal chalcogenides but they could be halides, pnictides and combination of them. The greatest advantage of the initial screening criterion is that it didn't require any additional calculation, since ICSD and Materials Project (MP) already contains all the necessary data already publicly available. However, for energy of exfoliation calculation, we had to compute the energetics of both the bulk and single layer material using optB88 functions, which were not commonly available. Moreover, we concentrated first on the structures with available ICSD IDs because it almost guarantees the feasibility of the experimental material fabrication. Please note that the "discovery" aspect of this work is related to the identification of possible 2D materials from known bulk compounds, and that is accomplished by our strategy, even if no new 3D materials can be identified using our approach. Although the data is provided for MP only, it is being extended to other databases such as AFLOW and OQMD, and a concise database, encompassing all these above-mentioned databases for predicting 2D materials, will be available soon. The easy to use web-interface with HTML features, and publicly available input data will facilitate the data reproducibility and enhances the machine learning capabilities; a topic which will be discussed in another paper. Our database can give experimentalist a way to compare for the collected XRD and phonon data for initial experimental comparison. HSE06 band structures have recently become very popular in DFT community for their accurate results and relatively cheap cost compared to many many-body calculations. We plan to calculate the HSE06 band-structures for all the materials  in our database. We have not analyzed the phonon data for many materials, however it is commonly known that finite size effect, as well as the method for phonon calculation such finite-difference and density functional perturbation theory (DFPT), play critical role in designating phonons for Raman measurements. Later we plan to carry out supercell phonon calculation for all the structures with both finite difference and DFPT to enhance quality of the phonon data. We also plan to compute thermoelectric, piezoelectric, defect energetics, interface energetics and multi-layer properties of materials with and without strain in our database. In a nutshell, C 33 being very low, high work function and bandgaps in the range 0.5 to 3.5 eV were identified to be characteristics of 2D materials. Once again, it is worth to mention that bandgaps are generally underestimated in DFT (unless HSE06, GW or other high-level calculations are carried out), however, if a material has non-zero bandgap in DFT calculation then it must be semiconducting or insulating. Our publicly available data will allow user to select material with required characteristics at absolutely no cost. A characteristic feature of this database is a presence of both the bulk and single layer data which enables easy comparison of properties as material changes from infinite to finite size. Calculation of data for the samples comprised of multiple, but small number of layers (i.e. 2, 3, 4) is more challenging, but it will also be available in the database in the future. Despite this database still being in its initial developmental stage, it serves as a comprehensive detailed source for finding properties of the single layer 2D-layered materials in conjunction with their corresponding 3D-bulk state. Our property calculations here are our first steps in characterizing layered materials in a systematic way, much more involved property analysis is subjected to future work. After accumulating a enough data in this project, we can employ machine learning methods to predict behavior of 3D bulk as well as 2D materials in a single framework.

Methods
Density functional theory. All the calculations are done with Vienna Ab initio Simulation Package (VASP) software 66 . Please note commercial software is identified to specify procedures. Such identification does not imply recommendation by the National Institute of Standards and Technology. 2D bulk material structures are taken from MP, OQMD, AFLOW and ICSD or from latest research papers. These DFT databases were made with crystal structure data obtained from ICSD database. We use the ICSD lattice constants and DFT lattice constants with PBE as a guide for identifying 2D materials as shown in Fig. 7. We use 5% expansion criterion between PBE and ICSD to identify 2D materials. The geometrical input structure for input were obtained from application program interface (API) of the DFT databases mainly. The calculation set-up and handling of data is primarily tackled with Pymatgen 41 and atomic simulation environment (ASE) 67 . The 1 L materials are made manually by separating single layers and adding vacuum padding in z-direction. At present, non-polar layered materials are computed only. We use LDA, PBE and optB88 functionals for many 2D materials to facilitate comparison of properties using different functionals. The K-point and Energy cut-off convergence is achieved with 0.001 eV tolerance with a series of DFT calculations. This convergence calculation is performed automatically. We keep increasing the Energy cut-off and K-point density of mesh until the difference between the previous and present pointy in the plot is within 0.001 eV tolerance. After achieving the tolerance criterion, we test five more points to ensure if there is some other minimum doesn't exist. If another minimum is found, the convergence calculation is started again from the last point automatically. The K-point is chosen in an automatic mesh way as provided in VASP, which uses the reciprocal lattice vectors of the system at an increment of 5, while the energy cut-off is incremented with a step size of 50 eV. The force and energy convergence for DFT self-consistent calculations are 10 −6 eV and 0.001 eV/Å respectively. It is to be noted that other DFT databases do not enforce force-convergence criterion which can result in incorrect description of elastic and phonon properties. After the structure optimization, we calculate X-ray diffraction (XRD), radial distribution function (RDF) plots, electronic band-structure, optical properties, and elastic constants in different sets of calculations in automated way. The elastic constants are calculated on conventional cell using space group information of the structures, while phonon-bandstructures are calculated with Phonopy package 62 . We have carried out phonon convergence for finite size effects for few systems only. The phonons obtained from conventional cell are also provided with visualization tools. We used the conventional K-points for electronic and phonon bandstructure, computational XRD as implemented in Pymatgen 41 and AFLOW 35 . The data presented in html format as well as available in JSON format. A flow-chart for the process mentioned above is given below in Fig. 8.
X-ray diffraction. MoTe 2 powders and crystals were finely ground in hexane using an agate mortar; θ-2θ scans were collected at ambient temperature on a powder X-ray diffractometer equipped with general area detector using CuKα radiation (λ = 0.15418 nm). MDI-JADE 6.5 software was used to analyze the experimental data (Jade 6.5, Materials Data, Inc. Livermore, CA., 2015).

Raman-measurement.
The polarization-resolved Raman measurements were acquired by sending linearly polarized light at 532 nm through a motorized achromatic halfwave plate to control the excitation polarization. The same microscope objective (0.75 NA) was used for excitation and collection. The resulting Raman signal passed through a longpass filter followed by a motorized linear polarizer. Measurements were acquired by rotating the excitation polarization and analyzer together in either co-polarized or cross-polarized configurations. The spectra were acquired using a 500 mm focal length spectrometer with a 2400 lines/mm grating and CCD. To avoid damaging the samples an excitation power of 0.35 mW at the microscope objective was used with a 2 min spectral acquisition. The spectra were calibrated using an Ar lamp, and the laser frequency was measured with a wavemeter.

Conclusions
In this work, we have created a new publicly available database specially dedicated to 2D materials. Our database provides detailed analysis of all the bulk and layered materials with appropriate vdW functional. Energy of exfoliation, bandgap, work-function and elastic constant values are presented and can be used to easily and consistently compare the materials of interest. All the database and the input files are provided in public domain for free to facilitate reproducibility of the DFT results. We also mention all the approximations used in the calculation on the webpage to enable the reader understand the limitations of the results. With the help of the data from our Figure 7. Illustration about initial screening of 2D materials is shown. DFT databases used ICSD crystallographic database to characterize bulk materials with PBE functional. PBE functional is known to overestimate lattice constant of materials which can be 2D. Based on API, materials with reasonably high PBE to Experimental lattice constant difference were screened from thousands of materials available in the DFT databases. Layered counterparts of these materials were made manually. Both bulk and layered materials were then subjected to a series of DFT calculations for their property calculation.
repository, we have demonstrated the effectiveness of a simple criterion to identify 2D materials based on lattice constants that can be used in predicting new layered materials. We have predicted hundreds of novel materials which satisfy both the lattice difference as well as exfoliation energy criterion to deem them as 2D materials, and the database is still expanding. Experimental synthesis and characterization of these 2D materials could reveal a vast range of material properties and phenomenon. All the bulk material calculation webpages are hyperlinked to Materials-Project and AFLOW databases to enable comparison of properties in different databases. We believe our new predicted 2D materials will be advantageous for many material science problems.