Menu

Dielectric tensor of a slab

Elk Users
Vincent C
2014-04-25
2014-05-19
  • Vincent C

    Vincent C - 2014-04-25

    Dear Elk users,

    Currently I'm trying to calculate the dielectric tensor of a slab of Bi2Te2Se. The calculation of the ground state density and bandstructure went fine but the dielectric tensor is giving a lot of trouble. In the first attempts, I still had spin-orbit coupling switched on but then I got a fortran error that told me I had insufficient memory. (I ran the calculation on a node with 23gb available memory.) As a test, I switched off spin-orbit coupling and then I didn't get any errors from Elk but then the node I used went down. I tried it a couple of times with the same result so I guess it is not just a coincidence.

    Could someone please point me out where things are going wrong? Also, how would I be able to guess up on front how much memory such calculations would take?

    Best regards,
    Vincent

    I'm using elk version 2.2.5. and the input file I used is given below:

    tasks
    120
    121

    wplot
    1000 100 1 :nwplot, ngrkf, nswplot
    0.0 1.0 :wplot

    optcomp
    1 1 1
    2 2 1
    3 3 1

    autoswidth
    .false.

    swidth
    0.01

    stype
    3

    ngridk
    24 24 1

    autolinengy
    .true.

    !Avoids the linearization energies to get lost.
    demaxbnd
    3.0

    mixtype
    3

    broydpm
    0.05 0.01

    !Use PBE
    xctype
    20

    ! Reduce k-mesh with full crystal symmetry group
    reducek
    1

    !Use the high quality parameters
    highq
    .true.

    autoswidth
    .false.

    stype
    3

    autokpt
    .false.

    vkloff
    0 0 0

    ! Angular momentum cut-off for muffin-tin density and potential
    lmaxvr
    9

    ! Spin-orbit coupling
    spinorb
    .false.

    ! Species file are located in the working folder
    sppath
    ' '

    avec
    8.1863927841 0.0000000000 0.0000000000
    -4.0931963921 7.0896241164 0.0000000000
    0.0000000000 0.0000000000 104.4739662636507

    scale
    1.0

    atoms
    3 : nspecies
    'Se.in' : spfname
    4 : natoms
    0.666666667 0.333333333 0.066290915428357 0.0 0.0 0.0 :atposl, bfcmt
    0.333333333 0.666666667 0.248109097008752 0.0 0.0 0.0
    0.000000000 0.000000000 0.429927261680057 0.0 0.0 0.0
    0.666666667 0.333333333 0.611745442714998 0.0 0.0 0.0
    'Bi.in' : spfname
    8 : natoms
    0.333333333 0.666666667 0.032218181472337 0.0 0.0 0.0 : atposl, bfcmt
    0.000000000 0.000000000 0.100363633020742 0.0 0.0 0.0
    0.000000000 0.000000000 0.214036363052732 0.0 0.0 0.0
    0.666666667 0.333333333 0.282181798237501 0.0 0.0 0.0
    0.666666667 0.333333333 0.395854527724037 0.0 0.0 0.0
    0.333333333 0.666666667 0.463999985272442 0.0 0.0 0.0
    0.333333333 0.666666667 0.577672708758977 0.0 0.0 0.0
    0.000000000 0.000000000 0.645818160307382 0.0 0.0 0.0
    'Te.in' : spfname
    8 : natoms
    0.000000000 0.000000000 0.000000000000000 0.0 0.0 0.0 : atposl, bfcmt
    0.333333333 0.666666667 0.132581823220351 0.0 0.0 0.0
    0.666666667 0.333333333 0.181818189216759 0.0 0.0 0.0
    0.000000000 0.000000000 0.314400004255292 0.0 0.0 0.0
    0.333333333 0.666666667 0.363636354433518 0.0 0.0 0.0
    0.666666667 0.333333333 0.496218173290233 0.0 0.0 0.0
    0.000000000 0.000000000 0.545454527286640 0.0 0.0 0.0
    0.333333333 0.666666667 0.678036350506992 0.0 0.0 0.0

     
  • Sangeeta Sharma

    Sangeeta Sharma - 2014-04-25

    you vaccum size (104.47a.u.) is way way way too high. this is an all-electron calculation, you are not using plane-wave code. Please reduce it to something of the order 14a.u. and then run the code again.

     
  • Rob N

    Rob N - 2014-05-02

    As Dr. Sangeeta Sharma pointed out, your vacuum is not calculable. It looks like your top atom to bottom atom distance in the c direction is ~67.5a.u (assuming this is correct for your structure), make the c vector 79.5a.u. and that gives you a vacuum of 12a.u. That should be sufficient.

    You will of course need to recalculate all your z coordinates based on the new basis. Simple ratios should work (i.e. z*104.47/79.5) If you're making the structure in crystal maker you should be able to get it to do this for you.

    Best,

    RN

     
  • Vincent C

    Vincent C - 2014-05-03

    Dear Sangeeta, Rob,

    Thanks for the response.

    I ran the calculations with 14 a.u. vacuum in the unit cell (the c axis is now ~ 84 bohr) this time but the result is still the same. The groundstate calculation goes fine but when I attempt to calculate the dielectric tensor, the node where the code runs on goes down. The code itself does not output any errors.

    I asked one of our system admins to take a look at what exactly goes wrong on the node and currently I'm waiting for reply.

    Best regards,

    Vincent

     

    Last edit: Vincent C 2014-05-03
  • Sangeeta Sharma

    Sangeeta Sharma - 2014-05-04

    Sorry, I do not understand. Is your c-axis 14 or 84? the avec block must look like:
    8.1863927841 0.0000000000 0.0000000000
    -4.0931963921 7.0896241164 0.0000000000
    0.0000000000 0.0000000000 14.00000000000

    What do you mean when you say crashes during "dielectric" calculation. Does it crash at task 120 or 121?

     
  • Vincent C

    Vincent C - 2014-05-05

    Dear Sangeeta,

    The system I am studying is the topological insulator Bi2Te2Se, which is a layered system in and each layer has the ordering Te-Bi-Se-Bi-Te. Before the edge state of the system becomes metallic I need 4 of these quintuple layers, which is why I need the large unit cell. The parameters of the cell I have now are given by (also correcting a mistake I made with the a and b axis but that aside):

    avec
    8.1220397949 0.0000000000 0.0000000000
    -4.0610198975 7.0338927930 0.0000000000
    0.0000000000 0.0000000000 83.940776024272424

    atoms
    3 : nspecies
    'Se.in' : spfname
    4 : natoms
    0.666666667 0.333333333 0.081462644710502 0.0 0.0 0.0 :atposl, bfcmt
    0.333333333 0.666666667 0.304892806027257 0.0 0.0 0.0
    0.000000000 0.000000000 0.528322946565008 0.0 0.0 0.0
    0.666666667 0.333333333 0.751753107211473 0.0 0.0 0.0
    'Bi.in' : spfname
    8 : natoms
    0.333333333 0.666666667 0.039591824212111 0.0 0.0 0.0 : atposl, bfcmt
    0.000000000 0.000000000 0.123333445100179 0.0 0.0 0.0
    0.000000000 0.000000000 0.263021985528867 0.0 0.0 0.0
    0.666666667 0.333333333 0.346763586308220 0.0 0.0 0.0
    0.666666667 0.333333333 0.486452126066618 0.0 0.0 0.0
    0.333333333 0.666666667 0.570193754327880 0.0 0.0 0.0
    0.333333333 0.666666667 0.709882286713083 0.0 0.0 0.0
    0.000000000 0.000000000 0.793623907601150 0.0 0.0 0.0
    'Te.in' : spfname
    8 : natoms
    0.000000000 0.000000000 0.000000000000000 0.0 0.0 0.0 : atposl, bfcmt
    0.333333333 0.666666667 0.162925280036937 0.0 0.0 0.0
    0.666666667 0.333333333 0.223430170700822 0.0 0.0 0.0
    0.000000000 0.000000000 0.386355440683402 0.0 0.0 0.0
    0.333333333 0.666666667 0.446860311908863 0.0 0.0 0.0
    0.666666667 0.333333333 0.609785586583477 0.0 0.0 0.0
    0.000000000 0.000000000 0.670290462500971 0.0 0.0 0.0
    0.333333333 0.666666667 0.833215742537909 0.0 0.0 0.0

    Such that the last atom is located at 69.9407760243 and with 14 bohr of vacuum above the last atom to avoid interaction with the images of the slab.

    When I try to run the calculation of the dielectric tensor, the last output I get from elk is:

    Info(writepmat): 1 of 109 k-points
    (as many of these lines as the amount of processors I execute Elk with)

    So it crashes at task 120. The PMAT.OUT is created but no data are written to it.

    As a test I also attempted to calculate the dielectric response of the same system with 1 and 2 quintuple layers. The first didn't give any problems but the the second has the exact same problem as the larger system above. The cell parameters of this system are given by:

    avec
    8.1220397949 0.0000000000 0.0000000000
    -4.0610198975 7.0338927930 0.0000000000
    0.0000000000 0.0000000000 46.430975512164544

    atoms
    3 : nspecies
    'Se.in' : spfname
    2 : natoms
    0.666666666 0.333333333 0.147273184303388 0.0 0.0 0.0 : atposl, bfcmt
    0.333333333 0.666666666 0.551203985267137 0.0 0.0 0.0
    'Bi.in' : spfname
    4 : natoms
    0.333333333 0.666666666 0.071576537256050 0.0 0.0 0.0 : atposl, bfcmt
    0.000000000 0.000000000 0.222969794996955 0.0 0.0 0.0
    0.000000000 0.000000000 0.475507338219799 0.0 0.0 0.0
    0.666666666 0.333333333 0.626900559606933 0.0 0.0 0.0
    'Te.in' : spfname
    4 : natoms
    0.000000000 0.000000000 0.000000000000000 0.0 0.0 0.0 : atposl, bfcmt
    0.333333333 0.666666666 0.294546351641684 0.0 0.0 0.0
    0.666666666 0.333333333 0.403930817928842 0.0 0.0 0.0
    0.000000000 0.000000000 0.698477151393640 0.0 0.0 0.0

    I am somewhat surprised that calculation of the dielectric tensor turns out to be a problem for these systems. Could you maybe point out why calculating this quantity is a lot more demanding in comparison to the ground state calculation?

    Best regards,
    Vincent

     
  • Sangeeta Sharma

    Sangeeta Sharma - 2014-05-05

    Hi,
    Can you please do the first calculation with

    avec
    8.1220397949 0.0000000000 0.0000000000
    -4.0610198975 7.0338927930 0.0000000000
    0.0000000000 0.0000000000 14.00000000000

    and changing the z of your atomic coordinates accordingly? This should be doable.

     

    Last edit: Sangeeta Sharma 2014-05-05
  • Vincent C

    Vincent C - 2014-05-07

    Hi Sangeeta,

    I am sorry but I am somewhat puzzled by what you ask me to check. I hope I did not miss anything obvious ...

    The system of the first calculation, with 4 quintuple layers, really does not fit in a box with a c-axis of 14 bohr. I tried either way to use the cell parameters you suggested and leave the fractional coordinates fixed. Not unsurprisingly, the ground state calculation ran into problems in this case.

    In fact, I can not even fit a single quintuple layer in a box with a c-axis of 14 bohr since the size of such a layer is ~18,755 bohr in this direction. If I carry out the calculation for 1 quintuple layer and also add 14 bohr of vacuum, the calculation runs fine. The cell parameters are then given by:

    avec
    8.1220397949 0.0000000000 0.0000000000
    -4.0610198975 7.0338927930 0.0000000000
    0.0000000000 0.0000000000 27.676074440272430

    If this is not what you wanted know, could you please clarify what test you do want me to run?

    Kind regards,
    Vincent

     
  • Sangeeta Sharma

    Sangeeta Sharma - 2014-05-09

    It looks like your system is too big (memory wise) for your computers. Hence your caculation keeps crashing. This is clear from you being able to run the same system but smaller in size.

     
  • Vincent C

    Vincent C - 2014-05-09

    Hi Sangeeta,

    Could you maybe point out what parameters determine the memory requirements for the calculation of linear response functions?

    Kind regards,
    Vincent

     
  • Rob N

    Rob N - 2014-05-09

    Vincent,

    Can you check memory usage while the job is running? I strongly suspect you are running out of memory.
    One parameter you might change is the angular momentum cut off
    set lmaxvr = 3.0 to start with and work your way up.

    Also, you should eliminate highq. It is going to redefine parameters you defined above it, such as the kpoint grid.

    I'd start with

    ngridk
    15 15 1
    rgkmax
    7
    lmaxvr
    3.0
    nempty
    8

    reducing lmaxvr should dramatically decrease memory usage. The other paramater will as well to some degree.

    I might also add a shift
    vkloff
    0.25 0.5 0.625

    and reduce nwplot to ~200. 1000 gives you an energy res of 0.02eV. nwplot scales linearly with computational expense (typically time, not memory)

    wplot
    200 100 1 :nwplot, ngrkf, nswplot
    0.0 1.0 :wplot

    I'm assuming you are running this in MPI? If so, you can distribute the work load over more noads to have more memory.

    For example, instead of

    +++++++++++++++++++++++++

    PBS -S /bin/bash
    PBS -l nodes=5:ppn=12
    PBS -l walltime=2:00:00
    PBS -N name
    PBS -j oe

    mpiexec -np 60 ~/elk/src/elk

    +++++++++++++++++++++++++

    try

    +++++++++++++++++++++++++

    PBS -S /bin/bash
    PBS -l nodes=15:ppn=12
    PBS -l walltime=2:00:00
    PBS -N name
    PBS -j oe

    mpiexec -npernode 4 ~/elk/src/elk

    +++++++++++++++++++++++++

    This way you will only have 4 jobs running on each 12 core node, each node then requiring 4/12 the memory required had you run an mpijob for each core.

    I'd try reducing the parameters I discussed above first and see where you can get. If you're still exhausting the memory available then you can distribute the job over more nodes. Also, when you begin increasing these parameters to test convergence, you may need to distribute the jobs as well.

    You definitely want a way to measure the memory usage while you're running these jobs so you will know how to correct accordingly.

    Hope this helps,

    RN

     

    Last edit: Rob N 2014-05-09
  • Vincent C

    Vincent C - 2014-05-14

    Hi Rob,

    Thanks for your reply, it sure helped. I managed to carry out the calculations by lowering the parameters as you suggested.

    Best regards,
    Vincent

     
    • Rob N

      Rob N - 2014-05-19

      Glad I could help. Make sure to convergence test those parameters as well!
      In my experience, a system that size shouldn't require too many kpoints. I suspect 15x15x1 is well converged. If you used the basis set size of 7 (rgkmax), that is also most likely sufficient. I have used rgkmax=5 for several systems, especially those mostly made of light elements, such as carbon.

      I have less experience varying lmaxvr and you may want to rerun the calculation with it equal to 4/5 to be sure 3 was sufficient.

      Also, you may need to check the empty states. For response functions in my systems, 8-12 is usually good. I'm not sure what is necessary for your system.

      Did you turn off "highq"?

      Did you run a spin polarized calc?

      Did you shift the kmesh?

      Also, I have found for broyden mixing the parameters "0.5 0.01" work at least as good and sometimes much better than the default values for broydpm. These parameters shaved off over an order of magnitude from some c60 based calculations.

      Best,

      RN

       

      Last edit: Rob N 2014-05-19

Log in to post a comment.