Error while using bisection command

[lesheng]

Dear Qbox developer,

I am trying to run a FPMD simulation with the PBE0 functional using the recursive subspace bisection algorithm. However, it keeps giving me the following error

"MPI_RECV : Invalid buffer pointer
Aborting Program!
Abort: MPI_Abort() code: 1, rank 4, Aborting program !: at line 97 in file mpid_init.c"

when I use the bisection command and the "blHF" "btHF" variables

Code: Select all

bisection 2 2 2 0.02
set btHF 0.02
set blHF 2 2 2

Could you give me some suggestions on this please?

Thanks,

[fgygi]

HI Lesheng,
Could please post your complete input file (or add it as attachment with extension ".i" )?
Also, which version of Qbox are you running?
Thanks
Francois

[lesheng]

Hi Francois,

I am using 1.63.5, here is my input file:

Code: Select all

load BN_Li.qb.xml
reset_vcm
set atoms_dyn MD
set xc PBE0
set alpha_PBE0 0.25
bisection 2 2 2 0.02
set btHF 0.02
set blHF 2 2 2
set dt 20
set nempty 36
set thermostat SCALING
set th_time 100
set th_temp 295
set scf_tol 1.E-8
run 2 100 2
save BN_Li.md_test.xml

I was trying to read this saved wave function file BN_Li.qb.xml and run a FPMD simulation, and it generated the above error.

Then I tried to first get the converged face function using this bisection algorithm by using the following input:

Code: Select all

set cell 14.235305872  0  0     7.117652936 12.328136516  0    0  0  32.125341806
species boron /p/lscratche/li44/BN_Li/PP/B_HSCV_PBE-1.0.xml
species nitrogen /p/lscratche/li44/BN_Li/PP/N_HSCV_PBE-1.0.xml
species lithium /p/lscratche/li44/BN_Li/PP/Li_HSCV_PBE-1.0.xml
atom B1 boron 0.00641352 0.00434673 -0.114416
atom B2 boron 4.74927 0.0232024 0.0478804
atom B3 boron 9.49084 0.00531087 -0.118418
atom N1 nitrogen 2.37662 1.37271 0.0282404
atom N2 nitrogen 7.12087 1.37298 0.0229718
atom N3 nitrogen 11.8669 1.37134 -0.136053
atom B4 boron 2.39404 4.1024 0.0478262
atom B5 boron 7.10568 4.10309 0.046615
atom B6 boron 11.8662 4.08375 -0.0559615
atom N4 nitrogen 4.74897 5.48128 0.0224294
atom N5 nitrogen 9.50574 5.48884 0.137381
atom N6 nitrogen 14.2289 5.48917 0.140786
atom B7 boron 4.74946 8.21746 -0.118873
atom B8 boron 9.46917 8.23534 -0.0568062
atom B9 boron 14.2651 8.23639 -0.0631715
atom N7 nitrogen 7.1206 9.59217 -0.136521
atom N8 nitrogen 11.8676 9.579 0.141165
atom N9 nitrogen 16.6133 9.59229 -0.140056
atom Li1 lithium 11.9147 6.88693 3.61816
set xc PBE0
set alpha_PBE0 0.25
bisection 2 2 2 0.02
set btHF 0.02
set blHF 2 2 2
set ecut 40
set ecutprec 10
set wf_dyn PSDA
set wf_diag T
randomize_wf
set nempty 36
set scf_tol 1.E-8
set net_charge 1
run 0 200 2
save BN_Li.qb_wf.xml

By using this above input, the code just hang there while the out put is like this:

...

Code: Select all

[qbox] <cmd>atom B1 boron 0.00641352 0.00434673 -0.114416</cmd>
[qbox] <cmd>atom B2 boron 4.74927 0.0232024 0.0478804</cmd>
[qbox] <cmd>atom B3 boron 9.49084 0.00531087 -0.118418</cmd>
[qbox] <cmd>atom N1 nitrogen 2.37662 1.37271 0.0282404</cmd>
[qbox] <cmd>atom N2 nitrogen 7.12087 1.37298 0.0229718</cmd>
[qbox] <cmd>atom N3 nitrogen 11.8669 1.37134 -0.136053</cmd>
[qbox] <cmd>atom B4 boron 2.39404 4.1024 0.0478262</cmd>
[qbox] <cmd>atom B5 boron 7.10568 4.10309 0.046615</cmd>
[qbox] <cmd>atom B6 boron 11.8662 4.08375 -0.0559615</cmd>
[qbox] <cmd>atom N4 nitrogen 4.74897 5.48128 0.0224294</cmd>
[qbox] <cmd>atom N5 nitrogen 9.50574 5.48884 0.137381</cmd>
[qbox] <cmd>atom N6 nitrogen 14.2289 5.48917 0.140786</cmd>
[qbox] <cmd>atom B7 boron 4.74946 8.21746 -0.118873</cmd>
[qbox] <cmd>atom B8 boron 9.46917 8.23534 -0.0568062</cmd>
[qbox] <cmd>atom B9 boron 14.2651 8.23639 -0.0631715</cmd>
[qbox] <cmd>atom N7 nitrogen 7.1206 9.59217 -0.136521</cmd>
[qbox] <cmd>atom N8 nitrogen 11.8676 9.579 0.141165</cmd>
[qbox] <cmd>atom N9 nitrogen 16.6133 9.59229 -0.140056</cmd>
[qbox] <cmd>atom Li1 lithium 11.9147 6.88693 3.61816</cmd>
[qbox] <cmd>set xc PBE0</cmd>
[qbox] <cmd>set alpha_PBE0 0.25</cmd>
[qbox] <cmd>bisection 2 2 2 0.02</cmd>

It just stopped there...

I attached the input and output files in this .tar file.

Thanks,

HI Lesheng,
Could please post your complete input file (or add it as attachment with extension ".i" )?
Also, which version of Qbox are you running?
Thanks
Francois

[fgygi]

Hi Lesheng,
The error comes from the fact that the bisection command is used when the wave function size is zero. At that point in the script, the plane wave energy cutoff has not yet been defined, so the wave function size is zero. Of course Qbox should do a better job and print an error message rather than crash. We will file a bug report and fix that.

Regarding the use of the bisection algorithm, it appears that bisection will not be very effective on such a small system, and it is not likely to bring a significant speedup of the PBE0 calculation.

For a PBE0 calculation of the ground state, I have the following advice:

• The HSCV potentials are not the best choice we have today. I recommend using the more recent ONCV pseudopotentials available at http://www.quantum-simulation.org/poten ... /index.htm. The ONCV potentials are more accurate and have been tested against all-electron calculations of solids (see M. Schlipf and F. Gygi, Computer Physics Communications 196, 36 (2015).http://dx.doi.org/10.1016/j.cpc.2015.05.011).
• A plane wave cutoff of 60 Ry is advisable for accurate calculations involving the N atom, although you may want to test accuracy with different Ecut values.
• The default value of alpha_PBE0 is 0.25, so there is no need to set it if the normal PBE0 approximation is used.
• Using the JD algorithm for wave function optimization leads to faster convergence.
• The ecutprec variable can be left to its default value (zero) which leads to automatic preconditioning.
• Using the -atomic_density option with the run command help converging faster by starting from a sum of atomic densities.
• When computing a PBE0 ground state, it is usually more efficient to first compute the PBE ground state and then switch to PBE0 and run a few additional iterations. This can all be done in the same script, as you can see below.

I suggest the following script for the PBE0 ground state calculation:

Code: Select all

set cell 14.235305872  0  0   7.117652936 12.328136516  0    0  0  32.125341806
species boron B_ONCV_PBE-1.0.xml
species nitrogen N_ONCV_PBE-1.0.xml
species lithium Li_ONCV_PBE-1.0.xml
atom B1 boron 0.00641352 0.00434673 -0.114416
atom B2 boron 4.74927 0.0232024 0.0478804
atom B3 boron 9.49084 0.00531087 -0.118418
atom N1 nitrogen 2.37662 1.37271 0.0282404
atom N2 nitrogen 7.12087 1.37298 0.0229718
atom N3 nitrogen 11.8669 1.37134 -0.136053
atom B4 boron 2.39404 4.1024 0.0478262
atom B5 boron 7.10568 4.10309 0.046615
atom B6 boron 11.8662 4.08375 -0.0559615
atom N4 nitrogen 4.74897 5.48128 0.0224294
atom N5 nitrogen 9.50574 5.48884 0.137381
atom N6 nitrogen 14.2289 5.48917 0.140786
atom B7 boron 4.74946 8.21746 -0.118873
atom B8 boron 9.46917 8.23534 -0.0568062
atom B9 boron 14.2651 8.23639 -0.0631715
atom N7 nitrogen 7.1206 9.59217 -0.136521
atom N8 nitrogen 11.8676 9.579 0.141165
atom N9 nitrogen 16.6133 9.59229 -0.140056
atom Li1 lithium 11.9147 6.88693 3.61816
set ecut 60
set xc PBE
set wf_dyn JD
set nempty 36
randomize_wf
set scf_tol 1.E-8
set net_charge 1
run -atomic_density 0 200 2
set xc PBE0
run 0 30 2
save gspbe0.xml

Since this system appears to have a sizeable band gap, you may consider running your MD simulation using occupied states only. This makes the calculation substantially faster. Of course, it only works if you can safely assume that the band gap will not close during the simulation.

For a calculation of the PBE0 ground state using occupied states only, I would suggest the following script:

Code: Select all

set cell 14.235305872  0  0   7.117652936 12.328136516  0    0  0  32.125341806
species boron B_ONCV_PBE-1.0.xml
species nitrogen N_ONCV_PBE-1.0.xml
species lithium Li_ONCV_PBE-1.0.xml
atom B1 boron 0.00641352 0.00434673 -0.114416
atom B2 boron 4.74927 0.0232024 0.0478804
atom B3 boron 9.49084 0.00531087 -0.118418
atom N1 nitrogen 2.37662 1.37271 0.0282404
atom N2 nitrogen 7.12087 1.37298 0.0229718
atom N3 nitrogen 11.8669 1.37134 -0.136053
atom B4 boron 2.39404 4.1024 0.0478262
atom B5 boron 7.10568 4.10309 0.046615
atom B6 boron 11.8662 4.08375 -0.0559615
atom N4 nitrogen 4.74897 5.48128 0.0224294
atom N5 nitrogen 9.50574 5.48884 0.137381
atom N6 nitrogen 14.2289 5.48917 0.140786
atom B7 boron 4.74946 8.21746 -0.118873
atom B8 boron 9.46917 8.23534 -0.0568062
atom B9 boron 14.2651 8.23639 -0.0631715
atom N7 nitrogen 7.1206 9.59217 -0.136521
atom N8 nitrogen 11.8676 9.579 0.141165
atom N9 nitrogen 16.6133 9.59229 -0.140056
atom Li1 lithium 11.9147 6.88693 3.61816
set ecut 60
set xc PBE
set wf_dyn JD
randomize_wf
set scf_tol 1.E-8
set net_charge 1
set wf_diag T
run -atomic_density 0 200 2
set xc PBE0
run 0 30 2
save gspbe0_occ.xml

I will post below other suggestions for the MD simulation scripts.

Francois

[fgygi]

For MD simulations using the PBE0 functional, I suggest the following script:

Code: Select all

load gspbe0_occ.xml
set xc PBE0
set wf_dyn JD
set atoms_dyn MD
set dt 20
randomize_v 300
set scf_tol 1.e-6
run 100 10
save mdpbe0_occ.xml

Note:

• This script uses the restart file gspbe0_occ.xml saved at the end of the ground state calculation discussed above.
• The MD simulation only includes occupied states.
• An initial random velocity is given to the atoms, corresponding to an initial temperature of approximately 300K. The temperature is expected to decrease to about half that value after equilibration.

A plot of the Kohn-Sham energy <etotal> compared to the constant of the motion <econst> can be generated using the econste.plt script provided in the util directory:

<etotal> and <econst> during a PBE0 MD simulation.

mdpbe0_occ_econste.jpg (30.12 KiB) Viewed 24207 times

[lesheng]

Hi Francois,

Thanks for you suggestions.

Actually, I was trying to use this BN and Li structure as a test system to run the FPMD simulation with the bisection command. The most recent study we had worked on with Qbox was this BN-Li interface, so I just use it instead of a large system. But your comments does help me to understand how to run the job more efficiently.

As respect to your comment on the error I had,

The error comes from the fact that the bisection command is used when the wave function size is zero. At that point in the script, the plane wave energy cutoff has not yet been defined, so the wave function size is zero.

by saying this, do you mean I can't use this bisection algorithm at this point? If we neglect the fact that this test system is too small to let the bisection works efficiently, do you have any suggestions so I could let this job start running?

Thanks,
Lesheng

[fgygi]

Lesheng,
You can use bisection for this problem using the following script:

Code: Select all

set cell 14.235305872  0  0   7.117652936 12.328136516  0    0  0  32.125341806
species boron B_ONCV_PBE-1.0.xml
species nitrogen N_ONCV_PBE-1.0.xml
species lithium Li_ONCV_PBE-1.0.xml
atom B1 boron 0.00641352 0.00434673 -0.114416
atom B2 boron 4.74927 0.0232024 0.0478804
atom B3 boron 9.49084 0.00531087 -0.118418
atom N1 nitrogen 2.37662 1.37271 0.0282404
atom N2 nitrogen 7.12087 1.37298 0.0229718
atom N3 nitrogen 11.8669 1.37134 -0.136053
atom B4 boron 2.39404 4.1024 0.0478262
atom B5 boron 7.10568 4.10309 0.046615
atom B6 boron 11.8662 4.08375 -0.0559615
atom N4 nitrogen 4.74897 5.48128 0.0224294
atom N5 nitrogen 9.50574 5.48884 0.137381
atom N6 nitrogen 14.2289 5.48917 0.140786
atom B7 boron 4.74946 8.21746 -0.118873
atom B8 boron 9.46917 8.23534 -0.0568062
atom B9 boron 14.2651 8.23639 -0.0631715
atom N7 nitrogen 7.1206 9.59217 -0.136521
atom N8 nitrogen 11.8676 9.579 0.141165
atom N9 nitrogen 16.6133 9.59229 -0.140056
atom Li1 lithium 11.9147 6.88693 3.61816
set ecut 60
set xc PBE
set wf_dyn JD
randomize_wf
set scf_tol 1.E-8
set net_charge 1
set wf_diag T
run -atomic_density 0 200 2
set xc PBE0
set btHF 0.02
set blHF 2 2 2
run 0 30 2
save gspbe0b_occ.xml

It works, but you don't gain much in this case. The code will print this information at each iteration:

Code: Select all

 ExchangeOperator: bisection size: ispin=0: 0.2571790541
 ExchangeOperator: pair fraction:  ispin=0: 0.9160739687
 ExchangeOperator: bisection time: 0.061 s

and this indicates that the bisection algorithm is only able to compress the wave functions by a factor of 0.25 for this tolerance (btHF=0.02). In turn, this translates into a reduction in the number of overlapping pairs by 0.92, so it is only accelerating the calculation of the exchange operator by 10%.

This comes from the fact that this system is small, and essentially two-dimensional, so there is not much freedom to create localized states (only 37 states total). On three-dimensional systems such as 64 water molecules, bisection brings a speedup of 3.5-4.0. On larger systems, it really starts paying off.

[lesheng]

Hi Francois,

Thanks for your reply and suggestion, it works at this time.

I really appreciate your help.

Best,
Lesheng

[lesheng]

Hi Francois,

Thanks again for your previous suggestions and help.

I have read your original PRL paper (PRL, 2009, 102, 166406) about the recursive subspace bisection algorithm and I have a question about this algorithm.

After the linear transformation that localizes wave functions on the arbitrary parts of the basis set, is it true that the wave functions generated by performing bisection algorithm (with PBE0 calculations) are always "smaller" (spatially more localized) than the wave functions obtained by ordinary PBE0 calculations (just like MLWF)?

Thanks,
Lesheng

Lesheng,
You can use bisection for this problem using the following script:

Code: Select all

set cell 14.235305872  0  0   7.117652936 12.328136516  0    0  0  32.125341806
species boron B_ONCV_PBE-1.0.xml
species nitrogen N_ONCV_PBE-1.0.xml
species lithium Li_ONCV_PBE-1.0.xml
atom B1 boron 0.00641352 0.00434673 -0.114416
atom B2 boron 4.74927 0.0232024 0.0478804
atom B3 boron 9.49084 0.00531087 -0.118418
atom N1 nitrogen 2.37662 1.37271 0.0282404
atom N2 nitrogen 7.12087 1.37298 0.0229718
atom N3 nitrogen 11.8669 1.37134 -0.136053
atom B4 boron 2.39404 4.1024 0.0478262
atom B5 boron 7.10568 4.10309 0.046615
atom B6 boron 11.8662 4.08375 -0.0559615
atom N4 nitrogen 4.74897 5.48128 0.0224294
atom N5 nitrogen 9.50574 5.48884 0.137381
atom N6 nitrogen 14.2289 5.48917 0.140786
atom B7 boron 4.74946 8.21746 -0.118873
atom B8 boron 9.46917 8.23534 -0.0568062
atom B9 boron 14.2651 8.23639 -0.0631715
atom N7 nitrogen 7.1206 9.59217 -0.136521
atom N8 nitrogen 11.8676 9.579 0.141165
atom N9 nitrogen 16.6133 9.59229 -0.140056
atom Li1 lithium 11.9147 6.88693 3.61816
set ecut 60
set xc PBE
set wf_dyn JD
randomize_wf
set scf_tol 1.E-8
set net_charge 1
set wf_diag T
run -atomic_density 0 200 2
set xc PBE0
set btHF 0.02
set blHF 2 2 2
run 0 30 2
save gspbe0b_occ.xml

It works, but you don't gain much in this case. The code will print this information at each iteration:

Code: Select all

 ExchangeOperator: bisection size: ispin=0: 0.2571790541
 ExchangeOperator: pair fraction:  ispin=0: 0.9160739687
 ExchangeOperator: bisection time: 0.061 s

and this indicates that the bisection algorithm is only able to compress the wave functions by a factor of 0.25 for this tolerance (btHF=0.02). In turn, this translates into a reduction in the number of overlapping pairs by 0.92, so it is only accelerating the calculation of the exchange operator by 10%.

This comes from the fact that this system is small, and essentially two-dimensional, so there is not much freedom to create localized states (only 37 states total). On three-dimensional systems such as 64 water molecules, bisection brings a speedup of 3.5-4.0. On larger systems, it really starts paying off.