Page 1 of 1

What is the best CUTBox?

Posted: Fri Apr 27, 2018 5:42 am
by wufeng
Hi all,
For low dimensional systems there is an options CUTBox to apply the box coulomb cutoff for faster convergence.
From wiki and tutorial, for 2D monolayer BN the CUTBox is selected to be (c - 1), e.g. CUTBox z = 32Bohr for a 33Bohr cell along z-axis.

I wonder how to select CUTBox for a 2D system with thickness, for example a bilayer BN. Assume the 2D system is 4 Bohr thick and the cell z is 30Bohr, what is best CUTBox?
From my tests the BSE eps fluctuates a lot from CUTBox 15 to 25 (30 cell - 4 thickness - 1 padding ) Bohr. As all of these are indeed larger than 2*thickness, I don't quite understand why results is sensitive to CUTBox and what is the CUTBox should be used here.
Thanks!



Best,
Feng

Re: What is the best CUTBox?

Posted: Fri Apr 27, 2018 6:32 am
by Daniele Varsano
Dear Feng,
actually, you need to consider a factor 2 in the box size, as it is built as -L/2 to L/2. A cutbox of 15 describes an interaction up to 7au and passing from 15 to 25 you are cutting interaction along with your system you do not want to.
The recipe is to consider a value just a bit smaller of your cell size: this means essentially an interaction 1/|r-r'| up to a distance that is half of the cell and having enough vacuum this ensure a Coulomb interaction inside your box avoiding spurious interaction between replica.
One important note is that the box-shaped coulomb potential build involves coulomb integral over the Brillouin zone: this means that the RIM has to be considered when building the coulomb cutoff potential.
So, something like:

Code: Select all

rim_cut                        # [R RIM CUT] Coulomb potential
RandQpts= 1000000              # [RIM] Number of random q-points in the BZ
RandGvec= 123          RL      # [RIM] Coulomb interaction RS components
CUTGeo= "bpx z"               # [CUT] Coulomb Cutoff geometry: box/cylinder/sphere/ws X/Y/Z/XY..
% CUTBox
 0.00     | 0.00     | 29.00     |        # [CUT] [au] Box sides
%
Very soon, probably in the next release, it is planned to release a new cutoff expression based on a cut on the Wigner-Seitz cell that will avoid setting all these parameters and it is appropriate for both 0D, 1D, and 2D systems.

Best,

Daniele

Re: What is the best CUTBox?

Posted: Fri Apr 27, 2018 3:50 pm
by wufeng
Daniele Varsano wrote:Dear Feng,
actually, you need to consider a factor 2 in the box size, as it is built as -L/2 to L/2. A cutbox of 15 describes an interaction up to 7au and passing from 15 to 25 you are cutting interaction along with your system you do not want to.
The recipe is to consider a value just a bit smaller of your cell size: this means essentially an interaction 1/|r-r'| up to a distance that is half of the cell and having enough vacuum this ensure a Coulomb interaction inside your box avoiding spurious interaction between replica.
One important note is that the box-shaped coulomb potential build involves coulomb integral over the Brillouin zone: this means that the RIM has to be considered when building the coulomb cutoff potential.
So, something like:

Code: Select all

rim_cut                        # [R RIM CUT] Coulomb potential
RandQpts= 1000000              # [RIM] Number of random q-points in the BZ
RandGvec= 123          RL      # [RIM] Coulomb interaction RS components
CUTGeo= "bpx z"               # [CUT] Coulomb Cutoff geometry: box/cylinder/sphere/ws X/Y/Z/XY..
% CUTBox
 0.00     | 0.00     | 29.00     |        # [CUT] [au] Box sides
%
Very soon, probably in the next release, it is planned to release a new cutoff expression based on a cut on the Wigner-Seitz cell that will avoid setting all these parameters and it is appropriate for both 0D, 1D, and 2D systems.

Best,

Daniele

Dear Daniele,
Thanks for your information. The idea is my system is only 4au thick, so I don't expect cut 15-25au (interaction 7.5 au to 12.5 au) make a large difference. Is it correct?

For RandOpts and RandGvec I used at 1000000 and 10 RL. What is the suggested value that?


Best,
Feng

Re: What is the best CUTBox?

Posted: Fri Apr 27, 2018 4:39 pm
by Daniele Varsano
Dear Feng,
The idea is my system is only 4au thick, so I don't expect cut 15-25au (interaction 7.5 au to 12.5 au) make a large difference. Is it correct?
I cannot say, depends a bit on the "spill out" of you charge I would say that 25 is safer. Note that small oscillation can appear.
For RandOpts and RandGvec I used at 1000000 and 10 RL. What is the suggested value that?
! Million qpts is usally enough to have accurate integrals, RandGvec it really depends on the dimensions/anysotrpy of your BZ.
In the report you can find:

Code: Select all

 Summary of Coulomb integrals for non-metallic bands |Q|[au] RIM/Bare:
It gives you the ratio between coulomb integrals calculate by Monte Carlo methods and simply as \int V(qi)dqi = V(qi)omega(i) where omega is the small Bz volume around each qi. You can have a look at the ratio for the q at the edge of your Bz. When the Coulomb potential become sufficiently smooth, then this is a good approximation (ratio of the order 0.99) otherwise it is recommendable to use MC also outside the BZ. IN order to evaluate how many G vectors to treat with MC integrals a possible test is a simple Sigma_x calculation for just one state varying RandGvec and see where it converges.

Best,
Daniele

Re: What is the best CUTBox?

Posted: Fri Apr 27, 2018 5:00 pm
by wufeng
Daniele Varsano wrote:Dear Feng,
The idea is my system is only 4au thick, so I don't expect cut 15-25au (interaction 7.5 au to 12.5 au) make a large difference. Is it correct?
I cannot say, depends a bit on the "spill out" of you charge I would say that 25 is safer. Note that small oscillation can appear.
For RandOpts and RandGvec I used at 1000000 and 10 RL. What is the suggested value that?
! Million qpts is usally enough to have accurate integrals, RandGvec it really depends on the dimensions/anysotrpy of your BZ.
In the report you can find:

Code: Select all

 Summary of Coulomb integrals for non-metallic bands |Q|[au] RIM/Bare:
It gives you the ratio between coulomb integrals calculate by Monte Carlo methods and simply as \int V(qi)dqi = V(qi)omega(i) where omega is the small Bz volume around each qi. You can have a look at the ratio for the q at the edge of your Bz. When the Coulomb potential become sufficiently smooth, then this is a good approximation (ratio of the order 0.99) otherwise it is recommendable to use MC also outside the BZ. IN order to evaluate how many G vectors to treat with MC integrals a possible test is a simple Sigma_x calculation for just one state varying RandGvec and see where it converges.

Best,
Daniele
Hi Daniele,

Just to be sure, so this is what I am thinking about: for a slab with 4 au thickness , in a cell with z=30 au , the maximum CutBox is 30 - 4 = 26 au, and to avoid the "spill out" of the charge, to include intralayer interaction CutBox should be as large as possible, and to exclude interlayer interaction CutBox should be a bit smaller than 26. Is the statement above correct?

This is the output in my case:

Code: Select all

  Summary of Coulomb integrals for non-metallic bands |Q|[au] RIM/Bare:

  Q [1]:0.1000E-40.8755 * Q [7]:  0.07190  0.71071
  Q [2]:  0.07192  0.71159 * Q [8]: 0.101696 0.819092
  Q [13]: 0.143799 0.888265 * Q [3]: 0.143841 0.888692
  Q [14]: 0.160782 0.908328 * Q [9]: 0.160810 0.908499
  Q [15]: 0.203393 0.939286 * Q [19]: 0.215699 0.944609
  Q [4]: 0.215762 0.944865 * Q [20]: 0.227373 0.949887
  Q [10]: 0.227427 0.950043 * Q [21]: 0.259261 0.960812
  Q [16]: 0.259290 0.960871 * Q [25]: 0.287598 0.967440
  Q [5]: 0.287683 0.967617 * Q [26]: 0.296455 0.969303
  Q [11]: 0.296532 0.969430 * Q [22]: 0.305089 0.971103
  Q [27]: 0.321564 0.973737 * Q [17]: 0.321620 0.973809
  Q [31]: 0.359498 0.978635 * Q [28]: 0.359536 0.978771
  Q [23]: 0.359565 0.978799 * Q [6]: 0.359604 0.978648
  Q [32]: 0.366621 0.979440 * Q [12]: 0.366721 0.979435
  Q [33]: 0.387207 0.981501 * Q [18]: 0.387289 0.981477
  Q [29]: 0.406785 0.983214 * Q [34]: 0.419276 0.984122
  Q [24]: 0.419334 0.984082 * Q [37]: 0.431397 0.984891
  Q [38]: 0.437352 0.985292 * Q [39]: 0.454746 0.986360
  Q [35]: 0.460435 0.986722 * Q [30]: 0.460464 0.986673
  Q [40]: 0.482345 0.987820 * Q [43]: 0.503297 0.988800
  Q [44]: 0.508410 0.989021 * Q [36]: 0.508482 0.988953
  Q [41]: 0.518522 0.989390 * Q [45]: 0.523448 0.989623
  Q [46]: 0.547596 0.990481 * Q [42]: 0.561621 0.990838
  Q [47]: 0.579715 0.991459 * Q [48]: 0.618565 0.992398
Last one approaches 0.99 and is it enough?

Thanks very much!

Best,
Feng

Re: What is the best CUTBox?

Posted: Fri Apr 27, 2018 5:16 pm
by Daniele Varsano
Dear Feng
Just to be sure, so this is what I am thinking about: for a slab with 4 au thickness , in a cell with z=30 au , the maximum CutBox is 30 - 4 = 26 au, and to avoid the "spill out" of the charge, to include intralayer interaction CutBox should be as large as possible, and to exclude interlayer interaction CutBox should be a bit smaller than 26. Is the statement above correct?
Yes, it is correct but I would not set CutBox > L which means Lcell/2.


[---------<--- S--->--------][---------<--- S--->--------]
<------------L------------>



L is the side of your cell. S is the dimension of your system. L-S vacuum
You want a Lcut > S+d where d takes into account the spilling of the charge, as you want to have full interaction inside your supercell.
and you want also Lcut < L-(S+d) in order to ensure your system do not interact with your replica.
If you have at least L>2(S+d) an Lcut ~L/2 should fulfill both conditions. Note that there is factor 2 so Lcut~L is ok, but it is recommendable to not set it exactly as the lenght of the cell.
Last one approaches 0.99 and is it enough?
Well, it is not perfect and it indicates that some MC integration is needed outside the first BZ. Probably the 10G you consdered are enough. You can consider to use more as they are calculated once for all as they are stored in a database (ndb.RIM) and not calculated anymore.

Best,
Daniele

Re: What is the best CUTBox?

Posted: Fri Apr 27, 2018 7:34 pm
by wufeng
Daniele Varsano wrote: Yes, it is correct but I would not set CutBox > L which means Lcell/2.


[---------<--- S--->--------][---------<--- S--->--------]
<------------L------------>



L is the side of your cell. S is the dimension of your system. L-S vacuum
You want a Lcut > S+d where d takes into account the spilling of the charge, as you want to have full interaction inside your supercell.
and you want also Lcut < L-(S+d) in order to ensure your system do not interact with your replica.
If you have at least L>2(S+d) an Lcut ~L/2 should fulfill both conditions. Note that there is factor 2 so Lcut~L is ok, but it is recommendable to not set it exactly as the lenght of the cell.
Sorry but I didn't get what the 'I would not set CutBox > L which means Lcell/2' exactly means as L and Lcell are the same in discussion below.

Based on discussion below, Lcut < Lcell-(S+d), And CutBox = 2Lcut, so CutBox < 2Lcell - 2(S+d) is enough, which means CutBox can reach ~2Lcell for a thin layer.

but it is recommendable to set CutBox a bit smaller than Lcell. Is that what you mean?

Best,
Feng

Re: What is the best CUTBox?

Posted: Sat Apr 28, 2018 12:09 am
by Daniele Varsano
Dear Feng,
sorry I did not explain myself in a good way.
Let's assume the spill out contained in S for simplicity.
We have S<Lcut<Lcell-S
This means you need to have at least enough vacuum as S (Lcell=2S). Then Lcut = Lcell/2 satisfy the condition.
Cutbox is defined as 2Lcut in the code so the suggestion is to set Cutbox slightly smaller than L.
Consider also that the potential in planewave it will have the periodicity of the cell (no k points in the direction you want isolate) so does not make
sense to consider Lcut > Lcell/2.
Hope now it is clearer,
Best,
Daniele

PS: I know this is a bit messy, I really need to include the Wigner Seitz based cutoff asap.

Re: What is the best CUTBox?

Posted: Sat Apr 28, 2018 7:34 am
by wufeng
Daniele Varsano wrote:Dear Feng,
sorry I did not explain myself in a good way.
Let's assume the spill out contained in S for simplicity.
We have S<Lcut<Lcell-S
This means you need to have at least enough vacuum as S (Lcell=2S). Then Lcut = Lcell/2 satisfy the condition.
Cutbox is defined as 2Lcut in the code so the suggestion is to set Cutbox slightly smaller than L.
Consider also that the potential in planewave it will have the periodicity of the cell (no k points in the direction you want isolate) so does not make
sense to consider Lcut > Lcell/2.
Hope now it is clearer,
Best,
Daniele

PS: I know this is a bit messy, I really need to include the Wigner Seitz based cutoff asap.
Dear Daniele,
Thanks for your patience, it is clear to me now.

Also, I found

Code: Select all

     Summary of Coulomb integrals for non-metallic bands |Q|[au] RIM/Bare:
changes slowly respect to RandGvec.

For a system with

With 1M RandQpts and 10 RandGvec:

Code: Select all

  Q [46]: 0.547596 0.989350 * Q [42]: 0.561621 0.989990
  Q [47]: 0.579715 0.990345 * Q [48]: 0.618565 0.991313
With 1M RandQpts and 1000 RandGvec:

Code: Select all

  Q [46]: 0.547596 0.991227 * Q [42]: 0.561621 0.991632
  Q [47]: 0.579715 0.992198 * Q [48]: 0.618565 0.993215
I am not sure if it is the expected behavior.

Best,
Feng

Re: What is the best CUTBox?

Posted: Sat Apr 28, 2018 7:44 am
by Daniele Varsano
Dear Feng ,
the differences you observe using different RandGvec is most probably due to the statistic error associated with the value of the integral evaluated in MC method.
Different runs, different seeds. The discrepancy (ie the statistical error) should be reduced by using more random points (1/sqrt(Nq)).

Daniele