[Wien] Finite temperature DFT: electronic free energy

[Sabry Moustafa]

(0) PROBLEM:
WIEN2k did not accept lage temperatures (>~4000 K). I always get this 
error "FERMI5 - Error" in the job error file. No complain in the 
init_lapw step; only after run_lapw finishes. Any help here?


(1) Thanks Blaha and Gavin for your help -- I tried it and worked fine. 
Gavin, please note that in Fermi-Dirac statistics the energy is 
temperature dependent. See for example Marzari's thesis (e.g., Fig 4.2):
http://theossrv1.epfl.ch/Main/Theses?action=download&upname=Marzari_thesis_1996.pdf 



(2) Actually, the "TOTAL ENERGY IN Ry" in case.scf is the free energy 
F(T); i.e., F(T)=E(T) - T*S(T). This is written explicitly in the UG, 
when talking about "TEMPS" (Sec. 7.8.3):
"The total energy is corrected by -TS ..."
So, E(T) can be obtained from F(T) and -T*S(T), i.e., E = F + TS.

(3) Since I am beginner and familiar more with VASP at this point, I 
checked (using hcp-Fe crystal) both F(T) and -T*S(T) from both and found 
them to be consistent!


Thanks;
Sabry

On 05/15/2018 03:34 AM, Peter Blaha wrote:
>
> Thanks for the bug report in init_lapw.
>
> Please note: when the temperature (broadening) is set to zero, it will 
> use automatically "room temperature" (0.0018 Ry).
>
> On 05/15/2018 06:40 AM, Gavin Abo wrote:
>> I think you found a bug in init_lapw.  The fix seems like it should 
>> be simple though.  On line 498 in the init_lapw script in the sed 
>> command, it looks like $fermit just needs changed to $fermits.
>>
>> I made a patch file called init_lapw.patch [ 
>> https://github.com/gsabo/WIEN2k-Patches/tree/master/17.1 ], which you 
>> should be able to apply and try simply using:
>>
>> username at computername:~/Desktop$ cd $WIENROOT
>> username at computername:~/WIEN2k$ wget 
>> https://raw.githubusercontent.com/gsabo/WIEN2k-Patches/master/17.1/init_lapw.patch
>> username at computername:~/WIEN2k$ patch -b init_lapw init_lapw.patch
>> patching file init_lapw
>>
>> Maybe F is the :ENE (TOTAL ENERGY) in case.scf.
>>
>> I tried following what a previous user did [ 
>> https://www.mail-archive.com/wien@zeus.theochem.tuwien.ac.at/msg10319.html 
>> ] but with a quick TiC calculation using run_lapw in WIEN2k 17.1:
>>
>> T = 1000 K = 0.00633 Ry
>> init_lapw -b -fermits 0.00633
>> TiC.scf::ENE  : ********** TOTAL ENERGY IN Ry = -1783.95041939
>> TiC.output2:          Kb * T            =   0.00633000
>> TiC.output2:          -T*Entr           =  -0.00054181
>>
>> T = 0 K = 0.0 Ry
>> init_lapw -b -fermits 0.0
>> TiC.scf::ENE  : ********** TOTAL ENERGY IN Ry = -1783.94928338
>> TiC.output2:          Kb * T            =   0.00180000
>> TiC.output2:          -T*Entr           =  -0.00005329
>>
>> F_1000K - F_0K = (E - T*S) - (E - 0) = -T*S = -1783.95041939 - 
>> (-1783.94928338) = -0.001136 => -T*S = -0.001136
>>
>> If -T*S divided by 2:
>>
>> -T*S = -0.001136/2 = - 0.000568 <- This seems reasonably close to the 
>> above -0.00054181 for T = 1000 K.
>>
>> As I recall, the -T*S term was doubled in WIEN2k versions after 2014 
>> [ 
>> https://www.mail-archive.com/wien@zeus.theochem.tuwien.ac.at/msg10319.html 
>> ].
>>
>> Though, there may be a flaw in my above calculation.
>>
>> On 5/14/2018 5:35 PM, Sabry Moustafa wrote:
>>>
>>> Hi;
>>>
>>> I am interested in the electronic free energy (F=E-TS) using the 
>>> finite temperature DFT using Fermi-Dirac statistics (i.e., Mermin's 
>>> extension to 0 K DFT) -- i.e., I am not just looking for Fermi-Dirac 
>>> as a broadening/smearing technique. So,how this can be done in WIEN2k?
>>>
>>> It looks like I have to define "-fermits X" in the init_lapw command 
>>> (where X= kT, in Ry). But when I did (T=1000K in my case):
>>>
>>> init_lapw -b -fermits 0.00633
>>>
>>> I got "fermit: Undefined variable." at the end. This is fixed though 
>>> when defined fermit as well:
>>>
>>> init_lapw -b -fermit 0.00633 -fermits 0.00633
>>>
>>> So, do I have to define both?
>>>
>>>
>>> Also, where can I find this "free energy F". It does not seem to be 
>>> the "TOTAL ENERGY" in case.scf. Seems like I need to add this "TOTAL 
>>> ENERGY" (E) to "-T*Entr" (-TS) in case.output2 file to get F?
>>>
>>>
>>> Thanks ;
>>>
>>> Sabry
>>>
>>>
>>> :-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-:-: 
>>>
>>>
>>> Sabry G. Moustafa, Ph.D.
>>> Department of Chemical and Biological Engineering,
>>> University at Buffalo, The State University of New York.
>>> 511 Furnas Hall
>>> Buffalo, NY 14260-4200
>>> 716-645-1186 (office)
>>> 716-239-8543 (cell)
>>>
>>> E-mail: sabrygad at buffalo.edu <mailto:sabrygad at buffalo.edu>
>>
>>
>> _______________________________________________
>> Wien mailing list
>> Wien at zeus.theochem.tuwien.ac.at
>> http://zeus.theochem.tuwien.ac.at/mailman/listinfo/wien
>> SEARCH the MAILING-LIST at: 
>> http://www.mail-archive.com/wien@zeus.theochem.tuwien.ac.at/index.html
>>
>