Skip to content

Commit 52b7e75

Browse files
ddhhssddhhss
committed
Add function to get totn(g).
1 parent 744fb1b commit 52b7e75

File tree

5 files changed

+69
-34
lines changed

5 files changed

+69
-34
lines changed

‎source/module_surchem/H_correction_pw.cpp‎

Lines changed: 21 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -26,7 +26,7 @@ ModuleBase::matrix surchem::v_correction(const UnitCell &cell,
2626
complex<double> *Porter_g = new complex<double>[rho_basis->npw];
2727
ModuleBase::GlobalFunc::ZEROS(Porter_g, rho_basis->npw);
2828

29-
GlobalC::UFFT.ToReciSpace(Porter, Porter_g, rho_basis);
29+
rho_basis->real2recip(Porter, Porter_g);
3030

3131
complex<double> *N = new complex<double>[rho_basis->npw];
3232
complex<double> *TOTN = new complex<double>[rho_basis->npw];
@@ -57,7 +57,8 @@ void surchem::add_comp_chg(const UnitCell &cell,
5757
double l,
5858
double center,
5959
complex<double> *NG,
60-
int dim)
60+
int dim,
61+
bool flag)
6162
{
6263
// x dim
6364
double tmp_q = 0.0;
@@ -70,7 +71,10 @@ void surchem::add_comp_chg(const UnitCell &cell,
7071
{
7172
if (ig == rho_basis->ig_gge0)
7273
{
73-
NG[ig] = complex<double>(tmp_q * l / L, 0.0);
74+
if(flag)
75+
{
76+
NG[ig] = complex<double>(tmp_q * l / L, 0.0);
77+
}
7478
continue;
7579
}
7680
double GX = rho_basis->gcar[ig][0];
@@ -94,7 +98,10 @@ void surchem::add_comp_chg(const UnitCell &cell,
9498
{
9599
if (ig == rho_basis->ig_gge0)
96100
{
97-
NG[ig] = complex<double>(tmp_q * l / L, 0.0);
101+
if(flag)
102+
{
103+
NG[ig] = complex<double>(tmp_q * l / L, 0.0);
104+
}
98105
continue;
99106
}
100107
double GX = rho_basis->gcar[ig][0];
@@ -118,7 +125,10 @@ void surchem::add_comp_chg(const UnitCell &cell,
118125
{
119126
if (ig == rho_basis->ig_gge0)
120127
{
121-
NG[ig] = complex<double>(tmp_q * l / L, 0.0);
128+
if(flag)
129+
{
130+
NG[ig] = complex<double>(tmp_q * l / L, 0.0);
131+
}
122132
continue;
123133
}
124134
double GX = rho_basis->gcar[ig][0];
@@ -154,7 +164,7 @@ ModuleBase::matrix surchem::v_compensating(const UnitCell &cell,
154164
complex<double> *Porter_g = new complex<double>[rho_basis->npw];
155165
ModuleBase::GlobalFunc::ZEROS(Porter_g, rho_basis->npw);
156166

157-
GlobalC::UFFT.ToReciSpace(Porter, Porter_g, rho_basis);
167+
rho_basis->real2recip(Porter, Porter_g);
158168

159169
this->Porter_g_anchor = Porter_g[rho_basis->ig_gge0];
160170

@@ -164,7 +174,7 @@ ModuleBase::matrix surchem::v_compensating(const UnitCell &cell,
164174
cal_totn(cell, rho_basis, Porter_g, N, TOTN);
165175

166176
// save TOTN in real space
167-
GlobalC::UFFT.ToRealSpace(TOTN, this->TOTN_real, rho_basis);
177+
rho_basis->recip2real(TOTN, this->TOTN_real);
168178

169179
complex<double> *comp_reci = new complex<double>[rho_basis->npw];
170180
complex<double> *phi_comp_G = new complex<double>[rho_basis->npw];
@@ -174,9 +184,9 @@ ModuleBase::matrix surchem::v_compensating(const UnitCell &cell,
174184
ModuleBase::GlobalFunc::ZEROS(phi_comp_G, rho_basis->npw);
175185
ModuleBase::GlobalFunc::ZEROS(phi_comp_R, rho_basis->nrxx);
176186
// get compensating charge in reci space
177-
add_comp_chg(cell, rho_basis, comp_q, comp_l, comp_center, comp_reci, comp_dim);
187+
add_comp_chg(cell, rho_basis, comp_q, comp_l, comp_center, comp_reci, comp_dim, true);
178188
// save compensating charge in real space
179-
GlobalC::UFFT.ToRealSpace(comp_reci, this->comp_real, rho_basis);
189+
rho_basis->recip2real(comp_reci, this->comp_real);
180190

181191
// test sum of comp_real -> 0
182192
// for (int i = 0; i < rho_basis->nz;i++)
@@ -186,7 +196,7 @@ ModuleBase::matrix surchem::v_compensating(const UnitCell &cell,
186196
// double sum = 0;
187197
// for (int i = 0; i < rho_basis->nxyz; i++)
188198
// {
189-
// sum += comp_real[i];
199+
// sum += TOTN_real[i];
190200
// }
191201
// sum = sum * cell.omega / rho_basis->nxyz;
192202
// cout << "sum:" << sum << endl;
@@ -240,7 +250,7 @@ void surchem::cal_comp_force(ModuleBase::matrix &force_comp, ModulePW::PW_Basis
240250
std::complex<double> *N = new std::complex<double>[rho_basis->npw];
241251
std::complex<double> *phi_comp_G = new complex<double>[rho_basis->npw];
242252
std::complex<double> *vloc_at = new std::complex<double>[rho_basis->npw];
243-
GlobalC::UFFT.ToReciSpace(phi_comp_R, phi_comp_G, rho_basis);
253+
rho_basis->real2recip(phi_comp_R, phi_comp_G);
244254

245255
for (int it = 0; it < GlobalC::ucell.ntype; it++)
246256
{

‎source/module_surchem/corrected_energy.cpp‎

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -41,7 +41,7 @@ void surchem::cal_Acomp(const UnitCell &cell,
4141
ModuleBase::GlobalFunc::ZEROS(phi_comp_R, rho_basis->nrxx);
4242

4343
// part1: comp & comp
44-
GlobalC::UFFT.ToReciSpace(comp_real, comp_reci, rho_basis);
44+
rho_basis->real2recip(comp_real, comp_reci);
4545
for (int ig = 0; ig < rho_basis->npw; ig++)
4646
{
4747
if (rho_basis->gg[ig] >= 1.0e-12) // LiuXh 20180410
@@ -72,7 +72,7 @@ void surchem::cal_Acomp(const UnitCell &cell,
7272
}
7373

7474
// part2: electrons
75-
GlobalC::UFFT.ToRealSpace(phi_comp_G, phi_comp_R, rho_basis);
75+
rho_basis->recip2real(phi_comp_G, phi_comp_R);
7676
for (int ir = 0; ir < rho_basis->nrxx; ir++)
7777
{
7878
Acomp2 += n_elec_R[ir] * phi_comp_R[ir];

‎source/module_surchem/surchem.cpp‎

Lines changed: 22 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -61,4 +61,26 @@ surchem::~surchem()
6161
delete[] epspot;
6262
delete[] comp_real;
6363
delete[] phi_comp_R;
64+
}
65+
66+
void surchem::get_totn_reci(const UnitCell &cell, ModulePW::PW_Basis *rho_basis, complex<double> *totn_reci)
67+
{
68+
double *tmp_totn_real = new double[rho_basis->nrxx];
69+
double *tmp_comp_real = new double[rho_basis->nrxx];
70+
complex<double> *comp_reci = new complex<double>[rho_basis->npw];
71+
ModuleBase::GlobalFunc::ZEROS(tmp_totn_real, rho_basis->nrxx);
72+
ModuleBase::GlobalFunc::ZEROS(tmp_comp_real, rho_basis->nrxx);
73+
ModuleBase::GlobalFunc::ZEROS(comp_reci, rho_basis->npw);
74+
add_comp_chg(cell, rho_basis, comp_q, comp_l, comp_center, comp_reci, comp_dim, false);
75+
rho_basis->recip2real(comp_reci, tmp_comp_real);
76+
77+
for (int ir = 0; ir < rho_basis->nrxx;ir++)
78+
{
79+
tmp_totn_real[ir] = TOTN_real[ir] + tmp_comp_real[ir];
80+
}
81+
82+
rho_basis->real2recip(tmp_totn_real, totn_reci);
83+
delete[] tmp_totn_real;
84+
delete[] tmp_comp_real;
85+
delete[] comp_reci;
6486
}

‎source/module_surchem/surchem.h‎

Lines changed: 18 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -53,7 +53,8 @@ class surchem
5353
double l,
5454
double center,
5555
complex<double> *NG,
56-
int dim);
56+
int dim,
57+
bool flag); // Set value of comp_reci[ig_gge0] when flag is true.
5758

5859
void cal_comp_force(ModuleBase::matrix &force_comp, ModulePW::PW_Basis *rho_basis);
5960

@@ -84,7 +85,10 @@ class surchem
8485

8586
double cal_Acav(const UnitCell &cell, ModulePW::PW_Basis *rho_basis);
8687

87-
void cal_Acomp(const UnitCell &cell, ModulePW::PW_Basis *rho_basis, const double *const *const rho, vector<double> &res);
88+
void cal_Acomp(const UnitCell &cell,
89+
ModulePW::PW_Basis *rho_basis,
90+
const double *const *const rho,
91+
vector<double> &res);
8892

8993
void minimize_cg(const UnitCell &ucell,
9094
ModulePW::PW_Basis *rho_basis,
@@ -104,19 +108,24 @@ class surchem
104108
complex<double> *lp);
105109

106110
ModuleBase::matrix v_correction(const UnitCell &cell,
107-
ModulePW::PW_Basis* rho_basis,
108-
const int &nspin,
109-
const double *const *const rho);
111+
ModulePW::PW_Basis *rho_basis,
112+
const int &nspin,
113+
const double *const *const rho);
110114

111115
ModuleBase::matrix v_compensating(const UnitCell &cell,
112116
ModulePW::PW_Basis *rho_basis,
113117
const int &nspin,
114118
const double *const *const rho);
115119

116-
void test_V_to_N(ModuleBase::matrix &v, const UnitCell &cell, ModulePW::PW_Basis *rho_basis, const double *const *const rho);
117-
118-
void cal_force_sol(const UnitCell &cell, ModulePW::PW_Basis* rho_basis , ModuleBase::matrix& forcesol);
119-
120+
void test_V_to_N(ModuleBase::matrix &v,
121+
const UnitCell &cell,
122+
ModulePW::PW_Basis *rho_basis,
123+
const double *const *const rho);
124+
125+
void cal_force_sol(const UnitCell &cell, ModulePW::PW_Basis *rho_basis, ModuleBase::matrix &forcesol);
126+
127+
void get_totn_reci(const UnitCell &cell, ModulePW::PW_Basis *rho_basis, complex<double> *totn_reci);
128+
120129
private:
121130
};
122131

‎source/src_pw/potential.cpp‎

Lines changed: 6 additions & 12 deletions
Original file line numberDiff line numberDiff line change
@@ -364,22 +364,10 @@ ModuleBase::matrix Potential::v_of_rho(const double *const *const rho_in, const
364364
if(GlobalV::comp_chg)
365365
{
366366
v += GlobalC::solvent_model.v_compensating(GlobalC::ucell, GlobalC::rhopw, GlobalV::NSPIN, rho_in);
367-
vector<double> tmp(3, 0);
368-
GlobalC::solvent_model.cal_Acomp(GlobalC::ucell, GlobalC::rhopw, rho_in, tmp);
369-
ModuleBase::matrix forcecomp(GlobalC::ucell.nat, 3);
370-
GlobalC::solvent_model.cal_comp_force(forcecomp, GlobalC::rhopw);
371-
int p;
372-
// cin >> p;
373367
}
374368
if (GlobalV::imp_sol)
375369
{
376370
v += GlobalC::solvent_model.v_correction(GlobalC::ucell, GlobalC::rhopw, GlobalV::NSPIN, rho_in);
377-
/*
378-
// test energy outside
379-
cout << "energy Outside: " << endl;
380-
GlobalC::solvent_model.cal_Ael(GlobalC::ucell, GlobalC::rhopw);
381-
GlobalC::solvent_model.cal_Acav(GlobalC::ucell, GlobalC::rhopw);
382-
*/
383371
}
384372
}
385373

@@ -391,6 +379,12 @@ ModuleBase::matrix Potential::v_of_rho(const double *const *const rho_in, const
391379
v += Efield::add_efield(GlobalC::ucell, GlobalC::rhopw, GlobalV::NSPIN, rho_in);
392380
}
393381

382+
// test get ntot_reci
383+
// complex<double> *tmpn = new complex<double>[GlobalC::rhopw->npw];
384+
// ModuleBase::GlobalFunc::ZEROS(tmpn, GlobalC::rhopw->npw);
385+
// GlobalC::solvent_model.get_totn_reci(GlobalC::ucell, GlobalC::rhopw, tmpn);
386+
// delete[] tmpn;
387+
394388
ModuleBase::timer::tick("Potential", "v_of_rho");
395389
return v;
396390
} // end subroutine v_of_rho

0 commit comments

Comments
 (0)