Skip to content

Commit 264a307

Browse files
PeizeLinPeizeLin
committed
Feature: add ABFS_ORBITAL and ABFS_JLES_ORBITAL in Exx_Opt_Orb
1 parent 5411878 commit 264a307

File tree

10 files changed

+135
-74
lines changed

10 files changed

+135
-74
lines changed

‎source/source_cell/bcast_cell.cpp‎

Lines changed: 8 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -112,8 +112,14 @@ namespace unitcell
112112

113113
#ifdef __EXX
114114
ModuleBase::bcast_data_cereal(GlobalC::exx_info.info_ri.files_abfs,
115-
MPI_COMM_WORLD,
116-
0);
115+
MPI_COMM_WORLD,
116+
0);
117+
ModuleBase::bcast_data_cereal(GlobalC::exx_info.info_opt_abfs.files_abfs,
118+
MPI_COMM_WORLD,
119+
0);
120+
ModuleBase::bcast_data_cereal(GlobalC::exx_info.info_opt_abfs.files_jles,
121+
MPI_COMM_WORLD,
122+
0);
117123
#endif
118124
return;
119125
#endif

‎source/source_cell/read_atom_species.cpp‎

Lines changed: 10 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -112,6 +112,16 @@ bool read_atom_species(std::ifstream& ifa,
112112
std::string ofile;
113113
ifa >> ofile;
114114
GlobalC::exx_info.info_ri.files_abfs.push_back(ofile);
115+
GlobalC::exx_info.info_opt_abfs.files_abfs.push_back(ofile);
116+
}
117+
}
118+
if( ModuleBase::GlobalFunc::SCAN_LINE_BEGIN(ifa, "ABFS_JLES_ORBITAL") )
119+
{
120+
for(int i=0; i<ntype; i++)
121+
{
122+
std::string ofile;
123+
ifa >> ofile;
124+
GlobalC::exx_info.info_opt_abfs.files_jles.push_back(ofile);
115125
}
116126
}
117127
}

‎source/source_hamilt/module_xc/exx_info.h‎

Lines changed: 5 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -82,6 +82,11 @@ struct Exx_Info
8282
int abfs_Lmax = 0; // tmp
8383
double ecut_exx = 60;
8484
double tolerence = 1E-12;
85+
std::vector<std::string> files_jles;
86+
87+
double pca_threshold = 0;
88+
std::vector<std::string> files_abfs;
89+
8590
double kmesh_times = 4;
8691
};
8792
Exx_Info_Opt_ABFs info_opt_abfs;

‎source/source_io/input_conv.cpp‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -473,6 +473,7 @@ void Input_Conv::Convert()
473473
GlobalC::exx_info.info_ri.ccp_rmesh_times = std::stod(PARAM.inp.exx_ccp_rmesh_times);
474474
GlobalC::exx_info.info_ri.exx_symmetry_realspace = PARAM.inp.exx_symmetry_realspace;
475475

476+
GlobalC::exx_info.info_opt_abfs.pca_threshold = PARAM.inp.exx_pca_threshold;
476477
GlobalC::exx_info.info_opt_abfs.abfs_Lmax = PARAM.inp.exx_opt_orb_lmax;
477478
GlobalC::exx_info.info_opt_abfs.ecut_exx = PARAM.inp.exx_opt_orb_ecut;
478479
GlobalC::exx_info.info_opt_abfs.tolerence = PARAM.inp.exx_opt_orb_tolerence;

‎source/source_io/read_input_item_exx_dftu.cpp‎

Lines changed: 0 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -337,12 +337,6 @@ void ReadInput::item_exx()
337337
Input_Item item("exx_opt_orb_lmax");
338338
item.annotation = "the maximum l of the spherical Bessel functions for opt ABFs";
339339
read_sync_int(input.exx_opt_orb_lmax);
340-
item.check_value = [](const Input_Item& item, const Parameter& para) {
341-
if (para.input.exx_opt_orb_lmax < 0)
342-
{
343-
ModuleBase::WARNING_QUIT("ReadInput", "exx_opt_orb_lmax must >= 0");
344-
}
345-
};
346340
this->add_item(item);
347341
}
348342
{

‎source/source_lcao/module_ri/exx_abfs-io.cpp‎

Lines changed: 17 additions & 13 deletions
Original file line numberDiff line numberDiff line change
@@ -36,7 +36,7 @@ std::vector<std::vector<std::vector<Numerical_Orbital_Lm>>> Exx_Abfs::IO::constr
3636
const double kmesh_times )
3737
{
3838
std::vector<std::vector<std::vector<Numerical_Orbital_Lm>>>
39-
&&abfs = construct_abfs( orbs, files_abfs, kmesh_times );
39+
abfs = construct_abfs( orbs, files_abfs, kmesh_times );
4040

4141
assert( abfs.size() == abfs_pre.size() );
4242
for( size_t T=0; T!=abfs.size(); ++T )
@@ -73,7 +73,7 @@ std::vector<std::vector<Numerical_Orbital_Lm>> Exx_Abfs::IO::construct_abfs_T(
7373

7474
std::ifstream ifs( file_name.c_str() );
7575
if(!ifs)
76-
throw std::runtime_error(" Can't find the abfs ORBITAL file.");
76+
throw std::runtime_error(" Can't find the abfs ORBITAL file " + file_name);
7777

7878
while( ifs.good() )
7979
{
@@ -130,10 +130,22 @@ std::vector<std::vector<Numerical_Orbital_Lm>> Exx_Abfs::IO::construct_abfs_T(
130130
{
131131
ModuleBase::GlobalFunc::READ_VALUE( ifs, N_size[8] );
132132
}
133+
else if ( "Lorbital-->"==word )
134+
{
135+
ModuleBase::GlobalFunc::READ_VALUE( ifs, N_size[9] );
136+
}
137+
else if ( "Morbital-->"==word )
138+
{
139+
ModuleBase::GlobalFunc::READ_VALUE( ifs, N_size[10] );
140+
}
141+
else if ( "Norbital-->"==word )
142+
{
143+
ModuleBase::GlobalFunc::READ_VALUE( ifs, N_size[11] );
144+
}
133145
else if ( "END"==word )
134146
{
135147
break;
136-
}
148+
}
137149
}
138150

139151
ModuleBase::CHECK_NAME(ifs, "Mesh");
@@ -169,19 +181,11 @@ std::vector<std::vector<Numerical_Orbital_Lm>> Exx_Abfs::IO::construct_abfs_T(
169181
----------------------*/
170182
for( size_t L=0; L<=L_size; ++L )
171183
if( N_size.find(L) == N_size.end() )
172-
{
173-
std::stringstream ss;
174-
ss<<"Can't find N of L="<<L<<" in "<<file_name;
175-
throw std::domain_error(ss.str());
176-
}
184+
{ throw std::domain_error("Can't find N of L="+std::to_string(L)+" in "+file_name); }
177185
for( size_t L=0; L<=L_size; ++L )
178186
for( size_t N=0; N!=N_size[L]; ++N )
179187
if( psis.find(L)==psis.end() || psis[L].find(N)==psis[L].end() )
180-
{
181-
std::stringstream ss;
182-
ss<<"Can't find abf of L="<<L<<" T="<<T<<" in "<<file_name;
183-
throw std::domain_error(ss.str());
184-
}
188+
{ throw std::domain_error("Can't find abf of L="+std::to_string(L)+" T="+std::to_string(T)+" in "+file_name); }
185189

186190

187191
/*----------------------

‎source/source_lcao/module_ri/exx_abfs-jle.cpp‎

Lines changed: 6 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -5,6 +5,7 @@
55
#include "../../source_cell/unitcell.h"
66
#include "../../source_base/mathzone.h"
77
#include "../../source_base/math_sphbes.h" // mohan add 2021-05-06
8+
#include "source_base/tool_title.h"
89

910
std::vector<std::vector<std::vector<Numerical_Orbital_Lm>>>
1011
Exx_Abfs::Jle::init_jle(
@@ -13,15 +14,19 @@ Exx_Abfs::Jle::init_jle(
1314
const UnitCell& ucell,
1415
const LCAO_Orbitals& orb)
1516
{
17+
ModuleBase::TITLE("Exx_Abfs::Jle","init_jle");
1618
std::vector<std::vector<std::vector<Numerical_Orbital_Lm>>> jle( ucell.ntype );
1719

1820
for(int T=0; T<ucell.ntype; ++T)
1921
{
20-
jle[T].resize( info.abfs_Lmax+1 );
22+
if(info.abfs_Lmax+1>0)
23+
{ jle[T].resize( info.abfs_Lmax+1 ); }
2124
for(int L=0; L<=info.abfs_Lmax; ++L)
2225
{
2326
const size_t ecut_number
2427
= static_cast<size_t>( std::sqrt( info.ecut_exx ) * orb.Phi[T].getRcut() / ModuleBase::PI ); // Rydberg Unit.
28+
if(ecut_number<=0)
29+
{ continue; }
2530

2631
jle[T][L].resize( ecut_number );
2732

‎source/source_lcao/module_ri/exx_opt_orb-print.cpp‎

Lines changed: 39 additions & 40 deletions
Original file line numberDiff line numberDiff line change
@@ -1,18 +1,21 @@
11
#include "exx_opt_orb.h"
22
#include "../../source_pw/module_pwdft/global.h"
33
#include "exx_abfs-jle.h"
4+
#include "source_base/tool_title.h"
45
#include <iomanip>
56

67
void Exx_Opt_Orb::print_matrix(
78
const Exx_Info::Exx_Info_Opt_ABFs &info,
89
const UnitCell& ucell,
9-
const K_Vectors &kv,
10-
const std::string& file_name,
10+
const K_Vectors &kv,
11+
const int Lmax,
12+
const std::vector<std::size_t> &ecut_number,
13+
const std::string &file_name,
1114
const std::vector<RI::Tensor<double>> &matrix_Q,
1215
const std::vector<std::vector<RI::Tensor<double>>> &matrix_S,
1316
const RI::Tensor<double> &matrix_V,
14-
const size_t TA, const size_t IA, const size_t TB, const size_t IB,
15-
const std::vector<double>& orb_cutoff,
17+
const std::size_t TA, const std::size_t IA, const std::size_t TB, const std::size_t IB,
18+
const std::vector<double>& orb_cutoff,
1619
const ModuleBase::Element_Basis_Index::Range &range_jles,
1720
const ModuleBase::Element_Basis_Index::IndexLNM &index_jles) const
1821
{
@@ -61,45 +64,40 @@ void Exx_Opt_Orb::print_matrix(
6164
<< ucell.atoms[TB].tau[IB].z << std::endl;
6265
}
6366

64-
// ecutwfc_jlq determine the jlq corresponding to plane wave calculation.
65-
ofs << info.ecut_exx << " ecutwfc" << std::endl; // mohan add 2009-09-08
67+
ofs << info.ecut_exx << " ecutwfc" << std::endl;
6668

6769
// this parameter determine the total number of jlq.
68-
ofs << info.ecut_exx << " ecutwfc_jlq" << std::endl;//mohan modify 2009-09-08
70+
ofs << info.ecut_exx << " ecutwfc_jlq" << std::endl;
6971

7072
if(TA==TB)
7173
{ ofs << orb_cutoff[TA] << " rcut_Jlq" << std::endl; }
7274
else
7375
{ ofs << orb_cutoff[TA] << " " << orb_cutoff[TB] << " rcut_Jlq" << std::endl; }
7476

75-
// mohan add 'smooth' and 'smearing_sigma' 2009-08-28
7677
ofs << 0 << " smooth" << std::endl;
7778
ofs << 0 << " smearing_sigma" << std::endl;
7879

7980
ofs << info.tolerence << " tolerence" << std::endl;
8081

81-
ofs << info.abfs_Lmax << " lmax" << std::endl;
82+
ofs << Lmax << " lmax" << std::endl;
8283

8384
ofs << kv.get_nkstot() << " nks" << std::endl;
8485
assert( matrix_V.shape[0]*matrix_V.shape[1] == matrix_V.shape[2]*matrix_V.shape[3] );
8586
ofs << matrix_V.shape[0]*matrix_V.shape[1] << " nbands" << std::endl;
8687

87-
auto cal_sum_M = [&range_jles](size_t T) -> size_t
88+
auto cal_sum_M = [&range_jles](std::size_t T) -> std::size_t
8889
{
89-
size_t sum_M = 0;
90-
for( size_t L = 0; L!=range_jles[T].size(); ++L )
90+
std::size_t sum_M = 0;
91+
for( std::size_t L = 0; L!=range_jles[T].size(); ++L )
9192
{ sum_M += range_jles[T][L].M; }
9293
return sum_M;
9394
};
94-
const size_t nwfc = (TA==TB && IA==IB) ? cal_sum_M(TA) : cal_sum_M(TA)+cal_sum_M(TB);
95+
const std::size_t nwfc = (TA==TB && IA==IB) ? cal_sum_M(TA) : cal_sum_M(TA)+cal_sum_M(TB);
9596
ofs << nwfc << " nwfc" << std::endl;
9697

97-
const size_t ecut_numberA = static_cast<size_t>( std::sqrt( info.ecut_exx ) * orb_cutoff[TA] / ModuleBase::PI ); // Rydberg Unit
98-
const size_t ecut_numberB = static_cast<size_t>( std::sqrt( info.ecut_exx ) * orb_cutoff[TB] / ModuleBase::PI ); // Rydberg Unit
99-
if(TA==TB)
100-
{ ofs << ecut_numberA << " ne" << std::endl; }
101-
else
102-
{ ofs << ecut_numberA << " " << ecut_numberB << " ne" << std::endl; }
98+
for(const std::size_t ne : ecut_number)
99+
{ ofs << ne << " "; }
100+
ofs << "ne" << std::endl;
103101

104102
ofs << "<WEIGHT_OF_KPOINTS>" << std::endl;
105103
for( int ik=0; ik!=kv.get_nkstot(); ++ik )
@@ -119,18 +117,18 @@ void Exx_Opt_Orb::print_matrix(
119117
// < Psi | jY >
120118
//---------------------
121119
ofs<< "<OVERLAP_Q>" << std::endl;
122-
for( size_t iw0=0; iw0!=matrix_V.shape[0]; ++iw0 )
120+
for( std::size_t iw0=0; iw0!=matrix_V.shape[0]; ++iw0 )
123121
{
124-
for( size_t iw1=0; iw1!=matrix_V.shape[1]; ++iw1 )
122+
for( std::size_t iw1=0; iw1!=matrix_V.shape[1]; ++iw1 )
125123
{
126-
for( size_t iat=0; iat!=matrix_Q.size(); ++iat )
124+
for( std::size_t iat=0; iat!=matrix_Q.size(); ++iat )
127125
{
128-
const size_t it = (iat==0) ? TA : TB;
129-
for( size_t il=0; il!=range_jles[it].size(); ++il )
126+
const std::size_t it = (iat==0) ? TA : TB;
127+
for( std::size_t il=0; il!=range_jles[it].size(); ++il )
130128
{
131-
for( size_t im=0; im!=range_jles[it][il].M; ++im )
129+
for( std::size_t im=0; im!=range_jles[it][il].M; ++im )
132130
{
133-
for( size_t iq=0; iq!=range_jles[it][il].N; ++iq )
131+
for( std::size_t iq=0; iq!=range_jles[it][il].N; ++iq )
134132
{
135133
ofs<<matrix_Q[iat]( iw0, iw1, index_jles[it][il][iq][im] )<<"\t"<<0<<std::endl;
136134
}
@@ -149,23 +147,23 @@ void Exx_Opt_Orb::print_matrix(
149147
// < jY | jY >
150148
//---------------------
151149
ofs<< "<OVERLAP_Sq>" <<std::endl;
152-
for( size_t iat1=0; iat1!=matrix_S.size(); ++iat1 )
150+
for( std::size_t iat1=0; iat1!=matrix_S.size(); ++iat1 )
153151
{
154-
const size_t it1 = (iat1==0) ? TA : TB;
155-
for( size_t il1=0; il1!=range_jles[it1].size(); ++il1 )
152+
const std::size_t it1 = (iat1==0) ? TA : TB;
153+
for( std::size_t il1=0; il1!=range_jles[it1].size(); ++il1 )
156154
{
157-
for( size_t im1=0; im1!=range_jles[it1][il1].M; ++im1 )
155+
for( std::size_t im1=0; im1!=range_jles[it1][il1].M; ++im1 )
158156
{
159-
for( size_t iat2=0; iat2!=matrix_S[iat1].size(); ++iat2 )
157+
for( std::size_t iat2=0; iat2!=matrix_S[iat1].size(); ++iat2 )
160158
{
161-
const size_t it2 = (iat2==0) ? TA : TB;
162-
for( size_t il2=0; il2!=range_jles[it2].size(); ++il2 )
159+
const std::size_t it2 = (iat2==0) ? TA : TB;
160+
for( std::size_t il2=0; il2!=range_jles[it2].size(); ++il2 )
163161
{
164-
for( size_t im2=0; im2!=range_jles[it2][il2].M; ++im2 )
162+
for( std::size_t im2=0; im2!=range_jles[it2][il2].M; ++im2 )
165163
{
166-
for( size_t iq1=0; iq1!=range_jles[it1][il1].N; ++iq1 )
164+
for( std::size_t iq1=0; iq1!=range_jles[it1][il1].N; ++iq1 )
167165
{
168-
for( size_t iq2=0; iq2!=range_jles[it2][il2].N; ++iq2 )
166+
for( std::size_t iq2=0; iq2!=range_jles[it2][il2].N; ++iq2 )
169167
{
170168
ofs<<matrix_S[iat1][iat2]( index_jles[it1][il1][iq1][im1], index_jles[it2][il2][iq2][im2] )*scale<<"\t"<<0<<std::endl;
171169
}
@@ -186,13 +184,13 @@ void Exx_Opt_Orb::print_matrix(
186184
// < Psi | Psi >
187185
//---------------------
188186
ofs << "<OVERLAP_V>" << std::endl;
189-
for( size_t iw0=0; iw0!=matrix_V.shape[0]; ++iw0 )
187+
for( std::size_t iw0=0; iw0!=matrix_V.shape[0]; ++iw0 )
190188
{
191-
for( size_t iw1=0; iw1!=matrix_V.shape[1]; ++iw1 )
189+
for( std::size_t iw1=0; iw1!=matrix_V.shape[1]; ++iw1 )
192190
{
193-
for( size_t iw2=0; iw2!=matrix_V.shape[2]; ++iw2 )
191+
for( std::size_t iw2=0; iw2!=matrix_V.shape[2]; ++iw2 )
194192
{
195-
for( size_t iw3=0; iw3!=matrix_V.shape[3]; ++iw3 )
193+
for( std::size_t iw3=0; iw3!=matrix_V.shape[3]; ++iw3 )
196194
{
197195
ofs<<matrix_V(iw0,iw1,iw2,iw3)*scale<<"\t";
198196
}
@@ -203,6 +201,7 @@ void Exx_Opt_Orb::print_matrix(
203201
ofs << "</OVERLAP_V>" << std::endl << std::endl;
204202
};
205203

204+
ModuleBase::TITLE("Exx_Opt_Orb","print_matrix");
206205
std::ofstream ofs(file_name+"_"+std::to_string(TA)+"_"+std::to_string(IA)+"_"+std::to_string(TB)+"_"+std::to_string(IB));
207206
print_header(ofs);
208207
ofs<<std::setprecision(15);

0 commit comments

Comments
 (0)