diff --git a/source/module_hamilt_lcao/module_deepks/deepks_force.cpp b/source/module_hamilt_lcao/module_deepks/deepks_force.cpp index 1457632a27b..bd392b6651d 100644 --- a/source/module_hamilt_lcao/module_deepks/deepks_force.cpp +++ b/source/module_hamilt_lcao/module_deepks/deepks_force.cpp @@ -102,9 +102,9 @@ void DeePKS_domain::cal_f_delta(const std::vector>& dm, int dRz = 0; if constexpr (std::is_same>::value) // for multi-k { - dRx = dR2.x - dR1.x; - dRy = dR2.y - dR1.y; - dRz = dR2.z - dR1.z; + dRx = dR1.x - dR2.x; + dRy = dR1.y - dR2.y; + dRz = dR1.z - dR2.z; } ModuleBase::Vector3 dR(dRx, dRy, dRz); diff --git a/source/module_hamilt_lcao/module_deepks/deepks_fpre.cpp b/source/module_hamilt_lcao/module_deepks/deepks_fpre.cpp index 2616fb5d1cb..61274c3d0d7 100644 --- a/source/module_hamilt_lcao/module_deepks/deepks_fpre.cpp +++ b/source/module_hamilt_lcao/module_deepks/deepks_fpre.cpp @@ -70,14 +70,15 @@ void DeePKS_domain::cal_gdmx(const int lmaxd, int dRz = 0; if constexpr (std::is_same>::value) { - dRx = (dR2 - dR1).x; - dRy = (dR2 - dR1).y; - dRz = (dR2 - dR1).z; + dRx = (dR1 - dR2).x; + dRy = (dR1 - dR2).y; + dRz = (dR1 - dR2).z; } ModuleBase::Vector3 dR(dRx, dRy, dRz); hamilt::AtomPair dm_pair(ibt1, ibt2, dRx, dRy, dRz, &pv); - dm_pair.allocate(nullptr, 1); + dm_pair.allocate(nullptr, true); + // not support nspin = 2 now for (int ik = 0; ik < nks; ik++) { TK kphase = TK(0); diff --git a/source/module_hamilt_lcao/module_deepks/deepks_orbpre.cpp b/source/module_hamilt_lcao/module_deepks/deepks_orbpre.cpp index d679177b5be..7d1631a92c1 100644 --- a/source/module_hamilt_lcao/module_deepks/deepks_orbpre.cpp +++ b/source/module_hamilt_lcao/module_deepks/deepks_orbpre.cpp @@ -197,15 +197,15 @@ void DeePKS_domain::cal_orbital_precalc(const std::vector& dm_hl, { hamilt::AtomPair dm_pair(ibt1, ibt2, - (dR2 - dR1).x, - (dR2 - dR1).y, - (dR2 - dR1).z, + (dR1 - dR2).x, + (dR1 - dR2).y, + (dR1 - dR2).z, &pv); dm_pair.allocate(&dm_array[ik * row_size * col_size], 0); const double arg - = -(kvec_d[ik] * ModuleBase::Vector3(dR2 - dR1)) * ModuleBase::TWO_PI; + = -(kvec_d[ik] * ModuleBase::Vector3(dR1 - dR2)) * ModuleBase::TWO_PI; double sinp, cosp; diff --git a/source/module_hamilt_lcao/module_deepks/deepks_pdm.cpp b/source/module_hamilt_lcao/module_deepks/deepks_pdm.cpp index 5374e82531e..54b23591da8 100644 --- a/source/module_hamilt_lcao/module_deepks/deepks_pdm.cpp +++ b/source/module_hamilt_lcao/module_deepks/deepks_pdm.cpp @@ -273,36 +273,51 @@ void DeePKS_domain::cal_pdm(bool& init_pdm, } } // prepare DM from DMR - std::vector dm_array(row_size * col_size, 0.0); - const double* dm_current = nullptr; int dRx = 0, dRy = 0, dRz = 0; if constexpr (std::is_same>::value) { - dRx = dR2.x - dR1.x; - dRy = dR2.y - dR1.y; - dRz = dR2.z - dR1.z; + dRx = dR1.x - dR2.x; + dRy = dR1.y - dR2.y; + dRz = dR1.z - dR2.z; } - // dm_k + ModuleBase::Vector3 dR(dRx, dRy, dRz); + + hamilt::AtomPair dm_pair(ibt1, ibt2, dRx, dRy, dRz, &pv); + dm_pair.allocate(nullptr, true); auto dm_k = dm->get_DMK_vector(); - const int nrow = pv.nrow; - for (int ir = 0; ir < row_size; ir++) + + if constexpr (std::is_same::value) // for gamma-only + { + for (int is = 0; is < dm_k.size(); is++) + { + if (ModuleBase::GlobalFunc::IS_COLUMN_MAJOR_KS_SOLVER(PARAM.inp.ks_solver)) + { + dm_pair.add_from_matrix(dm_k[is].data(), pv.get_row_size(), 1.0, 1); + } + else + { + dm_pair.add_from_matrix(dm_k[is].data(), pv.get_col_size(), 1.0, 0); + } + } + } + else // for multi-k { - for (int ic = 0; ic < col_size; ic++) + for (int ik = 0; ik < dm_k.size(); ik++) { - int iglob = (pv.atom_begin_row[ibt1] + ir) + nrow * (pv.atom_begin_col[ibt2] + ic); - int iloc = ir * col_size + ic; - std::complex tmp = 0.0; - for(int ik = 0; ik < dm_k.size(); ik++) // dm_k.size() == _nk * _nspin + const double arg = -(kvec_d[ik] * dR) * ModuleBase::TWO_PI; + const std::complex kphase = std::complex(cos(arg), sin(arg)); + if (ModuleBase::GlobalFunc::IS_COLUMN_MAJOR_KS_SOLVER(PARAM.inp.ks_solver)) + { + dm_pair.add_from_matrix(dm_k[ik].data(), pv.get_row_size(), kphase, 1); + } + else { - const double arg = (kvec_d[ik] * ModuleBase::Vector3(dR1 - dR2)) * ModuleBase::TWO_PI; - const std::complex kphase = std::complex(cos(arg), sin(arg)); - tmp += dm_k[ik][iglob] * kphase; + dm_pair.add_from_matrix(dm_k[ik].data(), pv.get_col_size(), kphase, 0); } - dm_array[iloc] += tmp.real(); } } - dm_current = dm_array.data(); + const double* dm_current = dm_pair.get_pointer(); // use s_2t and dm_current to get g_1dmt // dgemm_: C = alpha * A * B + beta * C // C = g_1dmt, A = dm_current, B = s_2t diff --git a/source/module_hamilt_lcao/module_deepks/deepks_spre.cpp b/source/module_hamilt_lcao/module_deepks/deepks_spre.cpp index dfc83cad1e2..5ec34735d79 100644 --- a/source/module_hamilt_lcao/module_deepks/deepks_spre.cpp +++ b/source/module_hamilt_lcao/module_deepks/deepks_spre.cpp @@ -78,9 +78,9 @@ void DeePKS_domain::cal_gdmepsl(const int lmaxd, int dRz = 0; if constexpr (std::is_same>::value) { - dRx = (dR2 - dR1).x; - dRy = (dR2 - dR1).y; - dRz = (dR2 - dR1).z; + dRx = (dR1 - dR2).x; + dRy = (dR1 - dR2).y; + dRz = (dR1 - dR2).z; } ModuleBase::Vector3 dR(dRx, dRy, dRz); diff --git a/source/module_hamilt_lcao/module_hcontainer/atom_pair.cpp b/source/module_hamilt_lcao/module_hcontainer/atom_pair.cpp index cef77a34acf..0a6d3162721 100644 --- a/source/module_hamilt_lcao/module_hcontainer/atom_pair.cpp +++ b/source/module_hamilt_lcao/module_hcontainer/atom_pair.cpp @@ -701,7 +701,7 @@ void AtomPair::add_from_matrix(const std::complex* hk, hk_real_pointer = (T*)hk_tmp; hk_imag_pointer = hk_real_pointer+1; BlasConnector::axpy(this->col_size, kphase.real(), hk_real_pointer, 2, hr_tmp, 1); - BlasConnector::axpy(this->col_size, kphase.imag(), hk_imag_pointer, 2, hr_tmp, 1); + BlasConnector::axpy(this->col_size, -kphase.imag(), hk_imag_pointer, 2, hr_tmp, 1); hk_tmp += ld_hk; hr_tmp += this->col_size; } @@ -715,7 +715,7 @@ void AtomPair::add_from_matrix(const std::complex* hk, hk_real_pointer = (T*)hk_tmp; hk_imag_pointer = hk_real_pointer+1; BlasConnector::axpy(this->col_size, kphase.real(), hk_real_pointer, ld_hk_2, hr_tmp, 1); - BlasConnector::axpy(this->col_size, kphase.imag(), hk_imag_pointer, ld_hk_2, hr_tmp, 1); + BlasConnector::axpy(this->col_size, -kphase.imag(), hk_imag_pointer, ld_hk_2, hr_tmp, 1); hk_tmp ++; hr_tmp += this->col_size; } diff --git a/tests/deepks/603_NO_deepks_H2O_multik/result.ref b/tests/deepks/603_NO_deepks_H2O_multik/result.ref index bbe886355d8..14532b6eea9 100644 --- a/tests/deepks/603_NO_deepks_H2O_multik/result.ref +++ b/tests/deepks/603_NO_deepks_H2O_multik/result.ref @@ -1,4 +1,4 @@ -etotref -466.8143178204656 -etotperatomref -155.6047726068 -totalforceref 10.244182 -totaltimeref 7.52 +etotref -466.8999964506086 +etotperatomref -155.6333321502 +totalforceref 10.047085 +totaltimeref 6.22