Update DOS codes and related changes with respect to 'magnetzation' variable and default results of 204_NO_KP_NC example, reconstruct esolver by adding a new file pw_setup.cpp in esolver

docs/advanced/input_files/input-main.md

@@ -1705,7 +1705,8 @@ These variables are used to control the output of properties.
 - **Description**: Whether to output the density of states (DOS). For more information, refer to the [dos.md](../elec_properties/dos.md).
   - 0: no output
   - 1: output the density of states (DOS)
-  - 2: (lcao-only) output the density of states (DOS) and the projected density of states (PDOS)
+  - 2: (LCAO) output the density of states (DOS) and the projected density of states (PDOS)
+  - 3: output the Fermi surface file (fermi.bxsf) in BXSF format that can be visualized by XCrySDen 
 - **Default**: 0
 
 ### out_ldos

source/Makefile.Objects

@@ -495,13 +495,13 @@ OBJS_IO=input_conv.o\
     bessel_basis.o\
     cal_test.o\
     write_dos_pw.o\
+    nscf_fermi_surf.o\
     nscf_band.o\
     cal_dos.o\
     cal_pdos_gamma.o\
     cal_pdos_multik.o\
     cal_ldos.o\
     cif_io.o\
-    dos_nao.o\
     numerical_descriptor.o\
     numerical_basis.o\
     numerical_basis_jyjy.o\
@@ -564,7 +564,6 @@ OBJS_IO_LCAO=cal_r_overlap_R.o\
       write_dos_lcao.o\
       write_proj_band_lcao.o\
       write_istate_info.o\
-      nscf_fermi_surf.o\
       get_pchg_lcao.o\
       get_wf_lcao.o\
       io_dmk.o\

source/module_esolver/esolver_ks_lcao.cpp

@@ -1,16 +1,19 @@
 #include "esolver_ks_lcao.h"
 
+#include "module_io/write_dos_lcao.h"       // write DOS and PDOS
+#include "module_io/write_proj_band_lcao.h" // projcted band structure
+
 #include "module_base/formatter.h"
 #include "module_base/global_variable.h"
 #include "module_base/tool_title.h"
 #include "module_elecstate/elecstate_tools.h"
+
 #include "module_elecstate/module_dm/cal_dm_psi.h"
 #include "module_hamilt_lcao/module_deltaspin/spin_constrain.h"
 #include "module_hamilt_lcao/module_dftu/dftu.h"
 #include "module_io/berryphase.h"
 #include "module_io/cal_ldos.h"
 #include "module_io/cube_io.h"
-#include "module_io/dos_nao.h"
 #include "module_io/io_dmk.h"
 #include "module_io/io_npz.h"
 #include "module_io/output_dmk.h"
@@ -36,7 +39,6 @@
 #include "module_io/write_vxc.hpp"
 #include "module_io/write_vxc_r.hpp"
 
-//--------------temporary----------------------------
 #include "module_base/global_function.h"
 #include "module_cell/module_neighbor/sltk_grid_driver.h"
 #include "module_elecstate/cal_ux.h"
@@ -387,6 +389,8 @@ void ESolver_KS_LCAO<TK, TR>::cal_stress(UnitCell& ucell, ModuleBase::matrix& st
 //! the 8th function of ESolver_KS_LCAO: after_all_runners
 //! mohan add 2024-05-11
 //------------------------------------------------------------------------------
+
+
 template <typename TK, typename TR>
 void ESolver_KS_LCAO<TK, TR>::after_all_runners(UnitCell& ucell)
 {
@@ -397,7 +401,7 @@ void ESolver_KS_LCAO<TK, TR>::after_all_runners(UnitCell& ucell)
 
     const int nspin0 = (PARAM.inp.nspin == 2) ? 2 : 1;
 
-    // 4) write projected band structure by jiyy-2022-4-20
+    // 4) write projected band structure
     if (PARAM.inp.out_proj_band)
     {
         ModuleIO::write_proj_band_lcao(this->psi, this->pv, this->pelec, this->kv, ucell, this->p_hamilt);
@@ -406,7 +410,7 @@ void ESolver_KS_LCAO<TK, TR>::after_all_runners(UnitCell& ucell)
     // 5) print out density of states (DOS)
     if (PARAM.inp.out_dos)
 	{
-		ModuleIO::out_dos_nao(this->psi,
+		ModuleIO::write_dos_lcao(this->psi,
 				this->p_hamilt,
 				this->pv,
 				ucell,
@@ -417,7 +421,8 @@ void ESolver_KS_LCAO<TK, TR>::after_all_runners(UnitCell& ucell)
 				this->pelec->wg,
 				PARAM.inp.dos_edelta_ev,
 				PARAM.inp.dos_scale,
-				PARAM.inp.dos_sigma);
+				PARAM.inp.dos_sigma,
+                GlobalV::ofs_running);
     }
 
     // out ldos

source/module_esolver/esolver_ks_lcaopw.cpp

@@ -4,7 +4,6 @@
 #include "module_io/input_conv.h"
 #include "module_io/nscf_band.h"
 #include "module_io/output_log.h"
-#include "module_io/write_dos_pw.h"
 #include "module_io/write_istate_info.h"
 #include "module_io/write_wfc_pw.h"
 

source/module_esolver/esolver_ks_pw.cpp

@@ -837,16 +837,17 @@ void ESolver_KS_PW<T, Device>::after_all_runners(UnitCell& ucell)
     //! 1) Compute density of states (DOS)
     if (PARAM.inp.out_dos)
     {
-        ModuleIO::write_dos_pw(this->pelec->ekb,
-                               this->pelec->wg,
-                               this->kv,
-                               PARAM.inp.nbands,
-                               this->pelec->eferm,
-                               PARAM.inp.dos_edelta_ev,
-                               PARAM.inp.dos_scale,
-                               PARAM.inp.dos_sigma,
-                               GlobalV::ofs_running);
-    }
+		ModuleIO::write_dos_pw(ucell, 
+				this->pelec->ekb,
+				this->pelec->wg,
+				this->kv,
+				PARAM.inp.nbands,
+				this->pelec->eferm,
+				PARAM.inp.dos_edelta_ev,
+				PARAM.inp.dos_scale,
+				PARAM.inp.dos_sigma,
+				GlobalV::ofs_running);
+	}
 
     // out ldos
     if (PARAM.inp.out_ldos[0])

source/module_esolver/lcao_after_scf.cpp

@@ -1,4 +1,7 @@
 #include "esolver_ks_lcao.h"
+#include "module_io/write_HS.h"         // write H and S matrices
+#include "module_rdmft/rdmft.h"         // RDMFT codes
+
 #include "module_base/formatter.h"
 #include "module_base/global_variable.h"
 #include "module_base/tool_title.h"
@@ -7,7 +10,6 @@
 #include "module_hamilt_lcao/module_dftu/dftu.h"
 #include "module_io/berryphase.h"
 #include "module_io/cube_io.h"
-#include "module_io/dos_nao.h"
 #include "module_io/io_dmk.h"
 #include "module_io/io_npz.h"
 #include "module_io/nscf_band.h"
@@ -19,16 +21,12 @@
 #include "module_io/to_qo.h"
 #include "module_io/to_wannier90_lcao.h"
 #include "module_io/to_wannier90_lcao_in_pw.h"
-#include "module_io/write_HS.h"
-#include "module_io/write_dmr.h"
+#include "module_io/write_dmr.h" 
 #include "module_io/write_elecstat_pot.h"
 #include "module_io/write_istate_info.h"
-#include "module_io/write_proj_band_lcao.h"
 #include "module_io/write_wfc_nao.h"
 #include "module_io/cal_pLpR.h"
 #include "module_parameter/parameter.h"
-
-//--------------temporary----------------------------
 #include "module_base/global_function.h"
 #include "module_cell/module_neighbor/sltk_grid_driver.h"
 #include "module_elecstate/cal_ux.h"
@@ -38,8 +36,6 @@
 #include "module_hamilt_lcao/module_dftu/dftu.h"
 #include "module_hamilt_pw/hamilt_pwdft/global.h"
 #include "module_io/print_info.h"
-
-// mohan add 20250302
 #include "module_hamilt_lcao/hamilt_lcaodft/operator_lcao/ekinetic_new.h"
 
 #include <memory>
@@ -59,12 +55,7 @@
 #include "module_elecstate/elecstate_lcao.h"
 #include "module_hamilt_lcao/hamilt_lcaodft/hamilt_lcao.h"
 #include "module_hsolver/hsolver_lcao.h"
-// function used by deepks
-// #include "module_elecstate/cal_dm.h"
-//---------------------------------------------------
 
-// test RDMFT
-#include "module_rdmft/rdmft.h"
 
 #include <iostream>
 

source/module_esolver/pw_others.cpp

@@ -5,7 +5,6 @@
 #include "module_io/input_conv.h"
 #include "module_io/nscf_band.h"
 #include "module_io/output_log.h"
-#include "module_io/write_dos_pw.h"
 #include "module_io/write_istate_info.h"
 #include "module_io/write_wfc_pw.h"
 

source/module_io/CMakeLists.txt

@@ -7,6 +7,7 @@ list(APPEND objects
     cif_io.cpp
     write_dos_pw.cpp
     nscf_band.cpp
+    nscf_fermi_surf.cpp
     write_istate_info.cpp
     numerical_basis.cpp
     numerical_basis_jyjy.cpp
@@ -52,15 +53,13 @@ if(ENABLE_LCAO)
       cal_pdos_multik.cpp
       write_orb_info.cpp
       write_proj_band_lcao.cpp
-      nscf_fermi_surf.cpp
       get_pchg_lcao.cpp
       get_wf_lcao.cpp
       read_wfc_nao.cpp
       read_wfc_lcao.cpp
       write_wfc_nao.cpp
       io_dmk.cpp
       write_dmr.cpp
-      dos_nao.cpp
       sparse_matrix.cpp
       file_reader.cpp
       csr_reader.cpp

source/module_io/nscf_fermi_surf.cpp

@@ -15,13 +15,13 @@ void ModuleIO::nscf_fermi_surface(const std::string &out_band_dir,
 	ModuleBase::timer::tick("ModuleIO", "nscf_fermi_surface");
 #ifdef __MPI
 
-	int start = 1;
-	int end = PARAM.inp.nbands;
+	const int start = 1;
+	const int end = PARAM.inp.nbands;
 
 	std::ofstream ofs;
 	if(GlobalV::MY_RANK==0)
 	{
-		ofs.open(out_band_dir.c_str());//make the file clear!!
+		ofs.open(out_band_dir.c_str());
 		ofs << std::setprecision(6);
 		ofs.close();	
 	}

source/module_io/nscf_fermi_surf.h

@@ -1,5 +1,6 @@
 #ifndef NSCF_FERMI_SURF_H
 #define NSCF_FERMI_SURF_H
+
 #include "module_base/matrix.h"
 #include "module_cell/klist.h"
 #include "module_cell/unitcell.h"

source/module_io/test/CMakeLists.txt

@@ -62,7 +62,7 @@ AddTest(
 AddTest(
   TARGET io_write_dos_pw
   LIBS parameter  ${math_libs} base device symmetry 
-  SOURCES write_dos_pw_test.cpp ../cal_dos.cpp ../write_dos_pw.cpp ../output.cpp ../../module_cell/parallel_kpoints.cpp ../../module_cell/klist.cpp
+  SOURCES write_dos_pw_test.cpp ../cal_dos.cpp ../write_dos_pw.cpp ../output.cpp ../../module_cell/parallel_kpoints.cpp ../../module_cell/klist.cpp ../nscf_fermi_surf.cpp
 )
 
 AddTest(

source/module_io/test/write_dos_pw_test.cpp

@@ -68,7 +68,12 @@ TEST_F(DosPWTest,Dos1)
     elecstate::efermi fermi_energy;
 
 	std::ofstream ofs("write_dos_pw.log");
-	ModuleIO::write_dos_pw(dosp.ekb,
+
+    UnitCell ucell;
+
+	ModuleIO::write_dos_pw(
+            ucell, // this should be unitcell, 2025-04-12
+            dosp.ekb,
 			dosp.wg,
 			*kv,
 			PARAM.inp.nbands,
@@ -134,7 +139,12 @@ TEST_F(DosPWTest,Dos2)
     elecstate::efermi fermi_energy;
 
 	std::ofstream ofs("write_dos_pw.log");
-	ModuleIO::write_dos_pw(dosp.ekb,
+
+    UnitCell ucell;
+
+	ModuleIO::write_dos_pw(
+			ucell,
+			dosp.ekb,
 			dosp.wg,
 			*kv,
             PARAM.inp.nbands,