VASP overshooting Cell Volume and Lattice constants

[sankh]

Dear VASP Masters!

I had been trying to obtain accurate cell volumes and lattice parameters, as well as, internal ion co-ordinates of 48-atoms super-cells of Anatase and Rutile containing substitutive and interstitial dopes. Since some of the doped spices do not have experimental data, precise prognosis of cell geometry (cell vol, parameters and internal co-ords) is of particular importance and a prime aim of my current project. In order to do so, I was trying to calibrate VASP's results for the known systems as anatase and rutiles following the suggestions for cell relaxation at an equilibrium volume. I was using, within LSDA formalism of XC, the PBE PAW potential for core electrons and PAW bases for valence electrons with increasing cutoffs, ranging between 400 to 850 ev. I used experimental cell geometry and followed two approaches as follows:
1.
a. ENCUT = Max(EMax) = 400 ev, ISIF = 2, Cell geom = experimental, reoptimize from CONTCAR ---> There is positive Pulay stress but internal co-ordinates are correct!
b. Read WAVECAR, converged K-points, ENCUT=520, ISIF = 7 (as well as 6), to remove Pulay stress and to obtain equilibrium vol/lattice parameter. It shows that volume converges at 6% higher than experimental value.

2.
a. Obtain equilibrium vol from an EOS fit of the E/V data using several ENCUT (from 400 to 850ev) using ISF = 7 and experimental cell geom input.
b. Reoptimize with ISIF =4 and obtain equilibrium vol, K-points and WAVECAR.

Both approaches produce about 6% higher volume of Anatase and Rutile. Therefore, I cannot trust the equilibrium vols of the doped Anatase/Rutile systems, produced by either of these two approaches. I was also playing around with ISIF = 3 with tight ENCUT=400 and input of experimental geometries. To my surprise, this approach gives excellent cell geometry for Anatase and Rutile at such a tighter cutoff while with higher cutoffs reproduces similar results as 1 and 2 approaches.

Could VASP masters suggest an unambiguous method/strategy with which VASP produces best cell volume (as well as cell geometry) that would always try to match the experimental counterpart as closely as possible? Thanks a lot for this in advance.

Best regards,
Sankh

[forsdan]

To start with: if you use the PAW potentials in the potpaw_PBE directory and do not explicitely set the GGA-tag in the INCAR file you will use the GGA-PBE exchange-correlation functional. You will only use LDA if you have explicitely set GGA = CA. If you want to use use LDA you should however use the potentials in the potpaw directory.

In the case you in fact use GGA-PBE you really should be aware of that GGA in general tend to overestimate the equilibirum volume. Please consult the literature of previous DFT studies of TiO2 to check what is to be expected for this compound.

When comparing your obtained results to experimental data you should be aware of the restrictions of DFT, in particular the impact of the exchange-correlation approximation. Furthermore, you e.g. also need to consider the convergence of the used basis set, temperature effects, and how well the pseudopotentials compare to full potential all electron methods. A discussion of basic considerations can e.g. be found in Modelling Simul. Mater. Sci. Eng. 13 (2005), R1-R31. In addition, the experimental values could been measured on samples which include different levels of impurities. So altogether, there is no garantee at all that the DFT results will accurately represent the experimental data. Consequently, for your study you should therefore use the equilibrium volume as obtained by the functional you are using. If the discrepancies are too large you should investigate if another functional would be more suitable.

Further, neither of your approaches 1 or 2 in general gives a fully consistent combination of shape, volume and internal degrees of freedom (DOF). Unless the internal DOF are identical for all volumes you have two options:

1. Use ISIF = 4 and perform calculations for a several volumes. Then fit a EOS to the values. Please see http://cms.mpi.univie.ac.at/vasp/vasp/A ... s_two.html

2. Use ISIF = 3 in combination with an increased ENCUT (and PREC = Accurate). Please note that 30% is a recommended starting point. You however need to make sure that the evaluation of the stress tensor is sufficiently converged wrt ENCUT to allow for accurate volumes and shapes. The fact that you get a different volume with ISIF=3 and 400 eV than with higher ENCUT is a clear indication that the stress tensor is not sufficiently converged at 400 eV.

For further discussions and considerations, I recommend you to go through the pdf-talks http://cms.mpi.univie.ac.at/vasp-worksh ... tation.htm to learn more.

Cheers,
/Dan






<span class='smallblacktext'>[ Edited Wed Jan 05 2011, 08:47PM ]</span>

[sankh]

Hi Dan,
Thanks a lot of your rapid answers and considerate thoughts which I found very useful.

[quote author= CA. If you want to use use LDA you should however use the potentials in the potpaw directory.

In the case you in fact use GGA-PBE you really should be aware of that GGA in general tend to overestimate the equilibirum volume. Please consult the literature of previous DFT studies of TiO2 to check what is to be expected for this compound.
[/quote]Your these two paragraphs confuses me a bit, because I would presume PAW pseudopotentials are to treat/replace the core electrons (an add-on to 1-particle part of KS hamiltonian by replacing the core electrons with an average potential), while any chosen XC functional describes the interaction among the valence electrons spanned over PAW bases. Certainly, a type of pseudopotentials for core would be more suitable over the others for a particular XC functional for valence electrons, and this is subject to a sizable research work, which I presume that VASP developers have already done. If I am not very incorrect (if I am, please correct me), PAW_PBE pseudopotentials were created over CI reference states and that Prof. Kresse in a conference encouraged to use them due to their versatility over the ones created on LDA/LSDA or GGA reference states. VASP 5.2 pdf manual cites that one has to set GGA tag explicitly in INCAR in order to obtain GGA correction to the chosen XC functionals of interest, irrespective of chosen pseudopotentials for the core; otherwise VASP switches to LDA/LSDA. Further VASP manual does not mention GGA=CA option in the discussion of the available GGA tags. Therefore, I got to ask the following questions:

1. Which set of PAW pseudopotentials (e.g. potpaw, paw_gga, paw_pbe) should be used in combination with which set of XC functionals (say LDA/LSDA [CA of Perdrew-Zugang], GGA types [PW91, revised PBE, PBEsol, AM05])?

2. Could PAW_PBE pseudopotentials be used in combination with AM05 GGA functional?

3. Which set of pseudopotentials should be used in combination with screened hybrid functionals of type a) LDA/LSDA + Fock, b) PBE + Fock, c) PBEsol + Fock, d) AM05+Fock ?
[/color]

[quote author= 4 and perform calculations for a several volumes. Then fit a EOS to the values. Please see http://cms.mpi.univie.ac.at/vasp/vasp/A ... s_two.html
[/quote]Thank you for your suggestion, but I beg to disagree with you, as this "1" approach is not very practical for an unknown lattice system especially wherein a cell expansion (as well as a dramatic change in lattice shape may occur) is expected due to a defect intruder. One has no idea about the guesses for volumes and lattice parameters to try out, as they might be strongly coupled with internal (ion) co-ordinates, giving rise to an overall very complicated potential energy hypersurface with tons of local minima. But this approach is excellent, as so far I found out, when the lattice geometry is close by or predictable. [/color]

[quote author= 3 in combination with an increased ENCUT (and PREC = Accurate). Please note that 30% is a recommended starting point. You however need to make sure that the evaluation of the stress tensor is sufficiently converged wrt ENCUT to allow for accurate volumes and shapes. The fact that you get a different volume with ISIF=3 and 400 eV than with higher ENCUT is a clear indication that the stress tensor is not sufficiently converged at 400 eV.
[/quote]Yes, I was aware of insufficient convergence in stress tensor at 400 ev cutoff. But what really surprised me is that with ISIF = 3, Prec=Acc, the higher cutoffs though removed the Pulay stresses, which is theoretically obvious, produced rather incorrect lattice and internal co-ordinates for known systems, as TiO2 anatase and rutile, while the lower cutoff 400ev produced the geometry in excellent agreement with the experimental as well as the previous theoretical data, albeit with significant Pulay stresses. Although I am not aware of the numerical details associated with ISIF=3, 4 and 2, I would expect from mathematical point of view that internal coordinates are the major contributers toward administering the lattice geometry rather than the numerical issues such as Pulay's stress in finite bases. However, I will discover more when I have sorted out the discrepancies associated use of pseudopotentials and XC funtionals once you have answered my above questions.

I am quite hopeful about an approach comprising ISIF = 2 followed by ISIF=3 for unknown system.

I am wondering what you think about rather constant density based approach than present constant volume based approaches towards geometry optimizations. Could VASP do constant density calculations at present?
[/color]

[quote="blue"]Thank you very much for the support and I look forward to hearing from you at your earliest convenience.

Best regards,
Sankh[/color]
<span class='smallblacktext'>[ Edited Fri Jan 07 2011, 12:31AM ]</span>

[admin]

please let me add one short note: please remind that 1 dopant atom in a 48-atoms cell refers to a macroscopic dopant concentration of more than 2%(at/at) due to the periodic boundary conditions. What you expect to be the equilibrium geometry of such a cell might have nothing to do with the actual probes. Please have a look whether you find any data on typical dopant concentrations of the systems you want to investigate. If they are significantly smaller, it is more reliable to choose the equilibrium structure of the undoped system for further studies, in order to avoid unreasonable shape/volume effects.

[sankh]

Deat Head Admin,

Thank you very much for recent message. I truely feel grateful to have had your response to my issue.

Our experimental collaborators report that they see notable amount of cell (well "grain size", in their terminology) at different temperatures and different methods of heat-treatment, through analyzing XRD data. The dopant concentration lies within 3-8% atomic concentrations. I am trying to theoretically model thier outcome, as one of the principal aims within this project, which has some novel features in the field. Eventually and luckily, my theoretical data, obtained through your great code VASP are in excellent agreement with our experimental observation, including the size expansion of the doped nano-crystallites at a particular set of conditions, and we are in a process to report these in a publication in very near future.

However, I will be very grateful to have a word with Prof. Kresse to discuss the underlying issues, if he kindly permits me, although I am well aware of the fact that he is overwhelmed with the research works at his hand.

I look forward to hear from you. And sincere thanks for your time for my sake.

Best regards,
Sankh