Feature: Add a new parameter out_freq_td to control output frequency in RT-TDDFT (#6796)

* Refactor in a more elegant way

* Refactor ctrl_output_td

* Refactor ctrl_scf_lcao

* Refactor ctrl_output_fp

* Reformat ctrl_scf_lcao.h

Author: AsTonyshment

docs/advanced/input_files/input-main.md

@@ -135,6 +135,7 @@
     - [cell\_factor](#cell_factor)
   - [Output Variables](#variables-related-to-output-information)
     - [out\_freq\_ion](#out_freq_ion)
+    - [out\_freq\_td](#out_freq_td)
     - [out\_freq\_elec](#out_freq_elec)
     - [out\_chg](#out_chg)
     - [out\_pot](#out_pot)
@@ -535,7 +536,7 @@ These variables are used to control general system parameters.
   - ofdft: orbital-free density functional theory
   - tdofdft: time-dependent orbital-free density functional theory
   - sdft: [stochastic density functional theory](#electronic-structure-sdft)
-  - tddft: real-time time-dependent density functional theory (TDDFT)
+  - tddft: real-time time-dependent density functional theory (RT-TDDFT)
   - lj: Leonard Jones potential
   - dp: DeeP potential, see details in [md.md](../md.md#dpmd)
   - nep: Neuroevolution Potential, see details in [md.md](../md.md#nep)
@@ -1688,15 +1689,22 @@ These variables are used to control the output of properties.
 ### out_freq_ion
 
 - **Type**: Integer
-- **Description**: Control the interval to print information every few ion steps. These properties cover charge density, local potential, electrostatic potential, Hamiltonian matrix, overlap matrix, density matrix, Mulliken population analysis and so on.
+- **Description**: Controls the output interval in **ionic steps**. When set to a positive integer $N$, information such as charge density, local potential, electrostatic potential, Hamiltonian matrix, overlap matrix, density matrix, and Mulliken population analysis is printed every $N$ ionic steps.
 - **Default**: 0
-- **Note**: The integer indicates to print information every 'out_freq_ion' ion steps. 
+- **Note**: In RT-TDDFT calculations, this parameter is inactive; output frequency is instead controlled by [`out_freq_td`](#out_freq_td)—see its description for details.
+
+### out_freq_td
+
+- **Type**: Integer
+- **Description**: Controls the output interval in **completed electronic evolution steps** during RT-TDDFT calculations. When set to a positive integer $N$, detailed information (see [`out_freq_ion`](#out_freq_ion)) is printed every $N$ electron time-evolution steps (i.e., every $N$ `STEP OF ELECTRON EVOLVE`). For example, if you wish to output information once per ionic step, you should set `out_freq_td` equal to [`estep_per_md`](#estep_per_md), since one ionic step corresponds to [`estep_per_md`](#estep_per_md) electronic evolution steps.
+- **Default**: 0
+- **Note**: This parameter is **only active in RT-TDDFT mode** (`esolver_type = tddft`). It has no effect in ground-state calculations.
 
 ### out_freq_elec
 
 - **Type**: Integer
-- **Description**: Output the charge density (only binary format, controlled by [out_chg](#out_chg)), wavefunction (controlled by [out_wfc_pw](#out_wfc_pw)) per `out_freq_elec` electronic iterations. Note that they are always output when converged or reach the maximum iterations [scf_nmax](#scf_nmax).
-- **Default**: [scf_nmax](#scf_nmax)
+- **Description**: Output the charge density (only binary format, controlled by [`out_chg`](#out_chg)), wavefunction (controlled by [`out_wfc_pw`](#out_wfc_pw)) per `out_freq_elec` electronic iterations. Note that they are always output when converged or reach the maximum iterations [`scf_nmax`](#scf_nmax).
+- **Default**: [`scf_nmax`](#scf_nmax)
 
 ### out_chg
 

source/source_esolver/esolver_ks_lcao_tddft.cpp

@@ -107,7 +107,9 @@ void ESolver_KS_LCAO_TDDFT<TR, Device>::runner(UnitCell& ucell, const int istep)
     {
         estep_max = PARAM.inp.estep_per_md + 1;
     }
-    // reset laststep matrix and wfc, if any atom cross the boundary
+
+    // Reset laststep matrix and wfc, if any atom cross the boundary
+    // Apply a phase correction to H, S, and psi to keep consistency when atoms cross periodic boundaries
     const size_t len_hs_ik = use_tensor && use_lapack ? PARAM.globalv.nlocal * PARAM.globalv.nlocal : this->pv.nloc;
     module_rt::reset_matrix_boundary(ucell,
                                      this->kv,
@@ -116,6 +118,7 @@ void ESolver_KS_LCAO_TDDFT<TR, Device>::runner(UnitCell& ucell, const int istep)
                                      this->Sk_laststep,
                                      this->psi_laststep,
                                      len_hs_ik);
+
     for (int estep = 0; estep < estep_max; estep++)
     {
         // calculate total time step

source/source_io/ctrl_output_fp.cpp

@@ -1,9 +1,10 @@
-#include "source_io/ctrl_output_fp.h" // use ctrl_output_fp() 
+#include "source_io/ctrl_output_fp.h" // use ctrl_output_fp()
+
+#include "cube_io.h"                                  // use write_vdata_palgrid
 #include "source_estate/module_charge/symmetry_rho.h" // use Symmetry_rho
-#include "source_io/write_elecstat_pot.h" // use write_elecstat_pot 
+#include "source_hamilt/module_xc/xc_functional.h"    // use XC_Functional
+#include "source_io/write_elecstat_pot.h"             // use write_elecstat_pot
 #include "source_io/write_elf.h"
-#include "cube_io.h"  // use write_vdata_palgrid
-#include "source_hamilt/module_xc/xc_functional.h" // use XC_Functional
 
 #ifdef USE_LIBXC
 #include "source_io/write_libxc_r.h"
@@ -12,14 +13,14 @@
 namespace ModuleIO
 {
 
-void ctrl_output_fp(UnitCell& ucell, 
-		elecstate::ElecState* pelec,	
-        ModulePW::PW_Basis_Big* pw_big,
-        ModulePW::PW_Basis* pw_rhod,
-        Charge &chr,
-        surchem &solvent,
-        Parallel_Grid &para_grid,
-		const int istep)
+void ctrl_output_fp(UnitCell& ucell,
+                    elecstate::ElecState* pelec,
+                    ModulePW::PW_Basis_Big* pw_big,
+                    ModulePW::PW_Basis* pw_rhod,
+                    Charge& chr,
+                    surchem& solvent,
+                    Parallel_Grid& para_grid,
+                    const int istep)
 {
     ModuleBase::TITLE("ModuleIO", "ctrl_output_fp");
     ModuleBase::timer::tick("ModuleIO", "ctrl_output_fp");
@@ -29,59 +30,64 @@ void ctrl_output_fp(UnitCell& ucell,
     const int nspin = PARAM.inp.nspin;
     const std::string global_out_dir = PARAM.globalv.global_out_dir;
 
-
-    // print out the 'g' index when istep_in!=-1
+    // print out the 'g' index when istep_in != -1
     int istep_in = -1;
-    if (PARAM.inp.out_freq_ion>0) // default value of out_freq_ion is 0
+    if (PARAM.inp.esolver_type != "tddft" && PARAM.inp.out_freq_ion > 0) // default value of out_freq_ion is 0
     {
         if (istep % PARAM.inp.out_freq_ion == 0)
         {
             istep_in = istep;
         }
     }
+    else if (PARAM.inp.esolver_type == "tddft" && PARAM.inp.out_freq_td > 0) // default value of out_freq_td is 0
+    {
+        if (istep % PARAM.inp.out_freq_td == 0)
+        {
+            istep_in = istep;
+        }
+    }
 
     std::string geom_block;
-    if(istep_in==-1)
+    if (istep_in == -1)
     {
         // do nothing
     }
-    else if(istep_in>=0)
+    else if (istep_in >= 0)
     {
         geom_block = "g" + std::to_string(istep + 1);
     }
 
-
     // 4) write charge density
     if (PARAM.inp.out_chg[0] > 0)
     {
         for (int is = 0; is < nspin; ++is)
         {
             pw_rhod->real2recip(chr.rho_save[is], chr.rhog_save[is]);
 
-            std::string fn =PARAM.globalv.global_out_dir + "chg";
+            std::string fn = PARAM.globalv.global_out_dir + "chg";
 
             std::string spin_block;
-            if(nspin == 2 || nspin == 4)
+            if (nspin == 2 || nspin == 4)
             {
-                spin_block= "s" + std::to_string(is + 1);
+                spin_block = "s" + std::to_string(is + 1);
             }
-            else if(nspin == 1)
+            else if (nspin == 1)
             {
                 // do nothing
             }
 
             fn += spin_block + geom_block + ".cube";
 
             ModuleIO::write_vdata_palgrid(para_grid,
-                    chr.rho_save[is],
-                    is,
-                    nspin,
-                    istep_in,
-                    fn,
-                    pelec->eferm.get_efval(is),
-                    &(ucell),
-                    PARAM.inp.out_chg[1],
-                    1);
+                                          chr.rho_save[is],
+                                          is,
+                                          nspin,
+                                          istep_in,
+                                          fn,
+                                          pelec->eferm.get_efval(is),
+                                          &(ucell),
+                                          PARAM.inp.out_chg[1],
+                                          1);
 
             if (XC_Functional::get_ked_flag())
             {
@@ -90,13 +96,13 @@ void ctrl_output_fp(UnitCell& ucell,
                 fn += spin_block + geom_block + ".cube";
 
                 ModuleIO::write_vdata_palgrid(para_grid,
-                        chr.kin_r_save[is],
-                        is,
-                        nspin,
-                        istep,
-                        fn,
-                        pelec->eferm.get_efval(is),
-                        &(ucell));
+                                              chr.kin_r_save[is],
+                                              is,
+                                              nspin,
+                                              istep,
+                                              fn,
+                                              pelec->eferm.get_efval(is),
+                                              &(ucell));
             }
         }
     }
@@ -106,49 +112,49 @@ void ctrl_output_fp(UnitCell& ucell,
     {
         for (int is = 0; is < nspin; is++)
         {
-            std::string fn =PARAM.globalv.global_out_dir + "pot";
+            std::string fn = PARAM.globalv.global_out_dir + "pot";
 
             std::string spin_block;
-            if(nspin == 2 || nspin == 4)
+            if (nspin == 2 || nspin == 4)
             {
-                spin_block= "s" + std::to_string(is + 1);
+                spin_block = "s" + std::to_string(is + 1);
             }
-            else if(nspin == 1)
+            else if (nspin == 1)
             {
                 // do nothing
             }
 
             fn += spin_block + geom_block + ".cube";
 
             ModuleIO::write_vdata_palgrid(para_grid,
-                    pelec->pot->get_eff_v(is),
-                    is,
-                    nspin,
-                    istep_in,
-                    fn,
-                    0.0, // efermi
-                    &(ucell),
-                    3,  // precision
-                    0); // out_fermi
+                                          pelec->pot->get_eff_v(is),
+                                          is,
+                                          nspin,
+                                          istep_in,
+                                          fn,
+                                          0.0, // efermi
+                                          &(ucell),
+                                          3,  // precision
+                                          0); // out_fermi
         }
     }
     else if (PARAM.inp.out_pot == 2)
     {
-        std::string fn =PARAM.globalv.global_out_dir + "potes";
+        std::string fn = PARAM.globalv.global_out_dir + "potes";
         fn += geom_block + ".cube";
 
         ModuleIO::write_elecstat_pot(
 #ifdef __MPI
-                pw_big->bz,
-                pw_big->nbz,
+            pw_big->bz,
+            pw_big->nbz,
 #endif
-                fn,
-                istep,
-                pw_rhod,
-                &chr,
-                &(ucell),
-                pelec->pot->get_fixed_v(),
-                solvent);
+            fn,
+            istep,
+            pw_rhod,
+            &chr,
+            &(ucell),
+            pelec->pot->get_fixed_v(),
+            solvent);
     }
 
     // 6) write ELF
@@ -161,40 +167,39 @@ void ctrl_output_fp(UnitCell& ucell,
             srho.begin(is, chr, pw_rhod, ucell.symm);
         }
 
-        std::string out_dir =PARAM.globalv.global_out_dir;
+        std::string out_dir = PARAM.globalv.global_out_dir;
         ModuleIO::write_elf(
 #ifdef __MPI
-                pw_big->bz,
-                pw_big->nbz,
+            pw_big->bz,
+            pw_big->nbz,
 #endif
-                out_dir,
-                istep,
-                nspin,
-                chr.rho,
-                chr.kin_r,
-                pw_rhod,
-                para_grid,
-                &(ucell),
-                PARAM.inp.out_elf[1]);
+            out_dir,
+            istep,
+            nspin,
+            chr.rho,
+            chr.kin_r,
+            pw_rhod,
+            para_grid,
+            &(ucell),
+            PARAM.inp.out_elf[1]);
     }
 
 #ifdef USE_LIBXC
     // 7) write xc(r)
-    if(PARAM.inp.out_xc_r[0]>=0)
+    if (PARAM.inp.out_xc_r[0] >= 0)
     {
-        ModuleIO::write_libxc_r(
-                PARAM.inp.out_xc_r[0],
-                XC_Functional::get_func_id(),
-                pw_rhod->nrxx, // number of real-space grid
-                ucell.omega, // volume of cell
-                ucell.tpiba,
-                chr,
-                *pw_big,
-                *pw_rhod);
+        ModuleIO::write_libxc_r(PARAM.inp.out_xc_r[0],
+                                XC_Functional::get_func_id(),
+                                pw_rhod->nrxx, // number of real-space grid
+                                ucell.omega,   // volume of cell
+                                ucell.tpiba,
+                                chr,
+                                *pw_big,
+                                *pw_rhod);
     }
 #endif
 
     ModuleBase::timer::tick("ModuleIO", "ctrl_output_fp");
 }
 
-} // End ModuleIO
+} // namespace ModuleIO

source/source_io/ctrl_output_fp.h

@@ -1,19 +1,19 @@
-#ifndef CTRL_OUTPUT_FP_H 
-#define CTRL_OUTPUT_FP_H 
+#ifndef CTRL_OUTPUT_FP_H
+#define CTRL_OUTPUT_FP_H
 
 #include "source_estate/elecstate_lcao.h"
 
 namespace ModuleIO
 {
 
-	void ctrl_output_fp(UnitCell& ucell, 
-			elecstate::ElecState* pelec,	
-			ModulePW::PW_Basis_Big* pw_big,
-			ModulePW::PW_Basis* pw_rhod,
-			Charge &chr,
-			surchem &solvent,
-			Parallel_Grid &para_grid,
-			const int istep);
+void ctrl_output_fp(UnitCell& ucell,
+                    elecstate::ElecState* pelec,
+                    ModulePW::PW_Basis_Big* pw_big,
+                    ModulePW::PW_Basis* pw_rhod,
+                    Charge& chr,
+                    surchem& solvent,
+                    Parallel_Grid& para_grid,
+                    const int istep);
 
 }
 #endif

source/source_io/ctrl_output_td.cpp

@@ -59,9 +59,9 @@ void ctrl_output_td(const UnitCell& ucell,
     }
 
     // (3) Output file for restart
-    if (PARAM.inp.out_freq_ion > 0) // default value of out_freq_ion is 0
+    if (PARAM.inp.out_freq_td > 0) // default value of out_freq_td is 0
     {
-        if (istep % PARAM.inp.out_freq_ion == 0)
+        if (istep % PARAM.inp.out_freq_td == 0)
         {
             if (td_p != nullptr)
             {