Merge branch 'develop' of https://github.com/deepmodeling/abacus-develop into develop

Author: YuLiu98

docs/examples/band-struc.md

@@ -27,6 +27,7 @@ pw_diag_thr 1.0e-7
 #Parameters (File)
 init_chg file
 out_band 1
+out_proj_band 1
 
 #Parameters (Smearing)
 smearing_method gaussian
@@ -58,4 +59,35 @@ points.
 Run the program, and you will see a file named BANDS_1.dat in the output directory. Plot it
 to get energy band structure.
 
-[back to top](#band-structure)
+If "out_proj_band" set 1, it will also produce the projected band structure in a file called PBAND_1 in xml format.
+
+The PBAND_1 file starts with number of atomic orbitals in the system, the text contents of element <band structure> is the same as data in the BANDS_1.dat file, such as:
+```
+<pband>
+<nspin>1</nspin>
+<norbitals>153</norbitals>
+<band_structure nkpoints="96" nbands="50" units="eV">
+...
+
+```
+
+The rest of the files arranged in sections, each section with a header such as below:
+    
+```
+<orbital
+ index="                                       1"
+ atom_index="                                       1"
+ species="Si"
+ l="                                       0"
+ m="                                       0"
+ z="                                       1"
+>
+<data>
+...
+</data>
+
+```
+
+The shape of text contents of element <data> is (Number of k-points, Number of bands) 
+ 
+[back to top](#band-structure)

docs/input-main.md

@@ -32,7 +32,7 @@
 
 - [Variables related to output information](#variables-related-to-output-information)
 
-    [out_force](#out_force) | [out_mul](#out_mul) | [out_freq_elec](#out_freq_elec) | [out_freq_ion](#out_freq_ion) | [out_chg](#out_chg) | [out_pot](#out_pot) | [out_dm](#out-dm) | [out_wfc_pw](#out_wfc_pw) | [out_wfc_r](#out_wfc_r) | [out_wfc_lcao](#out_wfc_lcao) | [out_dos](#out-dos) | [out_band](#out-band) | [out_stru](#out-stru) | [out_level](#out_level) | [out_alllog](#out-alllog) | [out_mat_hs](#out_mat_hs) | [out_mat_r](#out_mat_r) | [out_mat_hs2](#out_mat_hs2) | [out_element_info](#out-element-info) | [restart_save](#restart_save) | [restart_load](#restart_load)
+    [out_force](#out_force) | [out_mul](#out_mul) | [out_freq_elec](#out_freq_elec) | [out_freq_ion](#out_freq_ion) | [out_chg](#out_chg) | [out_pot](#out_pot) | [out_dm](#out-dm) | [out_wfc_pw](#out_wfc_pw) | [out_wfc_r](#out_wfc_r) | [out_wfc_lcao](#out_wfc_lcao) | [out_dos](#out-dos) | [out_band](#out-band) | [out_proj_band](#out-proj-band) | [out_stru](#out-stru) | [out_level](#out_level) | [out_alllog](#out-alllog) | [out_mat_hs](#out_mat_hs) | [out_mat_r](#out_mat_r) | [out_mat_hs2](#out_mat_hs2) | [out_element_info](#out-element-info) | [restart_save](#restart_save) | [restart_load](#restart_load)
 
 - [Density of states](#density-of-states)
 
@@ -735,6 +735,12 @@ This part of variables are used to control the output of properties.
 - **Description**: Controls whether to output the band structure. For mroe information, refer to the [worked example](examples/band-struc.md)
 - **Default**: 0
 
+#### out_proj_band
+
+- **Type**: Integer
+- **Description**: Controls whether to output the projected band structure. For mroe information, refer to the [worked example](examples/band-struc.md)
+- **Default**: 0
+
 #### out_stru
 
 - **Type**: Boolean

source/Makefile.Objects

@@ -121,6 +121,7 @@ gint_k_init.o\
 gint_k_vl.o\
 gint_k_rho.o\
 gint_k_fvl.o\
+gint_k_fvl_new.o\
 gint_k_env.o \
 ORB_control.o\
 ORB_read.o\

source/input.cpp

@@ -268,6 +268,7 @@ void Input::Default(void)
     out_wfc_r = 0;
     out_dos = 0;
     out_band = 0;
+    out_proj_band = 0;
     out_mat_hs = 0;
     out_mat_hs2 = 0; // LiuXh add 2019-07-15
     out_mat_r = 0; // jingan add 2019-8-14
@@ -995,6 +996,10 @@ bool Input::Read(const std::string &fn)
         {
             read_value(ifs, out_band);
         }
+        else if (strcmp("out_proj_band", word) == 0)
+        {
+            read_value(ifs, out_proj_band);
+        }
 
         else if (strcmp("out_mat_hs", word) == 0)
         {
@@ -1941,6 +1946,7 @@ void Input::Bcast()
     Parallel_Common::bcast_int(out_wfc_r);
     Parallel_Common::bcast_int(out_dos);
     Parallel_Common::bcast_int(out_band);
+    Parallel_Common::bcast_int(out_proj_band);
     Parallel_Common::bcast_int(out_mat_hs);
     Parallel_Common::bcast_int(out_mat_hs2); // LiuXh add 2019-07-15
     Parallel_Common::bcast_int(out_mat_r); // jingan add 2019-8-14
@@ -2227,6 +2233,7 @@ void Input::Check(void)
         out_stru = 0;
         out_dos = 0;
         out_band = 0;
+        out_proj_band = 0;
         cal_force = 0;
         init_wfc = "file";
         init_chg = "atomic"; // useless,
@@ -2248,6 +2255,7 @@ void Input::Check(void)
         out_stru = 0;
         out_dos = 0;
         out_band = 0;
+        out_proj_band = 0;
         cal_force = 0;
         init_wfc = "file";
         init_chg = "atomic";

source/input.h

@@ -206,6 +206,7 @@ class Input
     int out_wfc_r; // 0: no; 1: yes
     int out_dos; // dos calculation. mohan add 20090909
     int out_band; // band calculation pengfei 2014-10-13
+    int out_proj_band; // projected band structure calculation jiyy add 2022-05-11
     int out_mat_hs; // output H matrix and S matrix in local basis.
     int out_mat_hs2; // LiuXh add 2019-07-16, output H(R) matrix and S(R) matrix in local basis.
     int out_mat_r; // jingan add 2019-8-14, output r(R) matrix.

source/input_conv.cpp

@@ -434,6 +434,7 @@ void Input_Conv::Convert(void)
     GlobalC::wf.out_wfc_r = INPUT.out_wfc_r;
     GlobalC::en.out_dos = INPUT.out_dos;
     GlobalC::en.out_band = INPUT.out_band;
+    GlobalC::en.out_proj_band = INPUT.out_proj_band;
 #ifdef __LCAO
     Local_Orbital_Charge::out_dm = INPUT.out_dm;
     Pdiag_Double::out_mat_hs = INPUT.out_mat_hs;

source/module_elecstate/elecstate_lcao.cpp

@@ -94,7 +94,7 @@ void ElecStateLCAO::psiToRho(const psi::Psi<std::complex<double>>& psi)
     //------------------------------------------------------------
 
     ModuleBase::GlobalFunc::NOTE("Calculate the charge on real space grid!");
-    this->uhm->GK.cal_rho_k(this->loc->DM_R);
+    this->uhm->GK.cal_rho_k(this->loc->DM_R, this->charge);
 
     this->charge->renormalize_rho();
 
@@ -142,7 +142,7 @@ void ElecStateLCAO::psiToRho(const psi::Psi<double>& psi)
     // calculate the charge density on real space grid.
     //------------------------------------------------------------
     ModuleBase::GlobalFunc::NOTE("Calculate the charge on real space grid!");
-    this->uhm->GG.cal_rho(this->loc->DM);
+    this->uhm->GG.cal_rho(this->loc->DM, this->charge);
 
     this->charge->renormalize_rho();
 

source/module_esolver/esolver_ks_lcao.cpp

@@ -143,6 +143,19 @@ namespace ModuleESolver
         int ion_step = 0;
         GlobalC::pot.init_pot(ion_step, GlobalC::pw.strucFac);
 
+        //------------------init Basis_lcao----------------------
+        // Init Basis should be put outside of Ensolver.
+        // * reading the localized orbitals/projectors
+        // * construct the interpolation tables.
+        this->Init_Basis_lcao(this->orb_con, inp, ucell);
+        //------------------init Basis_lcao----------------------
+
+        //------------------init Hamilt_lcao----------------------
+        // * allocate H and S matrices according to computational resources
+        // * set the 'trace' between local H/S and global H/S
+        this->LM.divide_HS_in_frag(GlobalV::GAMMA_ONLY_LOCAL, orb_con.ParaV);
+        //------------------init Hamilt_lcao----------------------
+
 #ifdef __MPI  
         // PLEASE simplify the Exx_Global interface
         // mohan add 2021-03-25
@@ -175,19 +188,7 @@ namespace ModuleESolver
 
         if (INPUT.dft_plus_dmft) GlobalC::dmft.init(INPUT, ucell);
 
-        //------------------init Basis_lcao----------------------
-        // Init Basis should be put outside of Ensolver.
-        // * reading the localized orbitals/projectors
-        // * construct the interpolation tables.
-        this->Init_Basis_lcao(this->orb_con, inp, ucell);
-        //------------------init Basis_lcao----------------------
-
-        //------------------init Hamilt_lcao----------------------
-        // * allocate H and S matrices according to computational resources
-        // * set the 'trace' between local H/S and global H/S
-        this->LM.divide_HS_in_frag(GlobalV::GAMMA_ONLY_LOCAL, orb_con.ParaV);
-        //------------------init Hamilt_lcao----------------------
-      //init Psi
+        //init Psi
         if (GlobalV::GAMMA_ONLY_LOCAL)
             this->LOWF.wfc_gamma.resize(GlobalV::NSPIN);
         else

source/module_esolver/esolver_ks_lcao_elec.cpp

@@ -124,11 +124,11 @@ namespace ModuleESolver
             // calculate the charge density
             if (GlobalV::GAMMA_ONLY_LOCAL)
             {
-                this->UHM.GG.cal_rho(this->LOC.DM);
+                this->UHM.GG.cal_rho(this->LOC.DM, (Charge*)(&GlobalC::CHR));
             }
             else
             {
-                this->UHM.GK.cal_rho_k(this->LOC.DM_R);
+                this->UHM.GK.cal_rho_k(this->LOC.DM_R, (Charge*)(&GlobalC::CHR));
             }
 
             // renormalize the charge density

source/module_orbital/ORB_gen_tables.cpp

@@ -140,7 +140,18 @@ void ORB_gen_tables::snap_psibeta_half(
 
 	//find number of projectors on atom R0
 	const int nproj = infoNL_.nproj[T0];
-	assert(nproj>0); // mohan add 2021-04-25
+	if(nproj==0)
+	{
+		if(calc_deri)
+		{
+			nlm.resize(4);
+		}
+		else
+		{
+			nlm.resize(1);
+		}
+		return;	
+	}
 
 	std::vector<bool> calproj;
 	calproj.resize(nproj);