Merge: Update DFT+U & documentation for DFT+U and Wannier90 interface…

… from abacusmodeling/abacus-develop/develop (#1326) * Fix: SOC calculation in PW code would caused nan * Fix: setting of k-points in NSPIN=2 has conflict with Pkpoints in numerical_basis.cpp * Fix: nscf with hybrid functional 'hse' or 'pbe0' or 'scan0' * Fix: generating orbitals outputs with nspin=2 * Update the documents for DFT+U * Modify the help documents of wannier90 and berry phase Mainly modified Berry_phase.md and Wannier90.md in the docs/advanced/ folder. Adjusted the files in example/interface_wannier90. * Adjust the file in example/interface_wannier90. * Revise the error in the documents of DFT+U * Update docs/advanced/input_files/input-main.md Co-authored-by: Chun Cai <amoycaic@gmail.com> Co-authored-by: Qx80610702 <80610702@qq.com> Co-authored-by: jingan <jingan@mail.ustc.edu.cn> Co-authored-by: Chun Cai <amoycaic@gmail.com> Co-authored-by: Wenfei Li <38569667+wenfei-li@users.noreply.github.com>
deepmodeling · Sep 27, 2022 · f0fa34b · f0fa34b
1 parent 617bd3b
commit f0fa34b
Show file tree

Hide file tree

Showing 21 changed files with 3,993 additions and 103 deletions.
diff --git a/docs/advanced/elec_properties/Berry_phase.md b/docs/advanced/elec_properties/Berry_phase.md
@@ -15,7 +15,7 @@ cp KPT-nscf-c KPT
 mpirun -np 4 abacus
 ```
 
-In this example, we calculate the electric polarization along c axis for PbTiO3, and below are the INPUT file (nscf) and KPT file (nscf):
+In this example, we calculate the electric polarization along c axis for PbTiO~3~, and below are the INPUT file (nscf) and KPT file (nscf):
 
 ```
 INPUT_PARAMETERS
@@ -32,8 +32,8 @@ gdir          3 // calculate polarization along c axis
 ```
 
 Note: You need to turn off the symmetry when do Berry phase calculations. Currently, ABACUS support Berry phase calculation with nspin=1 and nspin=2. The Berry phase can be calculated in both pw and lcao bases.
-- berry_phase : 1, calculate berry phase; 0, no calculate berry phase.
-- gdir : 1, 2, 3, the lattice vector direction of the polarization you want to calculate.
+- [berry_phase](../input_files/input-main.md#berry_phase) : 1, calculate berry phase; 0, no calculate berry phase.
+- [gdir](../input_files/input-main.md#gdir) : 1, 2, 3, the lattice vector direction of the polarization you want to calculate.
 
 The KPT file need to be modified according to gdir in the INPUT file. Generally, you need denser k points along this direction. For example, in the following KPT file, 4 k-points are taken along the a and b axes, and 8 k-points are taken along the c-axis. You should check the convergence of the k points when calculating the polarization.
 
@@ -56,9 +56,6 @@ The results are shown as follows:
  |                                                                    |
  <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
 
-
-
-
  VALUES OF POLARIZATION
 
   The Ionic Phase:   -0.10600
@@ -73,4 +70,4 @@ The results are shown as follows:
  P =    0.8906925  (mod    2.1748536)  (   0.0000000,   0.0000000,   0.8906925) C/m^2
 ```
 
-The electric polarization <b>P</b> is multivalued, which modulo a quantum eR=<sub>Vcell</sub>. Note: the values in parentheses are the components of the <b>P</b> along the c axis in the x, y, z Cartesian coordinates when set gdir = 3 in INPUT file.
+The electric polarization **P** is multivalued, which modulo a quantum e**R**/V~cell~. Note: the values in parentheses are the components of the **P** along the c axis in the x, y, z Cartesian coordinates when set gdir = 3 in INPUT file.
diff --git a/docs/advanced/input_files/input-main.md b/docs/advanced/input_files/input-main.md
@@ -56,7 +56,7 @@
 
     [tddft](#tddft) | [td_scf_thr](#td_scf_thr) | [td_dt](#td_dt) | [td_force_dt](#td_force_dt) | [td_vext](#td_vext) | [td_vext_dire](#td_vext_dire) | [td_timescale](#td_timescale) | [td_vexttype](#td_vexttype) | [td_vextout](#td_vextout) | [td_dipoleout](#td_dipoleout) | [ocp](#ocp) | [ocp_set](#ocp_set)
 
-- [DFT+U correction](#DFTU-correction) (Under tests)
+- [DFT+*U* correction](#dft-u-correction) (Under development)
 
     [dft_plus_u](#dft_plus_u) | [orbital_corr](#orbital_corr) | [hubbard_u](#hubbard_u) | [hund_j](#hund_j) | [yukawa_potential](#yukawa_potential) | [omc](#omc)
 
@@ -1377,7 +1377,7 @@ temperature will fluctuate violently; if it is too small, the temperature will t
 
 [back to top](#full-list-of-input-keywords)
 
-### DFT+U correction
+### DFT+*U* correction
 
 This part of variables are used to control DFT+U correlated parameters
 

diff --git a/docs/advanced/interface/Wannier90.md b/docs/advanced/interface/Wannier90.md
@@ -1 +1,121 @@
-# Wannier90
+# Wannier90
+
+[Wannier90](http://www.wannier.org/) is a useful package to generating the maximally-localized Wannier functions (MLWFs), which can be used to compute advanced electronic properties. Some post-processing tools (such as WannierTools, etc.) will use MLWFs for further analysis and calculations. 
+
+Currently ABACUS provides an interface to Wannier90 package. The users are assumed to be familiar with the use of Wannier90. The ABACUS-Wannier90 interface is only suitable for nspin=1 or 2, not for nspin=4 or spin-orbit coupling (SOC). 
+
+To construct the MLWFs using the wave functions of ABACUS generally requires four steps. Here we use the diamond as an example which can be found in [examples/interface_wannier90/](https://github.com/abacusmodeling/abacus-develop/tree/develop/examples/interface_wannier90).
+
+1. Enter the `ABACUS_towannier90/` folder, prepare a Wannier90 input file `diamond.win`, which is the main input file for Wannier90. Then To generate `diamond.nnkp` file by running Wannier90,  which ABACUS will read later: 
+
+    ```
+    wannier90 -pp diamond.win
+    ```
+
+    The content of `diamond.win` is as follows:
+
+    ```
+    num_wann        =  4
+    num_iter        = 20
+    
+    wannier_plot=.true.
+    wannier_plot_supercell = 3
+    wvfn_formatted = .true.
+
+    begin atoms_frac
+    C   -0.12500  -0.1250    -0.125000
+    C    0.12500   0.1250     0.125000
+    end atoms_frac
+
+    begin projections
+    f=0.0,0.0,0.0:s
+    f=0.0,0.0,0.5:s
+    f=0.0,0.5,0.0:s
+    f=0.5,0.0,0.0:s
+    end projections
+
+    begin unit_cell_cart
+    -1.613990   0.000000   1.613990
+     0.000000   1.613990   1.613990
+    -1.613990   1.613990   0.000000
+    end unit_cell_cart
+
+    mp_grid : 4 4 4
+
+    begin kpoints
+    0.0000  0.0000     0.0000
+    0.0000  0.2500     0.0000
+    0.0000  0.5000     0.0000
+    0.0000  0.7500     0.0000
+    ...
+    end kpoints
+    ```
+
+2. Do a self-consistent calculation and get the converged charge density:
+
+    ```
+    cp INPUT-scf INPUT
+    cp KPT-scf KPT
+    mpirun -np 4 abacus
+    ```
+
+
+3. Do a non-self-consistent calculation: 
+
+    ```
+    cp INPUT-nscf INPUT
+    cp KPT-nscf KPT
+    mpirun -np 4 abacus
+    ```
+
+    below are the INPUT file (nscf):
+
+    ```
+    INPUT_PARAMETERS
+
+    ntype                   1
+    ecutwfc                 50
+    nbands                  4
+    calculation             nscf
+    scf_nmax                50
+    pw_diag_thr             1.0e-12
+    scf_thr                 1.0e-15
+    init_chg                file
+    symmetry                0
+    towannier90             1
+    nnkpfile                diamond.nnkp
+    ```
+
+    There are three interface-related parameters in the `INPUT` file:
+
+    - [towannier90](../input_files/input-main.md#towannier90): `1`, generate files for wannier90 code; `0`, do not generate.
+    - [nnkpfile](../input_files/input-main.md#nnkpfile) : the name of the file generated by running "wannier90 -pp ...".
+    - [wannier_spin](../input_files/input-main.md#wannier_spin): If you use nspin=2, `up`: calculate the Wannier functions for the spin up components ; `down`: calculate the Wannier functions spin down components.
+
+    Note: You need to turn off the symmetry during the entire nscf calculation.
+
+    To setup the `KPT` file according to the `diamond.win` file, which is similar to "begin kpoints ..." in the `diamond.win` file: 
+
+    ```
+    K_POINTS
+    64
+    Direct
+    0.0000 0.0000 0.0000 0.0156250
+    0.0000 0.2500 0.0000 0.0156250
+    0.0000 0.5000 0.0000 0.0156250
+    0.0000 0.7500 0.0000 0.0156250
+    ...
+    ```
+
+    After the nscf calculation, ABACUS will generate `diamond.amn`, `diamond.mmn`, `diamond.eig`, `UNK` files in the `OUT.` folder which are input files needed by Wannier90 code. 
+
+4. Copy `.amn`, `.mmn`, `.eig`, `UNK` file to `wannier/` folder, to get the MLWFs by running Wannier90: 
+
+    ```
+    wannier90 diamond.win
+    ```
+
+Notes:
+
+- The ABACUS-wannier90 interface can be used in both PW and LCAO basis.
+- If you want to plot the Wannier function, you must set `wvfn_formatted = .true.` in `diamond.win`, otherwise Wannier90 code cannot read files generated by ABACUS because these files are not binary files.
diff --git a/docs/advanced/scf/construct_H.md b/docs/advanced/scf/construct_H.md
@@ -80,4 +80,10 @@ Here, we use a simple [example calculation](https://github.com/deepmodeling/abac
     ```
     This is normal and it will not affect the results of calculation.
 
-## DFT+U
+## DFT+*U*
+
+Conventional functionals, e.g., L(S)DA and GGAs, encounter failures in strongly correlated systems, usually characterized by partially filled *d*/*f* shells. These include transition metals (TM) and their oxides, rare-earth compounds, and actinides, to name a few, where L(S)DA/GGAs typically yield quantitatively or even qualitatively wrong results. To address this failure, an efficient and successful method named DFT+*U*, which inherits the efficiency of L(S)DA/GGA but gains the strength of the Hubbard model in describing the physics of strongly correlatedsystems, has been developed.
+
+Now the DFT+*U* method is accessible in ABACUS. The details of the DFT+*U* method could be found in this [paper](https://doi.org/10.1063/5.0090122). It should be noted that the DFT+*U* works only within the NAO scheme, which means that the value of the keyword `basis_type` must be lcao when DFT+*U* is called. To turn on DFT+*U*, users need to set the value of the `dft_plus_u` keyword in the `INPUT` file to be 1. All relevant parmeters used in DFT+*U* calculations are listed in the [DFT+*U* correction](../input_files/input-main.md#DFTU-correction) part of the [list of keywords](../input_files/input-main.md).
+
+Examples of DFT+*U* calculations are provided in this [directory](https://github.com/deepmodeling/abacus-develop/tree/develop/examples/dft_plus_u).
diff --git a/examples/dft_plus_u/INPUT b/examples/dft_plus_u/INPUT
@@ -0,0 +1,36 @@
+INPUT_PARAMETERS
+suffix    NiO
+ntype    2
+nbands    40
+
+calculation    scf
+ecutwfc    100
+scf_thr    1.0e-6
+scf_nmax    200
+out_chg    1
+
+#init_chg    file
+#out_dos    1
+#dos_sigma    0.07
+#out_band    1
+
+smearing_method    gaussian
+smearing_sigma    0.05
+
+
+mixing_type    pulay
+mixing_beta    0.4
+#mixing_gg0    1.5
+
+ks_solver    genelpa
+basis_type    lcao
+gamma_only    0
+symmetry    0
+nspin    2
+
+#Parameter DFT+U
+dft_plus_u    1
+orbital_corr    2 2 -1
+hubbard_u    5.0 5.0 0.0
+hund_j    0.0 0.0 0.0
+
diff --git a/examples/dft_plus_u/KPT b/examples/dft_plus_u/KPT
@@ -0,0 +1,4 @@
+K_POINTS
+0
+Gamma
+7 7 7 0.0 0.0 0.0
diff --git a/examples/dft_plus_u/STRU b/examples/dft_plus_u/STRU
@@ -0,0 +1,35 @@
+ATOMIC_SPECIES
+Ni1 1.000 ../../tests/PP_ORB/Ni_ONCV_PBE-1.0.upf
+Ni2 1.000 ../../tests/PP_ORB/Ni_ONCV_PBE-1.0.upf
+O   1.000 ../../tests/PP_ORB/O_ONCV_PBE-1.0.upf
+
+NUMERICAL_ORBITAL
+../../tests/PP_ORB/Ni_gga_9au_100Ry_4s2p2d1f.orb
+../../tests/PP_ORB/Ni_gga_9au_100Ry_4s2p2d1f.orb
+../../tests/PP_ORB/O_gga_7au_100Ry_2s2p1d.orb
+
+LATTICE_CONSTANT
+7.884
+
+LATTICE_VECTORS
+ 1.00    0.50     0.50
+ 0.50    1.00     0.50
+ 0.50    0.50     1.00
+ATOMIC_POSITIONS
+Direct
+
+Ni1
+2.0
+1
+0.00            0.00            0.00            1        1       1
+
+Ni2
+-2.0
+1
+0.50            0.50            0.50            1        1       1
+
+O
+0.0
+2
+0.25            0.25            0.25            1        1       1
+0.75            0.75            0.75            1        1       1
diff --git a/examples/interface_wannier90/ABACUS_towannier90/INPUT-nscf b/examples/interface_wannier90/ABACUS_towannier90/INPUT-nscf
@@ -1,13 +1,15 @@
 INPUT_PARAMETERS
 
-ntype			1	
-ecutwfc			40
-nbands          	4
-calculation     	nscf
-scf_nmax			40
-pw_diag_thr                    1.0e-12
-scf_thr             	1.0e-15
-init_chg            file
-out_band                1
+pseudo_dir              ../../../tests/PP_ORB
+orbital_dir             ../../../tests/PP_ORB
+ntype                   1
+ecutwfc                 50
+nbands                  4
+calculation             nscf
+scf_nmax                50
+pw_diag_thr             1.0e-12
+scf_thr                 1.0e-15
+init_chg                file
+symmetry                0
 towannier90             1
-nnkpfile                seedname.nnkp
+nnkpfile                diamond.nnkp
diff --git a/examples/interface_wannier90/ABACUS_towannier90/INPUT-scf b/examples/interface_wannier90/ABACUS_towannier90/INPUT-scf
@@ -1,7 +1,9 @@
 INPUT_PARAMETERS
 
-ntype			1	
-ecutwfc			40
+pseudo_dir              ../../../tests/PP_ORB
+orbital_dir             ../../../tests/PP_ORB
+ntype                   1
+ecutwfc                 50
 calculation             scf
-scf_thr                     1e-13
-out_chg              1
+scf_thr                 1e-13
+out_chg                 1
diff --git a/examples/interface_wannier90/ABACUS_towannier90/clear.sh b/examples/interface_wannier90/ABACUS_towannier90/clear.sh
diff --git a/...annier90/ABACUS_towannier90/seedname.nnkp → ...wannier90/ABACUS_towannier90/diamond.nnkp b/...annier90/ABACUS_towannier90/seedname.nnkp → ...wannier90/ABACUS_towannier90/diamond.nnkp
diff --git a/.../interface_wannier90/wannier/seedname.win → ..._wannier90/ABACUS_towannier90/diamond.win b/.../interface_wannier90/wannier/seedname.win → ..._wannier90/ABACUS_towannier90/diamond.win
diff --git a/examples/interface_wannier90/ABACUS_towannier90/sub.sh b/examples/interface_wannier90/ABACUS_towannier90/sub.sh
diff --git a/examples/interface_wannier90/Read_before_test/note.txt b/examples/interface_wannier90/Read_before_test/note.txt
diff --git a/examples/interface_wannier90/Readme.md b/examples/interface_wannier90/Readme.md
@@ -0,0 +1,7 @@
+ABACUS to wannier90 step:
+
+1. You need a wannier90 input file `diamond.win`. then,you should run `wannier90 -pp diamond.win` and you will get `diamond.nnkp` file that is the necessary file for ABACUS.
+2. Prepare ABACUS's `INPUT`, `STRU`, `KPT`, pseudopotential file, `diamond.nnkp` file. Especially pay attention to kpt file in nscf calculation.
+3. First, you should do scf calculation, just like you always do. This step doesn't need `diamond.nnkp` file. Second, you should do nscf calculation, the kpt file is similar to "begin kpoints ..." in the `diamond.win` file.
+4. After you run ABACUS nscf calculation, you will get `diamond.amn`, `diamond.mmn` and `UNK*` file in OUT.* folder.
+5. Copy `diamond.amn`, `diamond.mmn`, `UNK*` file to wannier folder, then you can run `wannier90 diamond.win`, and you will get MLWF. There is a important things, you must set `wvfn_formatted = .true.` in `diamond.win`, otherwise the wannier90 code cannot read ABACUS file.