mohanchen · Oct 2, 2026 · Sep 25, 2026 · Sep 25, 2026 · Sep 25, 2026 · Sep 25, 2026
diff --git a/docs/advanced/elec_properties/dos.md b/docs/advanced/elec_properties/dos.md
@@ -9,7 +9,7 @@ We first, do a ground-state energy calculation ***with one additional keyword "[
 out_chg              1
 ```
 
-this will produce the converged charge density, which is contained in the file SPIN1_CHG.cube.
+this will produce the converged charge density, which is contained in the file chg.cube.
 Then, use the same `STRU` file, pseudopotential file and atomic orbital file (and the local density matrix file dm_onsite.txt if DFT+U is used) to do a non-self-consistent calculation. In this example, the potential is constructed from the ground-state charge density from the proceeding calculation. Now the INPUT file is like:
 
 ```
@@ -48,7 +48,7 @@ Some parameters in the INPUT file are explained:
   For LCAO calculations, this parameter will be neglected !
 - init_chg
 
-  the type of starting density. When doing scf calculation, this variable can be set ”atomic”. When doing nscf calculation, the charge density already exists(eg. in SPIN1_CHG.cube), and the variable should be set as ”file”. It means the density will be read from the existing file SPIN1_CHG.cube. For its more information please see the [here](https://abacus-rtd--1282.org.readthedocs.build/en/1282/advanced/input_files/input-main.html#init_chg).
+  the type of starting density. When doing scf calculation, this variable can be set ”atomic”. When doing nscf calculation, the charge density already exists(eg. in chg.cube), and the variable should be set as ”file”. It means the density will be read from the existing file chg.cube. For its more information please see the [here](https://abacus-rtd--1282.org.readthedocs.build/en/1282/advanced/input_files/input-main.html#init_chg).
 - out_dos
 
   output density of state(DOS). The unit of DOS is `(number of states)/(eV * unitcell)`. For its more information please see the [here](https://abacus-rtd--1282.org.readthedocs.build/en/1282/advanced/input_files/input-main.html#out_dos).
@@ -68,52 +68,30 @@ Gamma
 8 8 8 0 0 0
 ```
 
-Run the program, and you will see a file named DOS1_smearing.dat in the output directory. The first two columns in the file are the energy and DOS, respectively, and the third column is the sum of DOS. Plot file DOS1_smearing.dat with graphing software, and you’ll get the DOS.
+Run the program, and you will see a file named doss1g1_nao.txt in the output directory. The columns are: energy(eV), dos (number of states in each energy bin), dos_int (integrated DOS), dos_smear(1/eV) (Gaussian-smeared DOS), dos_smear_int (integrated smeared DOS). Plot the file with graphing software, and you'll get the DOS.
 
 ```
-            -5.49311           0.0518133           0.0518133
-            -5.48311           0.0641955            0.116009
-            -5.47311           0.0779299            0.193939
-            -5.46311           0.0926918             0.28663
-            -5.45311            0.108023            0.394653
-            -5.44311            0.123346            0.517999
-            ...
+#   energy(eV)           dos        dos_int  dos_smear(1/eV)  dos_smear_int
+      -5.49311       0.0518133       0.0518133       0.0518133       0.0518133
+      -5.48311       0.0641955        0.116009       0.0641955        0.116009
+      -5.47311       0.0779299        0.193939       0.0779299        0.193939
+      ...
 ```
 
 ## PDOS
 
-Along with the DOS1_smearing.dat file, we also produce the projected density of states (PDOS) in a file called PDOS.
+Along with the DOS files, we also produce the projected density of states (PDOS) in files named pdoss{spin}g{geom}_{basis}.txt (e.g., pdoss1g1_nao.txt).
 
-The PDOS file starts with number of atomic orbitals in the system, then a list of energy values, such as:
+The PDOS file uses a plain-text format. Each row corresponds to one energy point and one atom. Columns: energy(eV), atom (1-based), species, then pdos values ordered as s(1), p(3), d(5), f(7), etc. Zeta components are summed, and values below 1e-6 are zeroed out.
 
 ```
-<pdos>
-<nspin>1</nspin>
-<norbitals>26</norbitals>
-<energy_values units="eV">
-            -5.50311
-            -5.49311
-            -5.48311
-            -5.47311
-...
-
+# istep: 1
+# npoints: 2736
+# energy(eV)  atom  species  pdos(1/eV), columns: s(m=0) p(m=0,+1,-1) d(m=0,+1,-1,+2,-2) f(m=0,+1,-1,+2,-2,+3,-3)
+  -55.607730    1     Fe    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000
+  ...
 ```
 
-The rest of the fileis 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 `# istep:` line is only present for geometry-optimization (istep) output. Within each angular momentum `l`, the `2l+1` columns follow the magnetic-quantum-number order `m = 0, +1, -1, +2, -2, ...`.
 
-which tells the atom and symmetry of the current atomic orbital, and followed by the PDOS values. The values can thus be plotted against the energies. The unit of PDOS is also `(number of states)/(eV * unitcell)`.
+For nspin=2, two files are written (pdoss1* and pdoss2*), one per spin channel. For nspin=4, the two spinor components are summed into a single file. The unit of PDOS is also `(number of states)/(eV * unitcell)`.
diff --git a/docs/advanced/input_files/input-main.md b/docs/advanced/input_files/input-main.md
@@ -2038,13 +2038,12 @@
 ### out_dos
 
 - **Type**: Integer
-- **Description**: Whether to output the density of states (DOS). For more information, refer to the dos.md.
+- **Description**: Whether to output the density of states (DOS) and projected density of states (PDOS). For more information, refer to the dos.md.
   - 0: no output
-  - 1: output the density of states (DOS)
-   - nspin=1 or 4: doss1g{geom}_{basis}.txt, where geom is the geometry index when cell changes or ions move while basis is either pw or nao.
-   - nspin=2: doss1g{geom}_{basis}.txt and doss2g{geom}_{basis}.txt for two spin channles.
-  - 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
+  - 1: output the density of states (DOS) and projected density of states (PDOS, LCAO only)
+   - nspin=1 or 4: doss1g{geom}_{basis}.txt and pdoss1g{geom}_{basis}.txt, where geom is the geometry index when cell changes or ions move while basis is either pw or nao.
+   - nspin=2: doss1/doss2 and pdoss1/pdoss2 files for two spin channels.
+  Note: values 2 and 3 are no longer supported. Setting out_dos to 2 or 3 will raise an error.
 - **Default**: 0
 
 ### out_ldos
@@ -2523,8 +2522,8 @@
 ### dos_edelta_ev
 
 - **Type**: Real
-- **Description**: The step size in writing Density of States (DOS)
-- **Default**: 0.01
+- **Description**: The step size in writing Density of States (DOS). The default value was changed from 0.01 to 0.03 eV.
+- **Default**: 0.03
 - **Unit**: eV
 
 ### dos_sigma

diff --git a/docs/advanced/interface/pyabacus.md b/docs/advanced/interface/pyabacus.md
@@ -192,10 +192,10 @@ if bands_file:
 |------|-------------|
 | `running_scf.log` | Main calculation log |
 | `BANDS_1.dat` | Band structure data |
-| `PDOS` | Projected density of states |
-| `CHARGE.cube` | Charge density in cube format |
-| `SPIN1_CHG.cube` | Spin-up charge density |
-| `SPIN2_CHG.cube` | Spin-down charge density |
+| `pdoss1_*.txt` | Projected density of states (one file per spin channel) |
+| `chg.cube` | Charge density in cube format |
+| `chgs1.cube` | Spin-up charge density |
+| `chgs2.cube` | Spin-down charge density |
 | `istate.info` | Band eigenvalues and occupations |
 | `kpoints` | K-point information |
 

diff --git a/docs/advanced/scf/spin.md b/docs/advanced/scf/spin.md
@@ -49,9 +49,9 @@ Note: different settings for "noncolin" and "lspinorb" correspond to different c
 **Special case**: `noncolin=0, lspinorb=1` is commonly used for non-magnetic materials with SOC effects (e.g., topological insulators, semiconductors with spin-orbit splitting). In this case, the magnetization is NOT automatically set, implying no magnetic moments in the system.
 
 ## For the continuation job
-- Continuation job for "nspin 1" need file "SPIN1_CHG.cube" which is generated by setting "out_chg=1" in task before. By setting "init_chg file" in new job's INPUT file, charge density will start from file but not atomic. 
-- Continuation job for "nspin 2" need files "SPIN1_CHG.cube" and "SPIN2_CHG.cube" which are generated by "out_chg 1" with "nspin 2", and refer to spin-up and spin-down charge densities respectively. It should be note that reading "SPIN1_CHG.cube" only for the continuation target magnetic moment job is not supported now.
-- Continuation job for "nspin 4" need files "SPIN%s_CHG.cube", where %s in {1,2,3,4}, which are generated by "out_chg 1" with any variable setting leading to 'nspin'=4, and refer to charge densities in Pauli spin matrixes. It should be note that reading charge density files printing by 'nspin'=2 case is supported, which means only $\sigma_{tot}$ and $\sigma_{z}$ are read.
+- Continuation job for "nspin 1" need file "chg.cube" which is generated by setting "out_chg=1" in task before. By setting "init_chg file" in new job's INPUT file, charge density will start from file but not atomic. 
+- Continuation job for "nspin 2" need files "chgs1.cube" and "chgs2.cube" which are generated by "out_chg 1" with "nspin 2", and refer to spin-up and spin-down charge densities respectively. It should be note that reading "chgs1.cube" only for the continuation target magnetic moment job is not supported now.
+- Continuation job for "nspin 4" need files "chgs%s.cube", where %s in {1,2,3,4}, which are generated by "out_chg 1" with any variable setting leading to 'nspin'=4, and refer to charge densities in Pauli spin matrixes. It should be note that reading charge density files printing by 'nspin'=2 case is supported, which means only $\sigma_{tot}$ and $\sigma_{z}$ are read.
 
 # SOC Effects
 ## SOC

diff --git a/docs/community/faq.md b/docs/community/faq.md
@@ -108,7 +108,7 @@ This is a [known issue](https://github.com/open-mpi/ompi/issues/4948) of OpenMPI
 
 **1. How to visualize charge density file?**
 
-The output file SPIN1_CHG.cube can be visualized by using VESTA.
+The output file chg.cube can be visualized by using VESTA.
 
 **2. How to change cif file directly to STRU file?**
 

diff --git a/docs/parameters.yaml b/docs/parameters.yaml
@@ -3357,13 +3357,12 @@ parameters:
     category: Output information
     type: Integer
     description: |
-      Whether to output the density of states (DOS). For more information, refer to the dos.md.
+      Whether to output the density of states (DOS) and projected density of states (PDOS). For more information, refer to the dos.md.
       * 0: no output
-      * 1: output the density of states (DOS)
-       * nspin=1 or 4: doss1g{geom}_{basis}.txt, where geom is the geometry index when cell changes or ions move while basis is either pw or nao.
-       * nspin=2: doss1g{geom}_{basis}.txt and doss2g{geom}_{basis}.txt for two spin channles.
-      * 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
+      * 1: output the density of states (DOS) and projected density of states (PDOS, LCAO only)
+       * nspin=1 or 4: doss1g{geom}_{basis}.txt and pdoss1g{geom}_{basis}.txt, where geom is the geometry index when cell changes or ions move while basis is either pw or nao.
+       * nspin=2: doss1/doss2 and pdoss1/pdoss2 files for two spin channels.
+      Note: values 2 and 3 are no longer supported. Setting out_dos to 2 or 3 will raise an error.
     default_value: "0"
     unit: ""
     availability: ""
@@ -3940,8 +3939,8 @@ parameters:
     category: Density of states
     type: Real
     description: |
-      The step size in writing Density of States (DOS)
-    default_value: "0.01"
+      The step size in writing Density of States (DOS). The default value was changed from 0.01 to 0.03 eV.
+    default_value: "0.03"
     unit: eV
     availability: ""
   - name: dos_sigma

diff --git a/examples/13_dos/01_pw_Al/run.sh b/examples/13_dos/01_pw_Al/run.sh
@@ -18,9 +18,8 @@ if [[ ! -f scf.output ]] ||
    [[ ! -f nscf.output ]] ||
    [[ ! -f OUT.ABACUS/running_scf.log ]] ||
    [[ ! -f OUT.ABACUS/running_nscf.log ]] ||
-   [[ ! -f OUT.ABACUS/DOS1 ]] ||
-   [[ ! -f OUT.ABACUS/DOS1_smearing.dat ]] ||
-   [[ ! -f OUT.ABACUS/SPIN1_CHG.cube ]] ||
+   [[ ! -f OUT.ABACUS/dos.txt ]] ||
+   [[ ! -f OUT.ABACUS/chg.cube ]] ||
    [[ ! ( "$(tail -1 OUT.ABACUS/running_scf.log)" == " Total  Time  :"* ) ]] ||
    [[ ! ( "$(tail -1 OUT.ABACUS/running_nscf.log)" == " Total  Time  :"* ) ]]
 then

diff --git a/examples/13_dos/02_lcao_Si2/run.sh b/examples/13_dos/02_lcao_Si2/run.sh
@@ -18,9 +18,8 @@ if [[ ! -f scf.output ]] ||
    [[ ! -f nscf.output ]] ||
    [[ ! -f OUT.ABACUS/running_scf.log ]] ||
    [[ ! -f OUT.ABACUS/running_nscf.log ]] ||
-   [[ ! -f OUT.ABACUS/DOS1 ]] ||
-   [[ ! -f OUT.ABACUS/DOS1_smearing.dat ]] ||
-   [[ ! -f OUT.ABACUS/SPIN1_CHG.cube ]] ||
+   [[ ! -f OUT.ABACUS/doss1g1_nao.txt ]] ||
+   [[ ! -f OUT.ABACUS/chg.cube ]] ||
    [[ ! ( "$(tail -1 OUT.ABACUS/running_scf.log)" == " Total  Time  :"* ) ]] ||
    [[ ! ( "$(tail -1 OUT.ABACUS/running_nscf.log)" == " Total  Time  :"* ) ]]
 then

diff --git a/examples/13_dos/README b/examples/13_dos/README
@@ -16,11 +16,10 @@ It requires two steps of calculation:
    - set 'out_dos' to '1'
    - set 'pw_diag_thr' to a suitable value when basis_type is 'pw'
 
-Finish! Then you can see DOS1_smearing.dat in OUT.ABACUS. The first column is energy,
-and the second is density of states. Plot it!
+Finish! Then you can see the DOS file in OUT.ABACUS (dos.txt for PW basis, doss1g{geom}_nao.txt for LCAO basis). The columns are: energy(eV), dos(1/eV), dos_int (integrated DOS), dos_smear(1/eV), dos_smear_int (integrated smeared DOS). Plot it!
 
 Note:
-1)If the "out_dos" is used, the DOS1_smearing.dat will appear in the directory.
+1)If the "out_dos" is used, the DOS file will appear in the directory.
 2)The charge density needed to be output in scf calculations and be read
 in non-self consistent calculations.
 
diff --git a/interfaces/Wannier90_interface/abacusw90/interface.py b/interfaces/Wannier90_interface/abacusw90/interface.py
@@ -230,7 +230,7 @@ def step0_run_scf(
         for d in [self.scf_dir, self.scf_dir / "OUT.ABACUS"]:
             if not d.exists():
                 continue
-            for pat in ["SPIN*_CHG.cube", "SPIN*_CHG", "*CHARGE-DENSITY.restart"]:
+            for pat in ["chg*.cube", "chgs*.cube", "*CHARGE-DENSITY.restart"]:
                 if list(d.glob(pat)):
                     chg_found = True
                     break
@@ -283,11 +283,8 @@ def _prepare_scf_files(self):
 
         found = set()
         patterns = [
-            "SPIN*_CHG.cube",
-            "SPIN*_CHG",
             "chg*.cube",
-            "CHG*.cube",
-            "CHG*",
+            "chgs*.cube",
             "*CHARGE-DENSITY.restart",
         ]
         for pat in patterns:
@@ -468,7 +465,7 @@ def run(self, run_scf: bool = True):
             for d in [self.scf_dir, self.scf_dir / "OUT.ABACUS"]:
                 if not d.exists():
                     continue
-                for pat in ["SPIN*_CHG.cube", "SPIN*_CHG", "*CHARGE-DENSITY.restart"]:
+                for pat in ["chg*.cube", "chgs*.cube", "*CHARGE-DENSITY.restart"]:
                     if list(d.glob(pat)):
                         has_chg = True
                         break

diff --git a/source/Makefile.Objects b/source/Makefile.Objects
@@ -655,11 +655,11 @@ OBJS_IO=module_parameter/input_conv.o\
     module_bessel/bessel_basis.o\
     module_output/cal_test.o\
     module_dos/write_dos_pw.o\
-    module_energy/nscf_fermi_surf.o\
     module_energy/write_bands.o\
     module_dos/cal_dos.o\
     module_dos/cal_pdos_gamma.o\
     module_dos/cal_pdos_multik.o\
+    module_dos/write_pdos_text.o\
     module_dos/cal_ldos.o\
     module_bessel/numerical_descriptor.o\
     module_bessel/numerical_basis.o\

diff --git a/source/source_cell/unitcell.h b/source/source_cell/unitcell.h
@@ -52,6 +52,14 @@ class UnitCell : public BaseCell {
     /// @brief Get npol
     inline const int& get_npol() const { return npol; }
 
+    /// @brief Set iat2iwt (and npol) directly; intended for unit tests that
+    /// build a minimal UnitCell without running the full indexing setup.
+    inline void set_iat2iwt_for_test(const std::vector<int>& iat2iwt_in, const int npol_in)
+    {
+        this->iat2iwt = iat2iwt_in;
+        this->npol = npol_in;
+    }
+
     /// @brief Indexing tools for ia and it
     /// @return true if the last out is reset
     template <typename Tiat, typename Tiait>

diff --git a/source/source_io/CMakeLists.txt b/source/source_io/CMakeLists.txt
@@ -13,7 +13,6 @@ list(APPEND objects
     module_ml/cal_mlkedf_desc.cpp
     module_dos/write_dos_pw.cpp
     module_energy/write_bands.cpp
-    module_energy/nscf_fermi_surf.cpp
     module_energy/write_eig_occ.cpp
     module_bessel/numerical_basis.cpp
     module_bessel/numerical_basis_overlap.cpp
@@ -66,6 +65,7 @@ if(ENABLE_LCAO)
       module_dos/write_dos_lcao.cpp
       module_dos/cal_pdos_gamma.cpp
       module_dos/cal_pdos_multik.cpp
+      module_dos/write_pdos_text.cpp
       module_output/write_orb_info.cpp
       module_energy/write_proj_band_lcao.cpp
       module_chgpot/get_pchg_lcao.cpp
@@ -151,6 +151,7 @@ if(BUILD_TESTING)
   add_subdirectory(module_parameter/test)
   add_subdirectory(test_serial)
   add_subdirectory(test)
+  add_subdirectory(module_dos/unittests)
 endif()
 
 if(ENABLE_JSON)

diff --git a/source/source_io/module_ctrl/ctrl_output_pw.cpp b/source/source_io/module_ctrl/ctrl_output_pw.cpp
@@ -157,6 +157,16 @@ void ModuleIO::ctrl_scf_pw(const int istep,
                                    inp.dos_edelta_ev,
                                    inp.dos_scale,
                                    inp.dos_sigma,
+                                   inp.nspin,
+                                   inp.out_dos,
+                                   PARAM.globalv.dos_setemax,
+                                   inp.dos_emax_ev,
+                                   PARAM.globalv.dos_setemin,
+                                   inp.dos_emin_ev,
+                                   PARAM.globalv.two_fermi,
+                                   inp.out_app_flag,
+                                   inp.bndpar,
+                                   PARAM.globalv.global_out_dir,
                                    GlobalV::ofs_running);
         }
     }
@@ -271,7 +281,21 @@ void ModuleIO::ctrl_runner_pw(UnitCell& ucell,
                               *stp.template get_psi_d<T, Device>(),
                               ctx,
                               para_grid,
-                              ucell);
+                              ucell,
+                              inp.out_ldos,
+                              inp.stm_bias,
+                              inp.nspin,
+                              PARAM.globalv.global_out_dir,
+                              PARAM.globalv.two_fermi,
+                              inp.nbands,
+                              inp.dos_edelta_ev,
+                              inp.dos_scale,
+                              PARAM.globalv.dos_setemax,
+                              inp.dos_emax_ev,
+                              PARAM.globalv.dos_setemin,
+                              inp.dos_emin_ev,
+                              inp.dos_sigma,
+                              inp.ldos_line);
     }
 
     //----------------------------------------------------------

diff --git a/source/source_io/module_ctrl/ctrl_runner_lcao.cpp b/source/source_io/module_ctrl/ctrl_runner_lcao.cpp
@@ -51,7 +51,9 @@ void ctrl_runner_lcao(UnitCell& ucell,      // unitcell
 	if (inp.out_ldos[0])
     {
         ModuleIO::Cal_ldos<TK>::cal_ldos_lcao(pelec->eferm, chr, dmat, kv,
-          pelec->ekb, pelec->wg, psi[0], pgrid, pv, gd, ucell);
+          pelec->ekb, pelec->wg, psi[0], pgrid, pv, gd, ucell,
+          inp.stm_bias, inp.nspin, PARAM.globalv.global_out_dir,
+          PARAM.globalv.two_fermi, inp.out_ldos[1]);
     }
 
     // 3) print out exchange-correlation potential

diff --git a/source/source_io/module_ctrl/ctrl_scf_lcao.cpp b/source/source_io/module_ctrl/ctrl_scf_lcao.cpp
@@ -170,6 +170,15 @@ void ModuleIO::ctrl_scf_lcao(UnitCell& ucell,
                                  inp.dos_sigma,
                                  out_app_flag,
                                  istep,
+                                 inp.nspin,
+                                 PARAM.globalv.nlocal,
+                                 PARAM.globalv.dos_setemax,
+                                 inp.dos_emax_ev,
+                                 PARAM.globalv.dos_setemin,
+                                 inp.dos_emin_ev,
+                                 PARAM.globalv.two_fermi,
+                                 inp.bndpar,
+                                 PARAM.globalv.global_out_dir,
                                  GlobalV::ofs_running);
     }