Edge reconstruction in armchair phosphorene nanoribbons revealed by discontinuous Galerkin density functional theory

Edge reconstructions of black phosphorene: a global search

Despite reports of possible edge reconstructions of black phosphorene, the underlying mechanisms that determine the atomic configurations and appearance of black phosphorene edges have not been elucidated to date. In this study, the particle swarm optimization (PSO) algorithm is used to perform a global search of black phosphorene edge structures. In addition to the most stable edges, three databases of the typical black phosphorene zigzag edge, armchair edge, and skewed diagonal edge are constructed for the first time. The local phosphorus concentration plays an important role in determining the edge atomic configurations and the appearance of an edge. Variations in the local phosphorus concentration result in the rearrangement of sp3-hybrid bonds or the formation of double bonds that balance the dangling bonds at the edges and stabilize the black phosphorene edges. The black phosphorene edge databases provide a useful reference for edge studies of other 2D materials with puckered honeycomb structures.

High performance computing of DGDFT for tens of thousands of atoms using millions of cores on Sunway TaihuLight

High performance computing (HPC) is a powerful tool to accelerate the Kohn-Sham density functional theory (KS-DFT) calculations on modern heterogeneous supercomputers. Here, we describe a massively parallel implementation of discontinuous Galerkin density functional theory (DGDFT) method on the Sunway TaihuLight supercomputer. The DGDFT method uses the adaptive local basis (ALB) functions generated on-the-fly during the self-consistent field (SCF) iteration to solve the KS equations with high precision comparable to plane-wave basis set. In particular, the DGDFT method adopts a two-level parallelization strategy that deals with various types of data distribution, task scheduling, and data communication schemes, and combines with the master-slave multi-thread heterogeneous parallelism of SW26010 processor, resulting in large-scale HPC KS-DFT calculations on the Sunway TaihuLight supercomputer. We show that the DGDFT method can scale up to 8,519,680 processing cores (131,072 core groups) on the Sunway TaihuLight supercomputer for studying the electronic structures of two-dimensional (2D) metallic graphene systems that contain tens of thousands of carbon atoms.

Parallel Implementation of Large-Scale Linear Scaling Density Functional Theory Calculations With Numerical Atomic Orbitals in HONPAS

Linear-scaling density functional theory (DFT) is an efficient method to describe the electronic structures of molecules, semiconductors, and insulators to avoid the high cubic-scaling cost in conventional DFT calculations. Here, we present a parallel implementation of linear-scaling density matrix trace correcting (TC) purification algorithm to solve the Kohn-Sham (KS) equations with the numerical atomic orbitals in the HONPAS package. Such a linear-scaling density matrix purification algorithm is based on the Kohn's nearsightedness principle, resulting in a sparse Hamiltonian matrix with localized basis sets in the DFT calculations. Therefore, sparse matrix multiplication is the most time-consuming step in the density matrix purification algorithm for linear-scaling DFT calculations. We propose to use the MPI_Allgather function for parallel programming to deal with the sparse matrix multiplication within the compressed sparse row (CSR) format, which can scale up to hundreds of processing cores on modern heterogeneous supercomputers. We demonstrate the computational accuracy and efficiency of this parallel density matrix purification algorithm by performing large-scale DFT calculations on boron nitrogen nanotubes containing tens of thousands of atoms.

Two-Level Chebyshev Filter Based Complementary Subspace Method: Pushing the Envelope of Large-Scale Electronic Structure Calculations

We describe a novel iterative strategy for Kohn-Sham density functional theory calculations aimed at large systems (>1,000 electrons), applicable to metals and insulators alike. In lieu of explicit diagonalization of the Kohn-Sham Hamiltonian on every self-consistent field (SCF) iteration, we employ a two-level Chebyshev polynomial filter based complementary subspace strategy to (1) compute a set of vectors that span the occupied subspace of the Hamiltonian; (2) reduce subspace diagonalization to just partially occupied states; and (3) obtain those states in an efficient, scalable manner via an inner Chebyshev filter iteration. By reducing the necessary computation to just partially occupied states and obtaining these through an inner Chebyshev iteration, our approach reduces the cost of large metallic calculations significantly, while eliminating subspace diagonalization for insulating systems altogether. We describe the implementation of the method within the framework of the discontinuous Galerkin (DG) electronic structure method and show that this results in a computational scheme that can effectively tackle bulk and nano systems containing tens of thousands of electrons, with chemical accuracy, within a few minutes or less of wall clock time per SCF iteration on large-scale computing platforms. We anticipate that our method will be instrumental in pushing the envelope of large-scale ab initio molecular dynamics. As a demonstration of this, we simulate a bulk silicon system containing 8,000 atoms at finite temperature, and obtain an average SCF step wall time of 51 s on 34,560 processors; thus allowing us to carry out 1.0 ps of ab initio molecular dynamics in approximately 28 h (of wall time).

Interpolative Separable Density Fitting through Centroidal Voronoi Tessellation with Applications to Hybrid Functional Electronic Structure Calculations

The recently developed interpolative separable density fitting (ISDF) decomposition is a powerful way for compressing the redundant information in the set of orbital pairs and has been used to accelerate quantum chemistry calculations in a number of contexts. The key ingredient of the ISDF decomposition is to select a set of nonuniform grid points, so that the values of the orbital pairs evaluated at such grid points can be used to accurately interpolate those evaluated at all grid points. The set of nonuniform grid points, called the interpolation points, can be automatically selected by a QR factorization with column pivoting (QRCP) procedure. This is the computationally most expensive step in the construction of the ISDF decomposition. In this work, we propose a new approach to find the interpolation points based on the centroidal Voronoi tessellation (CVT) method, which offers a much less expensive alternative to the QRCP procedure when ISDF is used in the context of hybrid functional electronic structure calculations. The CVT method only uses information from the electron density and can be efficiently implemented using a K-Means algorithm. We find that this new method achieves comparable accuracy to the ISDF-QRCP method, at a cost that is negligible in the overall hybrid functional calculations. For instance, for a system containing 1000 silicon atoms simulated using the HSE06 hybrid functional on 2000 computational cores, the cost of the QRCP-based method for finding the interpolation points is 38.1 s, while the CVT procedure only takes 0.7 s. We also find that the ISDF-CVT method enhances the smoothness of the potential energy surface in the context of ab initio molecular dynamics (AIMD) simulations with hybrid functionals.

Interpolative Separable Density Fitting Decomposition for Accelerating Hybrid Density Functional Calculations with Applications to Defects in Silicon

We present a new efficient way to perform hybrid density functional theory (DFT)-based electronic structure calculations. The new method uses an interpolative separable density fitting (ISDF) procedure to construct a set of numerical auxiliary basis vectors and a compact approximation of the matrix consisting of products of occupied orbitals represented in a large basis set such as the planewave basis. Such an approximation allows us to reduce the number of Poisson solves from [Formula: see text] to [Formula: see text] when we apply the exchange operator to occupied orbitals in an iterative method for solving the Kohn-Sham equations, where N_e is the number of electrons in the system to be studied. We show that the ISDF procedure can be carried out in [Formula: see text] operations, with a much smaller preconstant compared to methods used in existing approaches. When combined with the recently developed adaptively compressed exchange (ACE) operator formalism, which reduces the number of times the exchange operator needs to be updated, the resulting ACE-ISDF method significantly reduces the computational cost associated with the exchange operator by nearly 2 orders of magnitude compared to existing approaches for a large silicon system with 1000 atoms. We demonstrate that the ACE-ISDF method can produce accurate energies and forces for insulating and metallic systems and that it is possible to obtain converged hybrid functional calculation results for a 1000-atom bulk silicon within 10 min on 2000 computational cores. We also show that ACE-ISDF can scale to 8192 computational cores for a 4096-atom bulk silicon system. We use the ACE-ISDF method to geometrically optimize a 1000-atom silicon system with a vacancy defect using the HSE06 functional and computes its electronic structure. We find that that the computed energy gap from the HSE06 functional is much closer to the experimental value compared to that produced by semilocal functionals in the DFT calculations.

Projected Commutator DIIS Method for Accelerating Hybrid Functional Electronic Structure Calculations

The commutator direct inversion of the iterative subspace (commutator DIIS or C-DIIS) method developed by Pulay is an efficient and the most widely used scheme in quantum chemistry to accelerate the convergence of self-consistent field (SCF) iterations in Hartree-Fock theory and Kohn-Sham density functional theory. The C-DIIS method requires the explicit storage of the density matrix, the Fock matrix, and the commutator matrix. Hence, the method can only be used for systems with a relatively small basis set, such as the Gaussian basis set. We develop a new method that enables the C-DIIS method to be efficiently employed in electronic structure calculations with a large basis set such as planewaves for the first time. The key ingredient is the projection of both the density matrix and the commutator matrix to an auxiliary matrix called the gauge-fixing matrix. The resulting projected commutator-DIIS method (PC-DIIS) only operates on matrices of the same dimension as that consists of Kohn-Sham orbitals. The cost of the method is comparable to that of standard charge mixing schemes used in large basis set calculations. The PC-DIIS method is gauge-invariant, which guarantees that its performance is invariant with respect to any unitary transformation of the Kohn-Sham orbitals. We demonstrate that the PC-DIIS method can be viewed as an extension of an iterative eigensolver for nonlinear problems. We use the PC-DIIS method for accelerating Kohn-Sham density functional theory calculations with hybrid exchange-correlation functionals, and demonstrate its superior performance compared to the commonly used nested two-level SCF iteration procedure. Furthermore, we demonstrate that in the context of ab initio molecular dynamics (MD) simulation with hybrid functionals one can extrapolate the gauge-fixing matrix to achieve the goal of extrapolating the entire density matrix implicitly along the MD trajectory. Numerical results indicate that the new method significantly reduces the number of SCF iterations per MD step, compared to the commonly used strategy of extrapolating the electron density.

Modulation of electronic transport properties in armchair phosphorene nanoribbons by doping and edge passivation

The electronic structures and transport properties of group IV atoms (C, Si, Ge)-doped armchair phosphorene nanoribbons (APNRs) are investigated using first-principles calculations, considering different edge passivation. The results show that the C, Si, Ge dopants can induce the transition occur from semiconductor to metal in the APNRs. The negative differential resistance (NDR) behavior in the doped APNR system is robust with respect to the doping concentration and edge passivation type. However, their current peak positions and peak-to-valley ratio (PVR) values are correlated with doping concentration and edge passivation type. In particular, for the C, Si-doped APNRs, the low bias NDR behavior with the PVR (10⁵–10⁸) can be observed when doping concentration is low in the APNRs with the F and H edge passivation. These results may play an important role for the fabrication of future low power consumption nano-electronic devices.

Quantum confinement in black phosphorus-based nanostructures

The modification of an idealized infinite bulk system by dimensional reduction or structural distortion results in quantum confinement effects (QCEs). For example, dimensional reduction of a black phosphorus structure leads to the realization of few-layer systems, creation of edges and surfaces, nanoribbons, quantum dots, and antidot lattices while structural distortion involves simple bending (including nanotubes) and rippling. Black phosphorus ('phosphorene' in the single-layer limit) has been of recent interest due to its relatively large charge carrier mobility and moderate semiconducting band gap, which remains direct irrespective of the number of layers. In this review the state-of-the-art properties of black phosphorus in its dimensionally reduced and structurally distorted forms are discussed, with emphasis on how quantum confinement impacts the material's properties.

Adaptively Compressed Exchange Operator for Large-Scale Hybrid Density Functional Calculations with Applications to the Adsorption of Water on Silicene

Density functional theory (DFT) calculations using hybrid exchange-correlation functionals have been shown to provide an accurate description of the electronic structures of nanosystems. However, such calculations are often limited to small system sizes due to the high computational cost associated with the construction and application of the Hartree-Fock (HF) exchange operator. In this paper, we demonstrate that the recently developed adaptively compressed exchange (ACE) operator formulation [J. Chem. Theory Comput. 2016, 12, 2242-2249] can enable hybrid functional DFT calculations for nanosystems with thousands of atoms. The cost of constructing the ACE operator is the same as that of applying the exchange operator to the occupied orbitals once, while the cost of applying the Hamiltonian operator with a hybrid functional (after construction of the ACE operator) is only marginally higher than that associated with applying a Hamiltonian constructed from local and semilocal exchange-correlation functionals. Therefore, this new development significantly lowers the computational barrier for using hybrid functionals in large-scale DFT calculations. We demonstrate that a parallel planewave implementation of this method can be used to compute the ground-state electronic structure of a 1000-atom bulk silicon system in less than 30 wall clock minutes and that this method scales beyond 8000 computational cores for a bulk silicon system containing about 4000 atoms. The efficiency of the present methodology in treating large systems enables us to investigate adsorption properties of water molecules on Ag-supported two-dimensional silicene. Our computational results show that water monomer, dimer, and trimer configurations exhibit distinct adsorption behaviors on silicene. In particular, the presence of additional water molecules in the dimer and trimer configurations induces a transition from physisorption to chemisorption, followed by dissociation on Ag-supported silicene. This is caused by the enhanced effect of hydrogen bonds on charge transfer and proton transfer processes. Such a hydrogen bond autocatalytic effect is expected to have broad applications for silicene as an efficient surface catalyst for oxygen reduction reactions and water dissociation.

Nano-sized graphene flakes: insights from experimental synthesis and first principles calculations

In this study, we proposed a cost-effective method for preparing graphene nano-flakes (GNFs) derived from carbon nanotubes (CNTs) via three steps (pressing, homogenization and sonication exfoliation processes). Scanning electron microscopy (SEM), transmission electron microscopy (TEM), atomic force microscopy (AFM), laser scattering, as well as ultraviolet-visible and photoluminescence (PL) measurements were carried out. The results indicated that the size of as-synthesized GNFs was approximately 40-50 nm. Furthermore, we also used first principles calculations to understand the transformation from CNTs to GNFs from the viewpoints of the edge formation energies of GNFs in different shapes and sizes. The corresponding photoluminescence measurements of GNFs were carried out in this work.

Chebyshev polynomial filtered subspace iteration in the discontinuous Galerkin method for large-scale electronic structure calculations

The Discontinuous Galerkin (DG) electronic structure method employs an adaptive local basis (ALB) set to solve the Kohn-Sham equations of density functional theory in a discontinuous Galerkin framework. The adaptive local basis is generated on-the-fly to capture the local material physics and can systematically attain chemical accuracy with only a few tens of degrees of freedom per atom. A central issue for large-scale calculations, however, is the computation of the electron density (and subsequently, ground state properties) from the discretized Hamiltonian in an efficient and scalable manner. We show in this work how Chebyshev polynomial filtered subspace iteration (CheFSI) can be used to address this issue and push the envelope in large-scale materials simulations in a discontinuous Galerkin framework. We describe how the subspace filtering steps can be performed in an efficient and scalable manner using a two-dimensional parallelization scheme, thanks to the orthogonality of the DG basis set and block-sparse structure of the DG Hamiltonian matrix. The on-the-fly nature of the ALB functions requires additional care in carrying out the subspace iterations. We demonstrate the parallel scalability of the DG-CheFSI approach in calculations of large-scale two-dimensional graphene sheets and bulk three-dimensional lithium-ion electrolyte systems. Employing 55 296 computational cores, the time per self-consistent field iteration for a sample of the bulk 3D electrolyte containing 8586 atoms is 90 s, and the time for a graphene sheet containing 11 520 atoms is 75 s.

Ultrafast charge transfer dynamics pathways in two-dimensional MoS2–graphene heterostructures: a core-hole clock approach

Ultrafast electron delocalization pathways on the MoS₂/graphene heterostructure were elucidated.

New Ferroelectric Phase in Atomic-Thick Phosphorene Nanoribbons: Existence of in-Plane Electric Polarization

Ferroelectrics have many significant applications in electric devices, such as capacitor or random-access memory, tuning the efficiency of solar cell. Although atomic-thick ferroelectrics are the necessary components for high-density electric devices or nanoscale devices, the development of such materials still faces a big challenge because of the limitation of intrinsic mechanism. Here, we reported that in-plane atomic-thick ferroelectricity can be induced by vertical electric field in phosphorene nanoribbons (PNRs). Through symmetry arguments, we predicted that ferroelectric direction is perpendicular to the direction of external electric field and lies in the plane. Further confirmed by the comprehensive first-principles calculations, we showed that such ferroelectricity is induced by the electron-polarization, which is different from the structural distortion in traditional ferroelectrics and the recent experimental discovery of in-plane atomic-thick ferroelectrics (Science 2016, 353, 274). Moreover, we found that the value of electronic polarization in bilayer is much larger than that in monolayer. Our results show that electron-polarization ferroelectricity maybe the most promising candidate for atomic-thick ferroelectrics.

Efficient Computation of the Hartree-Fock Exchange in Real-Space with Projection Operators

We describe an efficient projection-based real-space implementation of the nonlocal single-determinant exchange operator. Through a matrix representation of the projected operator, we show that this scheme works equally well for both occupied and virtual states. Our scheme reaches a speedup of 2 orders of magnitude and has no significant loss of accuracy compared to an implementation of the full nonlocal single-determinant exchange operator. We find excellent agreement upon comparing Hartree-Fock eigenvalues, dipoles, and polarizabilities of selected molecules calculated using our method to values in the literature. To illustrate the efficiency of this scheme we perform calculations on systems with up to 240 carbon atoms.

Ferromagnetism controlled by electric field in tilted phosphorene nanoribbon

Study on phosphorene nanoribbon was mostly focused on zigzag and armchair structures and no ferromagnetic ground state was observed in these systems. Here, we investigated the magnetic property of tilted black phosphorene nanoribbons (TPNRs) affected by an external electric field. We also studied the edge passivation effect on the magnetism and thermal stability of the nanoribbons. The pure TPNR displayed an edge magnetic state, but it disappeared in the edge reconstructed TPNR due to the self-passivation. In addition, we found that the bare TPNR was mechanically unstable because an imaginary vibration mode was obtained. However, the imaginary vibration mode disappeared in the edge passivated TPNRs. No edge magnetism was observed in hydrogen and fluorine passivated TPRNs. In contrast, the oxygen passivated TPNR was more stable than the pure TPNR and the edge-to-edge antiferromagntic (AFM) ground state was obtained. We found that the magnetic ground state could be tuned by the electric field from antiferromagnetic (AFM) to ferromagnetic (FM) ground state. Interestingly, the oxygen passivated TPNR displayed a half-metallic state at a proper electric field in both FM and AFM states. This finding may provoke an intriguing issue for potential spintronics application using the phosphorene nanoribbons.

Edge-Modified Phosphorene Nanoflake Heterojunctions as Highly Efficient Solar Cells

We propose to use edge-modified phosphorene nanoflakes (PNFs) as donor and acceptor materials for heterojunction solar cells. By using density functional theory based calculations, we show that heterojunctions consisting of hydrogen- and fluorine-passivated PNFs have a number of desired optoelectronic properties that are suitable for use in a solar cell. We explain why these properties hold for these types of heterojunctions. Our calculations also predict that the maximum energy conversion efficiency of these type of heterojunctions, which can be easily fabricated, can be as high as 20%, making them extremely competitive with other types of two-dimensional heterojunctions.

Tuning anisotropic electronic transport properties of phosphorene via substitutional doping

Using first-principles calculations, we studied the anisotropic electronic transport properties of pristine and X-doped phosphorene (X = B, Al, Ga, C, Si, Ge, N, As, O, S, and Se atoms).