PhononWork -> generating too many DFPT calculations

[cpashartis]

Hello,

Is anyone familiar with the PhononWork calculations? From what I understand, and from what I could decipher from the code, there seems to be an issue with generating the correct number of calculations (perhaps at the qpoint steps).

When looking at just qpt = (0,0,0):

1. I was trying a surface of Silicon in the 2x1 reconstruction and when I generate the PhononWork DFPT steps it generates 24 calculations, all of which have the same input...

2. I also had the same issue with crystalline cubic Si. It generated two tasks with the same input for DFPT, and when analyzed, each task had the same phonon frequencies (of course).

[cpashartis]

As an FYI, I noticed the 24 comes from the number of perturbations required on directions and atoms:

The list of irreducible perturbations for this q vector is:
   1)    idir= 1    ipert=   1
   2)    idir= 2    ipert=   1
   3)    idir= 3    ipert=   1
   4)    idir= 1    ipert=   3
   5)    idir= 2    ipert=   3
   6)    idir= 3    ipert=   3
   7)    idir= 1    ipert=   5
   8)    idir= 2    ipert=   5
   9)    idir= 3    ipert=   5
  10)    idir= 1    ipert=   7
  11)    idir= 2    ipert=   7
  12)    idir= 3    ipert=   7
  13)    idir= 1    ipert=   9
  14)    idir= 2    ipert=   9
  15)    idir= 3    ipert=   9
  16)    idir= 1    ipert=  11
  17)    idir= 2    ipert=  11
  18)    idir= 3    ipert=  11
  19)    idir= 1    ipert=  13
  20)    idir= 2    ipert=  13
  21)    idir= 3    ipert=  13
  22)    idir= 1    ipert=  15
  23)    idir= 2    ipert=  15
  24)    idir= 3    ipert=  15

[gmatteo]

I was trying a surface of Silicon in the 2x1 reconstruction and when I generate the PhononWork DFPT steps it generates 24 calculations, all of which have the same input...

These input files solve the DFPT equations for different perturbations.
You should see that each input file has a different value for (idir, ipert, qpt).
This point is discussed here

First of all, AbiPy generates PhononTasks only for the q-points in the irreducible wedge of the Brillouin zone corresponding to ph_ngqpt. Moreover, for a given q-point, only the irreducible atomic perturbations are explicitly computed since the other atomic perturbations can be reconstructed by symmetry

[cpashartis]

Yes so I see, now. It is a different direction for the response functions.

So the odd thing and I guess why I was confused: the input files are the exact same on startup. Then I think after they begin to run they change? I was using diff t0/run.abi t1/run.abi and they were the same until I started t0, which I found odd. (At least I'm pretty sure or I wouldn't have made this post ;)).

As an aside, I may also need the eigenvectors for the phonon frequencies, from what I saw in the code there is no easy connection via abipy to get this, correct me if I'm wrong.

Thanks!

[gmatteo]

Then I think after they begin to run they change

Not really. The input files with the different DFPT perturbations are generated at the beginning when you generate the flow.
This is what I get with e.g.:

vimdiff flow_phonons/w1/t0/run.abi flow_phonons/w1/t1/run.abi

As an aside, I may also need the eigenvectors for the phonon frequencies, from what I saw in the code there is no easy connection via abipy to get this, correct me if I'm wrong.

phonon frequencies and displacement vectors are stored in the PhononBands object
as documented here

To build a PhononBands instance, we usually start from a DdbFile object and then we call one of the anaget_ methods that are essentially wrappers around anaddb.

To get the ph displacement for a single q-point, I would use the following script:

#!/usr/bin/env python

from abipy.abilab import abiopen

with abiopen("trf2_3_DDB") as ddb:
    # ifcflag 0 means that no Fourier interpolation is performed in anaddb
    # hence the q-point must belong to the IBZ used in the DFPT calculation.
    # ifcflag 1 activates the Fourier interpolation of the dynamical matrix.
    # In this case, one can pass an arbitrary q-point to obtain interpolated quantities.

    phbands = ddb.anaget_phmodes_at_qpoint(qpoint=[0.0, 0, 0], ifcflag=0)

    print("Ph frequencies in eV", phbands.phfreqs)

    #phdispl_cart: [nqpt, 3*natom, 3*natom] array with displacement in Cartesian coordinates in Angstrom.
    #    The last dimension stores the cartesian components.
    #    This is an array of complex numbers

    mode = 0
    print(f"phdispl_cart for phonon mode {mode}:\n", phbands.phdispl_cart[0, mode])

[cpashartis]

Okay, cool. Thanks for helping out!