Skip to content

Commit ca6e24e

Browse files
committed
Merge branch 'develop' of https://github.com/deepmodeling/abacus-develop into mybfgs1
2 parents 73a8136 + 8040017 commit ca6e24e

39 files changed

+2089
-2039
lines changed

docs/advanced/input_files/input-main.md

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1389,7 +1389,7 @@ These variables are used to control the geometry relaxation.
13891389

13901390
- **Type**: Integer
13911391
- **Description**: The maximal number of ionic iteration steps, the minimum value is 1.
1392-
- **Default**: 1
1392+
- **Default**: 1 for SCF, 50 for relax and cell-relax calcualtions
13931393

13941394
### relax_cg_thr
13951395

source/Makefile.Objects

Lines changed: 4 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -255,7 +255,7 @@ OBJS_ESOLVER_LCAO=esolver_ks_lcao.o\
255255
dpks_cal_e_delta_band.o\
256256
set_matrix_grid.o\
257257
lcao_before_scf.o\
258-
lcao_gets.o\
258+
esolver_gets.o\
259259
lcao_others.o\
260260
lcao_init_after_vc.o\
261261

@@ -685,7 +685,9 @@ OBJS_VDW=vdw.o\
685685
vdwd3_parameters.o\
686686
vdwd2.o\
687687
vdwd3.o\
688-
vdwd3_parameters_tab.o
688+
vdwd3_parameters_tab.o\
689+
vdwd3_autoset_xcname.o\
690+
vdwd3_autoset_xcparam.o
689691

690692
OBJS_DFTU=dftu.o\
691693
dftu_force.o\

source/driver_run.cpp

Lines changed: 4 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -67,13 +67,16 @@ void Driver::driver_run() {
6767
Relax_Driver rl_driver;
6868
rl_driver.relax_driver(p_esolver);
6969
}
70+
else if (cal_type == "get_S")
71+
{
72+
p_esolver->runner(0, GlobalC::ucell);
73+
}
7074
else
7175
{
7276
//! supported "other" functions:
7377
//! get_pchg(LCAO),
7478
//! test_memory(PW,LCAO),
7579
//! test_neighbour(LCAO),
76-
//! get_S(LCAO),
7780
//! gen_bessel(PW), et al.
7881
const int istep = 0;
7982
p_esolver->others(istep);

source/module_basis/module_nao/two_center_bundle.cpp

Lines changed: 13 additions & 15 deletions
Original file line numberDiff line numberDiff line change
@@ -63,7 +63,7 @@ void TwoCenterBundle::tabulate()
6363
{
6464
ModuleBase::SphericalBesselTransformer sbt(true);
6565
orb_->set_transformer(sbt);
66-
beta_->set_transformer(sbt);
66+
if (beta_) { beta_->set_transformer(sbt); }
6767
if (alpha_) {
6868
alpha_->set_transformer(sbt);
6969
}
@@ -75,22 +75,17 @@ void TwoCenterBundle::tabulate()
7575
// build two-center integration tables
7676
//================================================================
7777
// set up a universal radial grid
78-
double rmax = std::max(orb_->rcut_max(), beta_->rcut_max());
79-
if (alpha_) {
80-
rmax = std::max(rmax, alpha_->rcut_max());
81-
}
78+
double rmax = orb_->rcut_max();
79+
if (beta_) { rmax = std::max(rmax, beta_->rcut_max()); }
80+
if (alpha_) { rmax = std::max(rmax, alpha_->rcut_max()); }
8281
double dr = 0.01;
8382
double cutoff = 2.0 * rmax;
8483
int nr = static_cast<int>(rmax / dr) + 1;
8584

8685
orb_->set_uniform_grid(true, nr, cutoff, 'i', true);
87-
beta_->set_uniform_grid(true, nr, cutoff, 'i', true);
88-
if (alpha_) {
89-
alpha_->set_uniform_grid(true, nr, cutoff, 'i', true);
90-
}
91-
if (orb_onsite_) {
92-
orb_onsite_->set_uniform_grid(true, nr, cutoff, 'i', true);
93-
}
86+
if (beta_) { beta_->set_uniform_grid(true, nr, cutoff, 'i', true); }
87+
if (alpha_) { alpha_->set_uniform_grid(true, nr, cutoff, 'i', true);}
88+
if (orb_onsite_) { orb_onsite_->set_uniform_grid(true, nr, cutoff, 'i', true);}
9489

9590
// build TwoCenterIntegrator objects
9691
kinetic_orb = std::unique_ptr<TwoCenterIntegrator>(new TwoCenterIntegrator);
@@ -101,9 +96,12 @@ void TwoCenterBundle::tabulate()
10196
overlap_orb->tabulate(*orb_, *orb_, 'S', nr, cutoff);
10297
ModuleBase::Memory::record("TwoCenterTable: Overlap", overlap_orb->table_memory());
10398

104-
overlap_orb_beta = std::unique_ptr<TwoCenterIntegrator>(new TwoCenterIntegrator);
105-
overlap_orb_beta->tabulate(*orb_, *beta_, 'S', nr, cutoff);
106-
ModuleBase::Memory::record("TwoCenterTable: Nonlocal", overlap_orb_beta->table_memory());
99+
if (beta_)
100+
{
101+
overlap_orb_beta = std::unique_ptr<TwoCenterIntegrator>(new TwoCenterIntegrator);
102+
overlap_orb_beta->tabulate(*orb_, *beta_, 'S', nr, cutoff);
103+
ModuleBase::Memory::record("TwoCenterTable: Nonlocal", overlap_orb_beta->table_memory());
104+
}
107105

108106
if (alpha_)
109107
{

source/module_cell/read_atoms.cpp

Lines changed: 18 additions & 50 deletions
Original file line numberDiff line numberDiff line change
@@ -1125,6 +1125,9 @@ void UnitCell::check_dtau() {
11251125

11261126
void UnitCell::read_orb_file(int it, std::string &orb_file, std::ofstream &ofs_running, Atom* atom)
11271127
{
1128+
// the maximum L is 7, according to the basissetexchange https://www.basissetexchange.org/
1129+
// there is no orbitals with L>7 presently
1130+
const std::string spectrum = "SPDFGHIK";
11281131
std::ifstream ifs(orb_file.c_str(), std::ios::in); // pengfei 2014-10-13
11291132
// mohan add return 2021-04-26
11301133
if (!ifs)
@@ -1133,71 +1136,36 @@ void UnitCell::read_orb_file(int it, std::string &orb_file, std::ofstream &ofs_r
11331136
std::cout << " orbital file: " << orb_file << std::endl;
11341137
ModuleBase::WARNING_QUIT("read_orb_file","ABACUS Cannot find the ORBITAL file (basis sets)");
11351138
}
1136-
char word[80];
1139+
std::string word;
11371140
atom->nw = 0;
1138-
int L =0;
11391141
while (ifs.good())
11401142
{
11411143
ifs >> word;
1142-
if (strcmp("Element", word) == 0) // pengfei Li 16-2-29
1144+
if (word == "Element") // pengfei Li 16-2-29
11431145
{
11441146
ModuleBase::GlobalFunc::READ_VALUE(ifs, atom->label_orb);
11451147
}
1146-
if (strcmp("Lmax", word) == 0)
1148+
if (word == "Lmax")
11471149
{
11481150
ModuleBase::GlobalFunc::READ_VALUE(ifs, atom->nwl);
11491151
delete[] atom->l_nchi;
1150-
atom->l_nchi = new int[ atom->nwl+1];
1152+
atom->l_nchi = new int[atom->nwl+1];
11511153
}
1152-
assert(atom->nwl<10);
1153-
if (strcmp("Cutoff(a.u.)", word) == 0) // pengfei Li 16-2-29
1154+
// assert(atom->nwl<10); // cannot understand why restrict the maximum value of atom->nwl
1155+
if (word == "Cutoff(a.u.)") // pengfei Li 16-2-29
11541156
{
11551157
ModuleBase::GlobalFunc::READ_VALUE(ifs, atom->Rcut);
11561158
}
1157-
if (strcmp("Sorbital-->", word) == 0)
1158-
{
1159-
ModuleBase::GlobalFunc::READ_VALUE(ifs, atom->l_nchi[L]);
1160-
atom->nw += (2*L + 1) * atom->l_nchi[L];
1161-
std::stringstream ss;
1162-
ss << "L=" << L << ", number of zeta";
1163-
ModuleBase::GlobalFunc::OUT(ofs_running,ss.str(),atom->l_nchi[L]);
1164-
L++;
1165-
}
1166-
if (strcmp("Porbital-->", word) == 0)
1167-
{
1168-
ModuleBase::GlobalFunc::READ_VALUE(ifs, atom->l_nchi[L]);
1169-
atom->nw += (2*L + 1) * atom->l_nchi[L];
1170-
std::stringstream ss;
1171-
ss << "L=" << L << ", number of zeta";
1172-
ModuleBase::GlobalFunc::OUT(ofs_running,ss.str(),atom->l_nchi[L]);
1173-
L++;
1174-
}
1175-
if (strcmp("Dorbital-->", word) == 0)
1176-
{
1177-
ModuleBase::GlobalFunc::READ_VALUE(ifs, atom->l_nchi[L]);
1178-
atom->nw += (2*L + 1) * atom->l_nchi[L];
1179-
std::stringstream ss;
1180-
ss << "L=" << L << ", number of zeta";
1181-
ModuleBase::GlobalFunc::OUT(ofs_running,ss.str(),atom->l_nchi[L]);
1182-
L++;
1183-
}
1184-
if (strcmp("Forbital-->", word) == 0)
1185-
{
1186-
ModuleBase::GlobalFunc::READ_VALUE(ifs, atom->l_nchi[L]);
1187-
atom->nw += (2*L + 1) * atom->l_nchi[L];
1188-
std::stringstream ss;
1189-
ss << "L=" << L << ", number of zeta";
1190-
ModuleBase::GlobalFunc::OUT(ofs_running,ss.str(),atom->l_nchi[L]);
1191-
L++;
1192-
}
1193-
if (strcmp("Gorbital-->", word) == 0)
1159+
for (int i = 0; i < spectrum.size(); i++)
11941160
{
1195-
ModuleBase::GlobalFunc::READ_VALUE(ifs, atom->l_nchi[L]);
1196-
atom->nw += (2*L + 1) * atom->l_nchi[L];
1197-
std::stringstream ss;
1198-
ss << "L=" << L << ", number of zeta";
1199-
ModuleBase::GlobalFunc::OUT(ofs_running,ss.str(),atom->l_nchi[L]);
1200-
L++;
1161+
if (word == spectrum.substr(i, 1) + "orbital-->")
1162+
{
1163+
ModuleBase::GlobalFunc::READ_VALUE(ifs, atom->l_nchi[i]);
1164+
atom->nw += (2*i + 1) * atom->l_nchi[i];
1165+
std::stringstream ss;
1166+
ss << "L=" << i << ", number of zeta";
1167+
ModuleBase::GlobalFunc::OUT(ofs_running,ss.str(),atom->l_nchi[i]);
1168+
}
12011169
}
12021170
}
12031171
ifs.close();

source/module_esolver/CMakeLists.txt

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -21,7 +21,7 @@ if(ENABLE_LCAO)
2121
dpks_cal_e_delta_band.cpp
2222
set_matrix_grid.cpp
2323
lcao_before_scf.cpp
24-
lcao_gets.cpp
24+
esolver_gets.cpp
2525
lcao_others.cpp
2626
lcao_init_after_vc.cpp
2727
)

source/module_esolver/esolver.cpp

Lines changed: 35 additions & 13 deletions
Original file line numberDiff line numberDiff line change
@@ -5,9 +5,10 @@
55
#include "module_base/module_device/device.h"
66
#include "module_parameter/parameter.h"
77
#ifdef __LCAO
8-
#include "esolver_ks_lcaopw.h"
8+
#include "esolver_gets.h"
99
#include "esolver_ks_lcao.h"
1010
#include "esolver_ks_lcao_tddft.h"
11+
#include "esolver_ks_lcaopw.h"
1112
#include "module_lr/esolver_lrtd_lcao.h"
1213
extern "C"
1314
{
@@ -188,18 +189,39 @@ ESolver* init_esolver(const Input_para& inp, UnitCell& ucell)
188189
{
189190
if (PARAM.globalv.gamma_only_local)
190191
{
191-
return new ESolver_KS_LCAO<double, double>();
192-
}
193-
else if (PARAM.inp.nspin < 4)
194-
{
195-
return new ESolver_KS_LCAO<std::complex<double>, double>();
196-
}
197-
else
198-
{
199-
return new ESolver_KS_LCAO<std::complex<double>, std::complex<double>>();
200-
}
201-
}
202-
else if (esolver_type == "ksdft_lcao_tddft")
192+
if (PARAM.inp.calculation == "get_S")
193+
{
194+
return new ESolver_GetS<double, double>();
195+
}
196+
else
197+
{
198+
return new ESolver_KS_LCAO<double, double>();
199+
}
200+
}
201+
else if (PARAM.inp.nspin < 4)
202+
{
203+
if (PARAM.inp.calculation == "get_S")
204+
{
205+
return new ESolver_GetS<std::complex<double>, double>();
206+
}
207+
else
208+
{
209+
return new ESolver_KS_LCAO<std::complex<double>, double>();
210+
}
211+
}
212+
else
213+
{
214+
if (PARAM.inp.calculation == "get_S")
215+
{
216+
return new ESolver_GetS<std::complex<double>, std::complex<double>>();
217+
}
218+
else
219+
{
220+
return new ESolver_KS_LCAO<std::complex<double>, std::complex<double>>();
221+
}
222+
}
223+
}
224+
else if (esolver_type == "ksdft_lcao_tddft")
203225
{
204226
return new ESolver_KS_LCAO_TDDFT();
205227
}

0 commit comments

Comments
 (0)