rename _nks as _nk to avoid ambiguity

maki49 · maki49 · commit da53d20ed7ad · 2024-10-13T14:52:20.000+08:00
diff --git a/source/module_elecstate/module_dm/density_matrix.cpp b/source/module_elecstate/module_dm/density_matrix.cpp
@@ -26,10 +26,10 @@ DensityMatrix<TK, TR>::~DensityMatrix()
 
 template <typename TK, typename TR>
 DensityMatrix<TK, TR>::DensityMatrix(const Parallel_Orbitals* paraV_in, const int nspin, const std::vector<ModuleBase::Vector3<double>>& kvec_d, const int nk)
-    : _paraV(paraV_in), _nspin(nspin), _kvec_d(kvec_d), _nks((nk > 0 && nk <= _kvec_d.size()) ? nk : _kvec_d.size())
+    : _paraV(paraV_in), _nspin(nspin), _kvec_d(kvec_d), _nk((nk > 0 && nk <= _kvec_d.size()) ? nk : _kvec_d.size())
 {
     ModuleBase::TITLE("DensityMatrix", "DensityMatrix-MK");
-    const int nks = _nks * _nspin;
+    const int nks = _nk * _nspin;
     this->_DMK.resize(nks);
     for (int ik = 0; ik < nks; ik++)
     {
@@ -39,7 +39,7 @@ DensityMatrix<TK, TR>::DensityMatrix(const Parallel_Orbitals* paraV_in, const in
 }
 
 template <typename TK, typename TR>
-DensityMatrix<TK, TR>::DensityMatrix(const Parallel_Orbitals* paraV_in, const int nspin) :_paraV(paraV_in), _nspin(nspin), _kvec_d({ ModuleBase::Vector3<double>(0,0,0) }), _nks(1)
+DensityMatrix<TK, TR>::DensityMatrix(const Parallel_Orbitals* paraV_in, const int nspin) :_paraV(paraV_in), _nspin(nspin), _kvec_d({ ModuleBase::Vector3<double>(0,0,0) }), _nk(1)
 {
     ModuleBase::TITLE("DensityMatrix", "DensityMatrix-GO");
     this->_DMK.resize(_nspin);
@@ -236,7 +236,7 @@ template <typename TK, typename TR>
 TK* DensityMatrix<TK, TR>::get_DMK_pointer(const int ik) const
 {
 #ifdef __DEBUG
-    assert(ik < this->_nks * this->_nspin);
+    assert(ik < this->_nk * this->_nspin);
 #endif
     return const_cast<TK*>(this->_DMK[ik].data());
 }
@@ -246,7 +246,7 @@ template <typename TK, typename TR>
 void DensityMatrix<TK, TR>::set_DMK_pointer(const int ik, TK* DMK_in)
 {
 #ifdef __DEBUG
-    assert(ik < this->_nks * this->_nspin);
+    assert(ik < this->_nk * this->_nspin);
 #endif
     this->_DMK[ik].assign(DMK_in, DMK_in + this->_paraV->nrow * this->_paraV->ncol);
 }
@@ -257,17 +257,17 @@ void DensityMatrix<TK, TR>::set_DMK(const int ispin, const int ik, const int i,
 {
 #ifdef __DEBUG
     assert(ispin > 0 && ispin <= this->_nspin);
-    assert(ik >= 0 && ik < this->_nks);
+    assert(ik >= 0 && ik < this->_nk);
 #endif
     // consider transpose col=>row
-    this->_DMK[ik + this->_nks * (ispin - 1)][i * this->_paraV->nrow + j] = value;
+    this->_DMK[ik + this->_nk * (ispin - 1)][i * this->_paraV->nrow + j] = value;
 }
 
 // set _DMK element
 template <typename TK, typename TR>
 void DensityMatrix<TK, TR>::set_DMK_zero()
 {
-    for (int ik = 0; ik < _nspin * _nks; ik++)
+    for (int ik = 0; ik < _nspin * _nk; ik++)
     {
         ModuleBase::GlobalFunc::ZEROS(this->_DMK[ik].data(),
                                       this->_paraV->get_row_size() * this->_paraV->get_col_size());
@@ -282,17 +282,17 @@ TK DensityMatrix<TK, TR>::get_DMK(const int ispin, const int ik, const int i, co
     assert(ispin > 0 && ispin <= this->_nspin);
 #endif
     // consider transpose col=>row
-    return this->_DMK[ik + this->_nks * (ispin - 1)][i * this->_paraV->nrow + j];
+    return this->_DMK[ik + this->_nk * (ispin - 1)][i * this->_paraV->nrow + j];
 }
 
 // get _DMK nks, nrow, ncol
 template <typename TK, typename TR>
 int DensityMatrix<TK, TR>::get_DMK_nks() const
 {
 #ifdef __DEBUG
-    assert(this->_DMK.size() == _nks * _nspin);
+    assert(this->_DMK.size() == _nk * _nspin);
 #endif
-    return _nks * _nspin;
+    return _nk * _nspin;
 }
 
 template <typename TK, typename TR>
@@ -364,7 +364,7 @@ void DensityMatrix<TK, TR>::cal_DMR_test()
 {
     for (int is = 1; is <= this->_nspin; ++is)
     {
-        int ik_begin = this->_nks * (is - 1); // jump this->_nks for spin_down if nspin==2
+        int ik_begin = this->_nk * (is - 1); // jump this->_nk for spin_down if nspin==2
         hamilt::HContainer<TR>* tmp_DMR = this->_DMR[is - 1];
         // set zero since this function is called in every scf step
         tmp_DMR->set_zero();
@@ -396,7 +396,7 @@ void DensityMatrix<TK, TR>::cal_DMR_test()
 #endif
                 std::complex<TR> tmp_res;
                 // loop over k-points
-                for (int ik = 0; ik < this->_nks; ++ik)
+                for (int ik = 0; ik < this->_nk; ++ik)
                 {
                     // cal k_phase
                     // if TK==std::complex<double>, kphase is e^{ikR}
@@ -438,7 +438,7 @@ void DensityMatrix<std::complex<double>, double>::cal_DMR()
     int ld_hk2 = 2 * ld_hk;
     for (int is = 1; is <= this->_nspin; ++is)
     {
-        int ik_begin = this->_nks * (is - 1); // jump this->_nks for spin_down if nspin==2
+        int ik_begin = this->_nk * (is - 1); // jump this->_nk for spin_down if nspin==2
         hamilt::HContainer<double>* tmp_DMR = this->_DMR[is - 1];
         // set zero since this function is called in every scf step
         tmp_DMR->set_zero();
@@ -471,7 +471,7 @@ void DensityMatrix<std::complex<double>, double>::cal_DMR()
                 // loop over k-points
                 if (PARAM.inp.nspin != 4)
                 {
-                    for (int ik = 0; ik < this->_nks; ++ik)
+                    for (int ik = 0; ik < this->_nk; ++ik)
                     {
                         // cal k_phase
                         // if TK==std::complex<double>, kphase is e^{ikR}
@@ -517,7 +517,7 @@ void DensityMatrix<std::complex<double>, double>::cal_DMR()
                     std::vector<std::complex<double>> tmp_DMR(this->_paraV->get_col_size()
                                                                   * this->_paraV->get_row_size(),
                                                               std::complex<double>(0.0, 0.0));
-                    for (int ik = 0; ik < this->_nks; ++ik)
+                    for (int ik = 0; ik < this->_nk; ++ik)
                     {
                         // cal k_phase
                         // if TK==std::complex<double>, kphase is e^{ikR}
@@ -603,7 +603,7 @@ void DensityMatrix<std::complex<double>, double>::cal_DMR(const int ik)
     int ld_hk2 = 2 * ld_hk;
     for (int is = 1; is <= this->_nspin; ++is)
     {
-        int ik_begin = this->_nks * (is - 1); // jump this->_nks for spin_down if nspin==2
+        int ik_begin = this->_nk * (is - 1); // jump this->_nk for spin_down if nspin==2
         hamilt::HContainer<double>* tmp_DMR = this->_DMR[is - 1];
         // set zero since this function is called in every scf step
         tmp_DMR->set_zero();
@@ -693,14 +693,14 @@ void DensityMatrix<double, double>::cal_DMR()
     int ld_hk = this->_paraV->nrow;
     for (int is = 1; is <= this->_nspin; ++is)
     {
-        int ik_begin = this->_nks * (is - 1); // jump this->_nks for spin_down if nspin==2
+        int ik_begin = this->_nk * (is - 1); // jump this->_nk for spin_down if nspin==2
         hamilt::HContainer<double>* tmp_DMR = this->_DMR[is - 1];
         // set zero since this function is called in every scf step
         tmp_DMR->set_zero();
 
 #ifdef __DEBUG
         // assert(tmp_DMR->is_gamma_only() == true);
-        assert(this->_nks == 1);
+        assert(this->_nk == 1);
 #endif
 #ifdef _OPENMP
 #pragma omp parallel for
@@ -828,7 +828,7 @@ void DensityMatrix<TK, TR>::read_DMK(const std::string directory, const int ispi
     {
         for (int j = 0; j < this->_paraV->ncol; ++j)
         {
-            ifs >> this->_DMK[ik + this->_nks * (ispin - 1)][i * this->_paraV->ncol + j];
+            ifs >> this->_DMK[ik + this->_nk * (ispin - 1)][i * this->_paraV->ncol + j];
         }
     }
     ifs.close();
@@ -865,7 +865,7 @@ void DensityMatrix<double, double>::write_DMK(const std::string directory, const
             {
                 ofs << "\n";
             }
-            ofs << " " << this->_DMK[ik + this->_nks * (ispin - 1)][i * this->_paraV->ncol + j];
+            ofs << " " << this->_DMK[ik + this->_nk * (ispin - 1)][i * this->_paraV->ncol + j];
         }
     }
 
@@ -902,7 +902,7 @@ void DensityMatrix<std::complex<double>, double>::write_DMK(const std::string di
             {
                 ofs << "\n";
             }
-            ofs << " " << this->_DMK[ik + this->_nks * (ispin - 1)][i * this->_paraV->ncol + j].real();
+            ofs << " " << this->_DMK[ik + this->_nk * (ispin - 1)][i * this->_paraV->ncol + j].real();
         }
     }
 
diff --git a/source/module_elecstate/module_dm/density_matrix.h b/source/module_elecstate/module_dm/density_matrix.h
@@ -243,8 +243,8 @@ class DensityMatrix
 
     /**
      * @brief density matrix in k space, which is a vector[ik]
-     * DMK should be a [_nspin][_nks][i][j] matrix,
-     * whose size is _nspin * _nks * _paraV->get_nrow() * _paraV->get_ncol()
+     * DMK should be a [_nspin][_nk][i][j] matrix,
+     * whose size is _nspin * _nk * _paraV->get_nrow() * _paraV->get_ncol()
      */
     // std::vector<ModuleBase::ComplexMatrix> _DMK;
     std::vector<std::vector<TK>> _DMK;
@@ -268,10 +268,10 @@ class DensityMatrix
 
     /**
      * @brief real number of k-points
-     * _nks is not equal to _kv->get_nks() when spin-polarization is considered
-     * _nks = kv->_nks / nspin
+     * _nk is not equal to _kv->get_nks() when spin-polarization is considered
+     * _nk = kv->get_nks() / nspin when nspin=2
      */
-    int _nks = 0;   // comment by lunasea: I want to rename it as nk
+    int _nk = 0;
 
 
 };
diff --git a/source/module_io/td_current_io.cpp b/source/module_io/td_current_io.cpp
@@ -28,7 +28,7 @@ void ModuleIO::cal_tmp_DM(elecstate::DensityMatrix<std::complex<double>, double>
     int ld_hk = DM_real.get_paraV_pointer()->nrow;
     int ld_hk2 = 2 * ld_hk;
     // tmp for is
-    int ik_begin = DM_real.get_DMK_nks() / nspin * (is - 1); // jump this->_nks for spin_down if nspin==2
+    int ik_begin = DM_real.get_DMK_nks() / nspin * (is - 1); // jump nk for spin_down if nspin==2
 
     hamilt::HContainer<double>* tmp_DMR_real = DM_real.get_DMR_vector()[is - 1];
     hamilt::HContainer<double>* tmp_DMR_imag = DM_imag.get_DMR_vector()[is - 1];
diff --git a/source/module_lr/dm_trans/dmr_complex.cpp b/source/module_lr/dm_trans/dmr_complex.cpp
@@ -11,7 +11,7 @@ namespace elecstate
         ModuleBase::timer::tick("DensityMatrix", "cal_DMR");
         for (int is = 1; is <= this->_nspin; ++is)
         {
-            int ik_begin = this->_nks * (is - 1); // jump this->_nks for spin_down if nspin==2
+            int ik_begin = this->_nk * (is - 1); // jump this->_nk for spin_down if nspin==2
             hamilt::HContainer<std::complex<double>>* tmp_DMR = this->_DMR[is - 1];
             // set zero since this function is called in every scf step
             tmp_DMR->set_zero();
@@ -43,7 +43,7 @@ namespace elecstate
 #endif
                     // loop over k-points
                     if (PARAM.inp.nspin != 4)
-                        for (int ik = 0; ik < this->_nks; ++ik)
+                        for (int ik = 0; ik < this->_nk; ++ik)
                         {
                             // cal k_phase
                             // if TK==std::complex<double>, kphase is e^{ikR}

Original file line number	Diff line number	Diff line change
`@@ -26,10 +26,10 @@ DensityMatrix<TK, TR>::~DensityMatrix()`
`26`	`26`
`27`	`27`	`template <typename TK, typename TR>`
`28`	`28`	`DensityMatrix<TK, TR>::DensityMatrix(const Parallel_Orbitals* paraV_in, const int nspin, const std::vector<ModuleBase::Vector3<double>>& kvec_d, const int nk)`
`29`		`- : _paraV(paraV_in), _nspin(nspin), _kvec_d(kvec_d), _nks((nk > 0 && nk <= _kvec_d.size()) ? nk : _kvec_d.size())`
	`29`	`+ : _paraV(paraV_in), _nspin(nspin), _kvec_d(kvec_d), _nk((nk > 0 && nk <= _kvec_d.size()) ? nk : _kvec_d.size())`
`30`	`30`	`{`
`31`	`31`	`ModuleBase::TITLE("DensityMatrix", "DensityMatrix-MK");`
`32`		`- const int nks = _nks * _nspin;`
	`32`	`+ const int nks = _nk * _nspin;`
`33`	`33`	`this->_DMK.resize(nks);`
`34`	`34`	`for (int ik = 0; ik < nks; ik++)`
`35`	`35`	`{`
`@@ -39,7 +39,7 @@ DensityMatrix<TK, TR>::DensityMatrix(const Parallel_Orbitals* paraV_in, const in`
`39`	`39`	`}`
`40`	`40`
`41`	`41`	`template <typename TK, typename TR>`
`42`		`-DensityMatrix<TK, TR>::DensityMatrix(const Parallel_Orbitals* paraV_in, const int nspin) :_paraV(paraV_in), _nspin(nspin), _kvec_d({ ModuleBase::Vector3<double>(0,0,0) }), _nks(1)`
	`42`	`+DensityMatrix<TK, TR>::DensityMatrix(const Parallel_Orbitals* paraV_in, const int nspin) :_paraV(paraV_in), _nspin(nspin), _kvec_d({ ModuleBase::Vector3<double>(0,0,0) }), _nk(1)`
`43`	`43`	`{`
`44`	`44`	`ModuleBase::TITLE("DensityMatrix", "DensityMatrix-GO");`
`45`	`45`	`this->_DMK.resize(_nspin);`
`@@ -236,7 +236,7 @@ template <typename TK, typename TR>`
`236`	`236`	`TK* DensityMatrix<TK, TR>::get_DMK_pointer(const int ik) const`
`237`	`237`	`{`
`238`	`238`	`#ifdef __DEBUG`
`239`		`- assert(ik < this->_nks * this->_nspin);`
	`239`	`+ assert(ik < this->_nk * this->_nspin);`
`240`	`240`	`#endif`
`241`	`241`	`return const_cast<TK*>(this->_DMK[ik].data());`
`242`	`242`	`}`
`@@ -246,7 +246,7 @@ template <typename TK, typename TR>`
`246`	`246`	`void DensityMatrix<TK, TR>::set_DMK_pointer(const int ik, TK* DMK_in)`
`247`	`247`	`{`
`248`	`248`	`#ifdef __DEBUG`
`249`		`- assert(ik < this->_nks * this->_nspin);`
	`249`	`+ assert(ik < this->_nk * this->_nspin);`
`250`	`250`	`#endif`
`251`	`251`	`this->_DMK[ik].assign(DMK_in, DMK_in + this->_paraV->nrow * this->_paraV->ncol);`
`252`	`252`	`}`
`@@ -257,17 +257,17 @@ void DensityMatrix<TK, TR>::set_DMK(const int ispin, const int ik, const int i,`
`257`	`257`	`{`
`258`	`258`	`#ifdef __DEBUG`
`259`	`259`	`assert(ispin > 0 && ispin <= this->_nspin);`
`260`		`- assert(ik >= 0 && ik < this->_nks);`
	`260`	`+ assert(ik >= 0 && ik < this->_nk);`
`261`	`261`	`#endif`
`262`	`262`	`// consider transpose col=>row`
`263`		`- this->_DMK[ik + this->_nks * (ispin - 1)][i * this->_paraV->nrow + j] = value;`
	`263`	`+ this->_DMK[ik + this->_nk * (ispin - 1)][i * this->_paraV->nrow + j] = value;`
`264`	`264`	`}`
`265`	`265`
`266`	`266`	`// set _DMK element`
`267`	`267`	`template <typename TK, typename TR>`
`268`	`268`	`void DensityMatrix<TK, TR>::set_DMK_zero()`
`269`	`269`	`{`
`270`		`- for (int ik = 0; ik < _nspin * _nks; ik++)`
	`270`	`+ for (int ik = 0; ik < _nspin * _nk; ik++)`
`271`	`271`	`{`
`272`	`272`	`ModuleBase::GlobalFunc::ZEROS(this->_DMK[ik].data(),`
`273`	`273`	`this->_paraV->get_row_size() * this->_paraV->get_col_size());`
`@@ -282,17 +282,17 @@ TK DensityMatrix<TK, TR>::get_DMK(const int ispin, const int ik, const int i, co`
`282`	`282`	`assert(ispin > 0 && ispin <= this->_nspin);`
`283`	`283`	`#endif`
`284`	`284`	`// consider transpose col=>row`
`285`		`- return this->_DMK[ik + this->_nks * (ispin - 1)][i * this->_paraV->nrow + j];`
	`285`	`+ return this->_DMK[ik + this->_nk * (ispin - 1)][i * this->_paraV->nrow + j];`
`286`	`286`	`}`
`287`	`287`
`288`	`288`	`// get _DMK nks, nrow, ncol`
`289`	`289`	`template <typename TK, typename TR>`
`290`	`290`	`int DensityMatrix<TK, TR>::get_DMK_nks() const`
`291`	`291`	`{`
`292`	`292`	`#ifdef __DEBUG`
`293`		`- assert(this->_DMK.size() == _nks * _nspin);`
	`293`	`+ assert(this->_DMK.size() == _nk * _nspin);`
`294`	`294`	`#endif`
`295`		`- return _nks * _nspin;`
	`295`	`+ return _nk * _nspin;`
`296`	`296`	`}`
`297`	`297`
`298`	`298`	`template <typename TK, typename TR>`
`@@ -364,7 +364,7 @@ void DensityMatrix<TK, TR>::cal_DMR_test()`
`364`	`364`	`{`
`365`	`365`	`for (int is = 1; is <= this->_nspin; ++is)`
`366`	`366`	`{`
`367`		`- int ik_begin = this->_nks * (is - 1); // jump this->_nks for spin_down if nspin==2`
	`367`	`+ int ik_begin = this->_nk * (is - 1); // jump this->_nk for spin_down if nspin==2`
`368`	`368`	`hamilt::HContainer<TR>* tmp_DMR = this->_DMR[is - 1];`
`369`	`369`	`// set zero since this function is called in every scf step`
`370`	`370`	`tmp_DMR->set_zero();`
`@@ -396,7 +396,7 @@ void DensityMatrix<TK, TR>::cal_DMR_test()`
`396`	`396`	`#endif`
`397`	`397`	`std::complex<TR> tmp_res;`
`398`	`398`	`// loop over k-points`
`399`		`- for (int ik = 0; ik < this->_nks; ++ik)`
	`399`	`+ for (int ik = 0; ik < this->_nk; ++ik)`
`400`	`400`	`{`
`401`	`401`	`// cal k_phase`
`402`	`402`	`// if TK==std::complex<double>, kphase is e^{ikR}`
`@@ -438,7 +438,7 @@ void DensityMatrix<std::complex<double>, double>::cal_DMR()`
`438`	`438`	`int ld_hk2 = 2 * ld_hk;`
`439`	`439`	`for (int is = 1; is <= this->_nspin; ++is)`
`440`	`440`	`{`
`441`		`- int ik_begin = this->_nks * (is - 1); // jump this->_nks for spin_down if nspin==2`
	`441`	`+ int ik_begin = this->_nk * (is - 1); // jump this->_nk for spin_down if nspin==2`
`442`	`442`	`hamilt::HContainer<double>* tmp_DMR = this->_DMR[is - 1];`
`443`	`443`	`// set zero since this function is called in every scf step`
`444`	`444`	`tmp_DMR->set_zero();`
`@@ -471,7 +471,7 @@ void DensityMatrix<std::complex<double>, double>::cal_DMR()`
`471`	`471`	`// loop over k-points`
`472`	`472`	`if (PARAM.inp.nspin != 4)`
`473`	`473`	`{`
`474`		`- for (int ik = 0; ik < this->_nks; ++ik)`
	`474`	`+ for (int ik = 0; ik < this->_nk; ++ik)`
`475`	`475`	`{`
`476`	`476`	`// cal k_phase`
`477`	`477`	`// if TK==std::complex<double>, kphase is e^{ikR}`
`@@ -517,7 +517,7 @@ void DensityMatrix<std::complex<double>, double>::cal_DMR()`
`517`	`517`	`std::vector<std::complex<double>> tmp_DMR(this->_paraV->get_col_size()`
`518`	`518`	`* this->_paraV->get_row_size(),`
`519`	`519`	`std::complex<double>(0.0, 0.0));`
`520`		`- for (int ik = 0; ik < this->_nks; ++ik)`
	`520`	`+ for (int ik = 0; ik < this->_nk; ++ik)`
`521`	`521`	`{`
`522`	`522`	`// cal k_phase`
`523`	`523`	`// if TK==std::complex<double>, kphase is e^{ikR}`
`@@ -603,7 +603,7 @@ void DensityMatrix<std::complex<double>, double>::cal_DMR(const int ik)`
`603`	`603`	`int ld_hk2 = 2 * ld_hk;`
`604`	`604`	`for (int is = 1; is <= this->_nspin; ++is)`
`605`	`605`	`{`
`606`		`- int ik_begin = this->_nks * (is - 1); // jump this->_nks for spin_down if nspin==2`
	`606`	`+ int ik_begin = this->_nk * (is - 1); // jump this->_nk for spin_down if nspin==2`
`607`	`607`	`hamilt::HContainer<double>* tmp_DMR = this->_DMR[is - 1];`
`608`	`608`	`// set zero since this function is called in every scf step`
`609`	`609`	`tmp_DMR->set_zero();`
`@@ -693,14 +693,14 @@ void DensityMatrix<double, double>::cal_DMR()`
`693`	`693`	`int ld_hk = this->_paraV->nrow;`
`694`	`694`	`for (int is = 1; is <= this->_nspin; ++is)`
`695`	`695`	`{`
`696`		`- int ik_begin = this->_nks * (is - 1); // jump this->_nks for spin_down if nspin==2`
	`696`	`+ int ik_begin = this->_nk * (is - 1); // jump this->_nk for spin_down if nspin==2`
`697`	`697`	`hamilt::HContainer<double>* tmp_DMR = this->_DMR[is - 1];`
`698`	`698`	`// set zero since this function is called in every scf step`
`699`	`699`	`tmp_DMR->set_zero();`
`700`	`700`
`701`	`701`	`#ifdef __DEBUG`
`702`	`702`	`// assert(tmp_DMR->is_gamma_only() == true);`
`703`		`- assert(this->_nks == 1);`
	`703`	`+ assert(this->_nk == 1);`
`704`	`704`	`#endif`
`705`	`705`	`#ifdef _OPENMP`
`706`	`706`	`#pragma omp parallel for`
`@@ -828,7 +828,7 @@ void DensityMatrix<TK, TR>::read_DMK(const std::string directory, const int ispi`
`828`	`828`	`{`
`829`	`829`	`for (int j = 0; j < this->_paraV->ncol; ++j)`
`830`	`830`	`{`
`831`		`- ifs >> this->_DMK[ik + this->_nks * (ispin - 1)][i * this->_paraV->ncol + j];`
	`831`	`+ ifs >> this->_DMK[ik + this->_nk * (ispin - 1)][i * this->_paraV->ncol + j];`
`832`	`832`	`}`
`833`	`833`	`}`
`834`	`834`	`ifs.close();`
`@@ -865,7 +865,7 @@ void DensityMatrix<double, double>::write_DMK(const std::string directory, const`
`865`	`865`	`{`
`866`	`866`	`ofs << "\n";`
`867`	`867`	`}`
`868`		`- ofs << " " << this->_DMK[ik + this->_nks * (ispin - 1)][i * this->_paraV->ncol + j];`
	`868`	`+ ofs << " " << this->_DMK[ik + this->_nk * (ispin - 1)][i * this->_paraV->ncol + j];`
`869`	`869`	`}`
`870`	`870`	`}`
`871`	`871`
`@@ -902,7 +902,7 @@ void DensityMatrix<std::complex<double>, double>::write_DMK(const std::string di`
`902`	`902`	`{`
`903`	`903`	`ofs << "\n";`
`904`	`904`	`}`
`905`		`- ofs << " " << this->_DMK[ik + this->_nks * (ispin - 1)][i * this->_paraV->ncol + j].real();`
	`905`	`+ ofs << " " << this->_DMK[ik + this->_nk * (ispin - 1)][i * this->_paraV->ncol + j].real();`
`906`	`906`	`}`
`907`	`907`	`}`
`908`	`908`