Skip to content

Commit 289f36c

Browse files
committed
Fix: Remove nonzero tails beyond spherical Bessel projector cutoffs
1 parent 70f7ed6 commit 289f36c

90 files changed

Lines changed: 6647 additions & 6871 deletions

File tree

Some content is hidden

Large Commits have some content hidden by default. Use the searchbox below for content that may be hidden.

examples/21_deepks/02_lcao_H2O/jle.orb

Lines changed: 7 additions & 13 deletions
Original file line numberDiff line numberDiff line change
@@ -8,7 +8,7 @@ Number of Dorbitals--> 2
88
---------------------------------------------------------------------------
99
SUMMARY END
1010

11-
Mesh 205
11+
Mesh 201
1212
dr 0.01
1313
Type L N
1414
0 0 0
@@ -62,8 +62,7 @@ dr 0.01
6262
6.345247331791e-02 5.791188607588e-02 5.241540711817e-02 4.696361790712e-02
6363
4.155709094895e-02 3.619638971392e-02 3.088206855820e-02 2.561467264749e-02
6464
2.039473788254e-02 1.522279082639e-02 1.009934863353e-02 5.024918980770e-03
65-
8.824636488425e-14 -4.974919786756e-03 -9.899361531700e-03 -1.477285612199e-02
66-
-1.959494423991e-02
65+
8.824636488425e-14
6766
Type L N
6867
0 0 1
6968
1.000000000000e+00 9.998355147105e-01 9.993421562398e-01 9.985202167122e-01
@@ -116,8 +115,7 @@ dr 0.01
116115
-6.232855556732e-02 -5.704953906849e-02 -5.177008647810e-02 -4.649509277994e-02
117116
-4.122940087424e-02 -3.597779810401e-02 -3.074501284474e-02 -2.553571116011e-02
118117
-2.035449352599e-02 -1.520589162508e-02 -1.009436521457e-02 -5.024299068919e-03
119-
-2.180811153812e-14 4.974306043314e-03 9.894476794432e-03 1.475645640459e-02
120-
1.955627809356e-02
118+
-2.180811153812e-14
121119
Type L N
122120
0 1 0
123121
0.000000000000e+00 7.488637748270e-03 1.497500757073e-02 2.245684235920e-02
@@ -170,8 +168,7 @@ dr 0.01
170168
6.163241347987e-02 5.629489953858e-02 5.098844044591e-02 4.571475520896e-02
171169
4.047554347042e-02 3.527248491362e-02 3.010723867693e-02 2.498144277778e-02
172170
1.989671354647e-02 1.485464507002e-02 9.856808646282e-03 4.904752248416e-03
173-
8.084377910810e-14 -4.855948338613e-03 -9.661617874115e-03 -1.441555907126e-02
174-
-1.911634823371e-02
171+
8.084377910810e-14
175172
Type L N
176173
0 1 1
177174
0.000000000000e+00 1.287349883354e-02 2.573547475531e-02 3.857441713205e-02
@@ -224,8 +221,7 @@ dr 0.01
224221
-6.113827004858e-02 -5.605889221235e-02 -5.095291481399e-02 -4.582759735266e-02
225222
-4.069015817980e-02 -3.554776565875e-02 -3.040752943992e-02 -2.527649186198e-02
226223
-2.016161948939e-02 -1.506979479638e-02 -1.000780800721e-02 -4.982349102379e-03
227-
-6.437579797176e-15 4.932773078295e-03 9.809627048423e-03 1.462434920303e-02
228-
1.937086424077e-02
224+
-6.437579797176e-15
229225
Type L N
230226
0 2 0
231227
0.000000000000e+00 5.535915211195e-05 2.213972080154e-04 4.979959858820e-04
@@ -278,8 +274,7 @@ dr 0.01
278274
5.992574581769e-02 5.478191562228e-02 4.965613429040e-02 4.455121891443e-02
279275
3.946996773718e-02 3.441515831428e-02 2.938954569414e-02 2.439586061649e-02
280276
1.943680773090e-02 1.451506383641e-02 9.633276143556e-03 4.794060559853e-03
281-
4.131818525851e-15 -4.746357278055e-03 -9.442499310155e-03 -1.408595203483e-02
282-
-1.867428087307e-02
277+
4.131818525851e-15
283278
Type L N
284279
0 2 1
285280
0.000000000000e+00 1.378450216078e-04 5.511357836611e-04 1.239139794773e-03
@@ -332,5 +327,4 @@ dr 0.01
332327
-5.983213508908e-02 -5.496252350283e-02 -5.004030907345e-02 -4.507498543086e-02
333328
-4.007604502571e-02 -3.505296252697e-02 -3.001517833094e-02 -2.497208221092e-02
334329
-1.993299713613e-02 -1.490716328841e-02 -9.903722304457e-03 -4.931701771028e-03
335-
4.773444701199e-14 4.882628890872e-03 9.707589564902e-03 1.446645969794e-02
336-
1.915100382853e-02
330+
4.773444701199e-14

examples/21_deepks/03_lcao_CsPbI3/jle.orb

Lines changed: 46 additions & 91 deletions
Large diffs are not rendered by default.

source/source_basis/module_ao/ORB_read.cpp

Lines changed: 24 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -7,7 +7,9 @@
77

88
#include <algorithm>
99
#include <cassert>
10+
#include <cmath>
1011
#include <cstring> // Peize Lin fix bug about strcmp 2016-08-02
12+
#include <iostream>
1113

1214
//==============================
1315
// Define an object here!
@@ -378,6 +380,7 @@ void LCAO_Orbitals::read_orb_file(std::ofstream& ofs_in, // GlobalV::ofs_running
378380
ModuleBase::TITLE("LCAO_Orbitals", "read_orb_file");
379381
char word[80];
380382
std::string orb_label;
383+
double declared_rcut = -1.0;
381384
if (my_rank == 0)
382385
{
383386
while (ifs.good())
@@ -388,6 +391,15 @@ void LCAO_Orbitals::read_orb_file(std::ofstream& ofs_in, // GlobalV::ofs_running
388391
ifs >> orb_label;
389392
continue;
390393
}
394+
if (std::strcmp(word, "Radius") == 0)
395+
{
396+
ifs >> word;
397+
if (std::strcmp(word, "Cutoff(a.u.)") == 0)
398+
{
399+
ifs >> declared_rcut;
400+
}
401+
continue;
402+
}
391403
if (std::strcmp(word, "Lmax") == 0)
392404
{
393405
ifs >> lmax;
@@ -448,6 +460,18 @@ void LCAO_Orbitals::read_orb_file(std::ofstream& ofs_in, // GlobalV::ofs_running
448460
}
449461
ModuleBase::CHECK_NAME(ifs, "dr");
450462
ifs >> dr;
463+
464+
if (declared_rcut >= 0.0 && meshr_read > 0)
465+
{
466+
const double mesh_rcut = (meshr_read - 1) * dr;
467+
const double tolerance = 1.0e-10 * std::max(1.0, std::abs(declared_rcut));
468+
if (std::abs(mesh_rcut - declared_rcut) > tolerance)
469+
{
470+
std::cout << " WARNING: The orbital file declares a cutoff radius of " << declared_rcut
471+
<< " Bohr, but (Mesh - 1) * dr is " << mesh_rcut
472+
<< " Bohr. The file will be read without modification." << std::endl;
473+
}
474+
}
451475
}
452476

453477
#ifdef __MPI

source/source_basis/module_ao/test/ORB_read_test.cpp

Lines changed: 2 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -286,12 +286,12 @@ TEST_F(LcaoOrbitalsTest, ReadOrbitals) {
286286
EXPECT_EQ(aod.PhiLN(L,N).getType(), 0);
287287
EXPECT_EQ(aod.PhiLN(L,N).getL(), L);
288288
EXPECT_EQ(aod.PhiLN(L,N).getChi(), N);
289-
EXPECT_EQ(aod.PhiLN(L,N).getNr(), 205);
289+
EXPECT_EQ(aod.PhiLN(L,N).getNr(), 201);
290290
EXPECT_EQ(aod.PhiLN(L,N).getNk(), lcao_.kmesh);
291291
EXPECT_EQ(aod.PhiLN(L,N).getDk(), lcao_.dk);
292292
EXPECT_EQ(aod.PhiLN(L,N).getDruniform(), lcao_.dr_uniform);
293293

294-
for (int ir = 0; ir != 205; ++ir) {
294+
for (int ir = 0; ir != 201; ++ir) {
295295
EXPECT_DOUBLE_EQ(aod.PhiLN(L,N).getRab(ir), 0.01);
296296
EXPECT_DOUBLE_EQ(aod.PhiLN(L,N).getRadial(ir), 0.01*ir);
297297
}
@@ -340,5 +340,3 @@ int main(int argc, char **argv)
340340

341341
return result;
342342
}
343-
344-

source/source_basis/module_ao/test/lcao_H2O/jle.orb

Lines changed: 7 additions & 13 deletions
Original file line numberDiff line numberDiff line change
@@ -8,7 +8,7 @@ Number of Dorbitals--> 2
88
---------------------------------------------------------------------------
99
SUMMARY END
1010

11-
Mesh 205
11+
Mesh 201
1212
dr 0.01
1313
Type L N
1414
0 0 0
@@ -62,8 +62,7 @@ dr 0.01
6262
6.345247331791e-02 5.791188607588e-02 5.241540711817e-02 4.696361790712e-02
6363
4.155709094895e-02 3.619638971392e-02 3.088206855820e-02 2.561467264749e-02
6464
2.039473788254e-02 1.522279082639e-02 1.009934863353e-02 5.024918980770e-03
65-
8.824636488425e-14 -4.974919786756e-03 -9.899361531700e-03 -1.477285612199e-02
66-
-1.959494423991e-02
65+
8.824636488425e-14
6766
Type L N
6867
0 0 1
6968
1.000000000000e+00 9.998355147105e-01 9.993421562398e-01 9.985202167122e-01
@@ -116,8 +115,7 @@ dr 0.01
116115
-6.232855556732e-02 -5.704953906849e-02 -5.177008647810e-02 -4.649509277994e-02
117116
-4.122940087424e-02 -3.597779810401e-02 -3.074501284474e-02 -2.553571116011e-02
118117
-2.035449352599e-02 -1.520589162508e-02 -1.009436521457e-02 -5.024299068919e-03
119-
-2.180811153812e-14 4.974306043314e-03 9.894476794432e-03 1.475645640459e-02
120-
1.955627809356e-02
118+
-2.180811153812e-14
121119
Type L N
122120
0 1 0
123121
0.000000000000e+00 7.488637748270e-03 1.497500757073e-02 2.245684235920e-02
@@ -170,8 +168,7 @@ dr 0.01
170168
6.163241347987e-02 5.629489953858e-02 5.098844044591e-02 4.571475520896e-02
171169
4.047554347042e-02 3.527248491362e-02 3.010723867693e-02 2.498144277778e-02
172170
1.989671354647e-02 1.485464507002e-02 9.856808646282e-03 4.904752248416e-03
173-
8.084377910810e-14 -4.855948338613e-03 -9.661617874115e-03 -1.441555907126e-02
174-
-1.911634823371e-02
171+
8.084377910810e-14
175172
Type L N
176173
0 1 1
177174
0.000000000000e+00 1.287349883354e-02 2.573547475531e-02 3.857441713205e-02
@@ -224,8 +221,7 @@ dr 0.01
224221
-6.113827004858e-02 -5.605889221235e-02 -5.095291481399e-02 -4.582759735266e-02
225222
-4.069015817980e-02 -3.554776565875e-02 -3.040752943992e-02 -2.527649186198e-02
226223
-2.016161948939e-02 -1.506979479638e-02 -1.000780800721e-02 -4.982349102379e-03
227-
-6.437579797176e-15 4.932773078295e-03 9.809627048423e-03 1.462434920303e-02
228-
1.937086424077e-02
224+
-6.437579797176e-15
229225
Type L N
230226
0 2 0
231227
0.000000000000e+00 5.535915211195e-05 2.213972080154e-04 4.979959858820e-04
@@ -278,8 +274,7 @@ dr 0.01
278274
5.992574581769e-02 5.478191562228e-02 4.965613429040e-02 4.455121891443e-02
279275
3.946996773718e-02 3.441515831428e-02 2.938954569414e-02 2.439586061649e-02
280276
1.943680773090e-02 1.451506383641e-02 9.633276143556e-03 4.794060559853e-03
281-
4.131818525851e-15 -4.746357278055e-03 -9.442499310155e-03 -1.408595203483e-02
282-
-1.867428087307e-02
277+
4.131818525851e-15
283278
Type L N
284279
0 2 1
285280
0.000000000000e+00 1.378450216078e-04 5.511357836611e-04 1.239139794773e-03
@@ -332,5 +327,4 @@ dr 0.01
332327
-5.983213508908e-02 -5.496252350283e-02 -5.004030907345e-02 -4.507498543086e-02
333328
-4.007604502571e-02 -3.505296252697e-02 -3.001517833094e-02 -2.497208221092e-02
334329
-1.993299713613e-02 -1.490716328841e-02 -9.903722304457e-03 -4.931701771028e-03
335-
4.773444701199e-14 4.882628890872e-03 9.707589564902e-03 1.446645969794e-02
336-
1.915100382853e-02
330+
4.773444701199e-14

source/source_basis/module_nao/atomic_radials.cpp

Lines changed: 22 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -9,6 +9,7 @@
99

1010
#include "source_base/projgen.h"
1111

12+
#include <cmath>
1213
#include <fstream>
1314
#include <iostream>
1415
#include <string>
@@ -119,6 +120,7 @@ void AtomicRadials::read_abacus_orb(std::ifstream& ifs, std::ofstream* ptr_log,
119120
* */
120121
int ngrid = 0; // number of grid points
121122
double dr = 0; // grid spacing
123+
double declared_rcut = -1.0;
122124
std::string tmp;
123125

124126
if (rank == 0)
@@ -143,6 +145,14 @@ void AtomicRadials::read_abacus_orb(std::ifstream& ifs, std::ofstream* ptr_log,
143145
{
144146
ifs >> orb_ecut_;
145147
}
148+
else if (tmp == "Radius")
149+
{
150+
ifs >> tmp;
151+
if (tmp == "Cutoff(a.u.)")
152+
{
153+
ifs >> declared_rcut;
154+
}
155+
}
146156
else if (tmp == "Lmax")
147157
{
148158
ifs >> lmax_;
@@ -167,6 +177,18 @@ void AtomicRadials::read_abacus_orb(std::ifstream& ifs, std::ofstream* ptr_log,
167177
}
168178
}
169179

180+
if (declared_rcut >= 0.0 && ngrid > 0)
181+
{
182+
const double mesh_rcut = (ngrid - 1) * dr;
183+
const double tolerance = 1.0e-10 * std::max(1.0, std::abs(declared_rcut));
184+
if (std::abs(mesh_rcut - declared_rcut) > tolerance)
185+
{
186+
std::cout << " WARNING: The orbital file declares a cutoff radius of " << declared_rcut
187+
<< " Bohr, but (Mesh - 1) * dr is " << mesh_rcut
188+
<< " Bohr. The file will be read without modification." << std::endl;
189+
}
190+
}
191+
170192
/*
171193
* calculate:
172194
*

source/source_io/module_bessel/bessel_basis.cpp

Lines changed: 21 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -5,6 +5,9 @@
55
#include "source_base/math_sphbes.h"
66
#include "source_base/parallel_common.h"
77
#include "source_base/timer.h"
8+
9+
#include <algorithm>
10+
#include <cmath>
811
#include <vector>
912

1013
Bessel_Basis::Bessel_Basis()
@@ -210,6 +213,15 @@ void Bessel_Basis::init_TableOne(
210213
std::vector<double> function(rmesh);
211214
std::vector<double> en(ecut_number);
212215

216+
const int cutoff_intervals = static_cast<int>(std::round(rcut / dr));
217+
const double cutoff_tolerance = 1.0e-10 * std::max(1.0, std::abs(rcut));
218+
const bool cutoff_is_on_grid = std::abs(cutoff_intervals * dr - rcut) <= cutoff_tolerance;
219+
220+
// Preserve the requested radial spacing. When the cutoff is on the grid,
221+
// omit the auxiliary integration points from the orbital file. Otherwise,
222+
// retain the original mesh and write zeros beyond the cutoff.
223+
const int output_rmesh = cutoff_is_on_grid ? cutoff_intervals + 1 : rmesh;
224+
213225
for(int ir=0; ir<rmesh; ir++)
214226
{
215227
r[ir] = static_cast<double>(ir) * dr;
@@ -251,8 +263,8 @@ void Bessel_Basis::init_TableOne(
251263
}
252264
ofs << "---------------------------------------------------------------------------"<< std::endl;
253265
ofs << "SUMMARY END" << std::endl << std::endl;
254-
ofs << std::setiosflags(std::ios::left) << std::setw(28) << "Mesh" << rmesh << std::endl;
255-
ofs << std::setiosflags(std::ios::left) << std::setw(28) << "dr" << dr << std::endl ;
266+
ofs << std::setiosflags(std::ios::left) << std::setw(28) << "Mesh" << output_rmesh << std::endl;
267+
ofs << std::setiosflags(std::ios::left) << std::setw(28) << "dr" << dr << std::endl;
256268
//=========output .orb format=============
257269

258270
// init eigenvalue of Jl
@@ -279,10 +291,14 @@ void Bessel_Basis::init_TableOne(
279291
//=========output .orb format=============
280292
ofs << std::setiosflags(std::ios::right) << std::setw(20) << "Type"<< std::setw(20) << "L" << std::setw(20) << "N" << std::endl;
281293
ofs << std::setiosflags(std::ios::right) << std::setw(20) << "0"<< std::setw(20) << l << std::setw(20) << ie << std::endl;
282-
for (int ir = 0; ir < rmesh; ir++)
294+
for (int ir = 0; ir < output_rmesh; ir++)
283295
{
284-
ofs << std::setiosflags(std::ios::scientific)
285-
<< std::setprecision(12) << jle[ir]<< " "; if ((ir+1) % 4 == 0) ofs << std::endl;
296+
const double output_value = r[ir] > rcut + cutoff_tolerance ? 0.0 : jle[ir];
297+
ofs << std::setiosflags(std::ios::scientific) << std::setprecision(12) << output_value << " ";
298+
if ((ir + 1) % 4 == 0)
299+
{
300+
ofs << std::endl;
301+
}
286302
}
287303
ofs << std::endl;
288304
//=========output .orb format=============
Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1 +1 @@
1-
-0.08058091803
1+
-0.08076909443
Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1 +1 @@
1-
-0.3135463017
1+
-0.313773958
Lines changed: 5 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -1,5 +1,5 @@
1-
-0.001966690223 0.0003863484866 -0.0002775679931
2-
0.0001490387822 0.0001336996716 -0.010866736
3-
-0.009536649342 -0.0001507865891 0.002944084409
4-
0.005665599528 0.00917496669 0.003914146719
5-
0.005688701256 -0.009544228259 0.004286072869
1+
-0.00197050246 0.0003857169325 -0.000286430865
2+
0.000147912611 0.0001343902978 -0.01086809874
3+
-0.009542180958 -0.000151654385 0.002945191164
4+
0.005671544221 0.009179414195 0.003918275884
5+
0.005693226587 -0.00954786704 0.004291062558

0 commit comments

Comments
 (0)