the derivative of the cell-periodic part of the orbitals w.r.t. k

[xiaoming_wang]

Hello,

As explained in the wiki wiki/index.php/LPEAD, there are three methods implemented in vasp to calculate the derivative of the cell-periodic part of the orbitals w.r.t k. One is the sum-over-states method used in LOPTICS=T. The second one is the DFPT method used in LEPSILON=T. Another method is PEAD method. Now, I'd like to evaluate the derivative of the cell-periodic part of the high lying empty orbitals w.r.t k. The PEAD method is not suitable since it only deals with occupied states. Then, for the other two methods, which one is more accurate for this purpose?

Best,

Xiaoming

[henrique_miranda]

Hi,

I think what you are referring to the WAVEDER file:
wiki/index.php/WAVEDER

The relevant information is here:
wiki/index.php/LPEAD

In short, you can use:

Code: Select all

LOPTICS=.TRUE.

to get the derivatives of the WF w.r.t. k using the Sternheimer equation and analytic derivatives of the hamiltonian and overlap operators or

Code: Select all

LOPTICS=.TRUE.
LPEAD=.TRUE.

to use the finite differences approach.
To use one or the other often depends on whether the analytic derivatives of the Hamiltonian are implemented or not for the level of theory you are using.
In some cases like for example self-consistent GW, only the finite differences approach is possible because the analytic derivatives of the hamiltonian and overlap are not implemented.

Also, in the limit of very high k-point sampling the approaches are in principle equivalent, however, in practice there are numeric differences for typically used k-point densities (i.e. k-point densities that would give an accurate value for the total energy).

Kind regards

[xiaoming_wang]

Thanks for your information.

In short, you can use:

Code: Select all

LOPTICS=.TRUE.

to get the derivatives of the WF w.r.t. k using the Sternheimer equation and analytic derivatives of the hamiltonian and overlap operators or

It seems LOPTICS=.TRUE. doesn't use the Sternheimer equation but the first equation (perturbation theory) in this page wiki/index.php/LPEAD. The Sternheimer equation (the second equation in the same page) is used when LEPSILON=.TRUE. is set.

For the derivatives of the WF w.r.t k only (I mean d u_nk / dk not the whole dipole matrix), either the perturbation theory or the LPEAD method might only valid for the occupied states. So is the Sternheimer equation valid for the empty states?

Best,
Xiaoming

[henrique_miranda]

It seems LOPTICS=.TRUE. doesn't use the Sternheimer equation but the first equation (perturbation theory) in this page wiki/index.php/LPEAD. The Sternheimer equation (the second equation in the same page) is used when LEPSILON=.TRUE. is set.

Yes, what I wrote is wrong.
What I meant to say is that one needs to know the analytic derivatives of the Hamiltonian that are also some of the ingredients you need for the Sternheimer equation.
With LOPTICS=.TRUE. VASP does not solve the Sternheimer equation!
I should not have mentioned it.

My point was:

• method 1 : LOPTICS=.TRUE. you obtain the derivatives of the wavefunction w.r.t k using analytic derivatives

• method 2 : LOPTICS=.TRUE.; LPEAD=.TRUE. you obtain the derivatives of the wavefunction w.r.t k using finite differences

For the derivatives of the WF w.r.t k only (I mean d u_nk / dk not the whole dipole matrix), either the perturbation theory or the LPEAD method might only valid for the occupied states. So is the Sternheimer equation valid for the empty states?

Either using method 1 or 2 you shoud get the derivatives of the WFs w.r.t k both for occupied and empty states.
These are indeed ingredients you need for GW and BSE, more precisely for the computation of chi and epsilon (the microscopic dielectric function).
Where did you read that PEAD is only valid for occupied states?

[xiaoming_wang]

Hi Henrique,

Thanks for your reply.

My point was:

• method 1 : LOPTICS=.TRUE. you obtain the derivatives of the wavefunction w.r.t k using analytic derivatives

• method 2 : LOPTICS=.TRUE.; LPEAD=.TRUE. you obtain the derivatives of the wavefunction w.r.t k using finite differences

Either using method 1 or 2 you shoud get the derivatives of the WFs w.r.t k both for occupied and empty states.
These are indeed ingredients you need for GW and BSE, more precisely for the computation of chi and epsilon (the microscopic dielectric function).
Where did you read that PEAD is only valid for occupied states?

It seems that you are talking about the WAVEDER file, which contains basically < psi_n| d/dk| psi_m >, but I'm interested in the derivative d/dk |psi_m>.

For method 2, the PEAD method, the overlap matrix S is defined on the occupied manifold, which can be seen from the pead.F file (PEAD_DPSI_DK subroutine), thus the derivative is only calculated for the occupied states. This doesn't restrict the calculation of the WAVEDER only within the occupied states because the Hermitian nature of the matrix can be employed to evaluate the empty orbital derivatives

< psi_occupied| d/dk| psi_empty > = (< psi_empty| d/dk| psi_occupied >)*

For method 1, the equation 1 in wiki/index.php/LPEAD is in principle exact only for a complete set of conduction bands. The summation would converge very slowly for plane wave basis.

Therefore, it seems that the Sternheimer equation should be the correct way to accurately calculate the k derivative of the empty states (if only d/dk |psi_m> is needed).

[henrique_miranda]

Sorry for my delayed answer.
I think I understand what you mean.

Using the < psi_n| d/dk| psi_m > matrix elements + the GS wavefunctions | psi_m > you should be able to evaluate d/dk |psi_m> for all bands using the sum-over-states approach (Eq. 1 from wiki/index.php/LPEAD).
There you have the problem that the convergence might be slow with respect to how the number of states you sum over.

Using the Sternheimer equation you solve the equations self-consistently directly for d/dk |psi_m> so you overcome the problem of summing over states.
These WFs are computed in linear_response.F after the comment:

Code: Select all

!=======================================================================
! determine G [H, r] |phi> = d/ dk | phi(k)>
!=======================================================================

and stored in the RPHI and RPHI_CPROJ arrays.
Perhaps in the meantime you've figured this out yourself.

PS: For PEAD the overlap is defined for the full number of bands, the sum over m (Eq. 3 in wiki/index.php/LPEAD) runs only over occupied states.